| Revision: | 1.1 |
| Committed: | Sun Apr 28 07:46:58 2013 UTC (10 years, 4 months ago) by tplarson |
| Content type: | text/plain |
| Branch: | MAIN |
| CVS Tags: | globalhs_version_5, Ver_8-7, Ver_8-5, globalhs_version_23, globalhs_version_22, globalhs_version_21, globalhs_version_20, Ver_LATEST, globalhs_version_24, Ver_8-3, globalhs_version_8, globalhs_version_9, globalhs_version_0, globalhs_version_1, globalhs_version_2, globalhs_version_3, globalhs_version_4, Ver_9-41, globalhs_version_6, globalhs_version_7, Ver_9-5, Ver_8-8, globalhs_version_19, Ver_8-2, Ver_8-10, Ver_8-1, Ver_8-6, Ver_9-1, Ver_8-4, Ver_9-2, globalhs_version_12, globalhs_version_13, globalhs_version_10, globalhs_version_11, globalhs_version_16, globalhs_version_17, globalhs_version_14, globalhs_version_15, globalhs_version_18, Ver_9-4, Ver_9-3, Ver_8-11, Ver_8-12, Ver_9-0, HEAD |
| Log Message: | functions needed for detrending and gapfilling |
| # | Content |
|---|---|
| 1 | // $Header: $ |
| 2 | |
| 3 | #include <stdio.h> |
| 4 | #include <stdlib.h> |
| 5 | #include <string.h> |
| 6 | #include <complex.h> |
| 7 | #include <float.h> |
| 8 | #include <math.h> |
| 9 | #include <assert.h> |
| 10 | #include <time.h> |
| 11 | |
| 12 | #include "ctypes.h" |
| 13 | #include "fahlman_ulrych_C99.h" |
| 14 | #include "multi_burg_C99.h" |
| 15 | #include "pcg_C99.h" |
| 16 | #include "levinson_C99.h" |
| 17 | //#include "xmem.h" |
| 18 | |
| 19 | |
| 20 | static inline int max(int a, int b) { return (a>b? a : b); } |
| 21 | static inline int min(int a, int b) { return (a<b? a : b); } |
| 22 | |
| 23 | typedef struct interval_struct |
| 24 | { int first; |
| 25 | int last; |
| 26 | } interval; |
| 27 | |
| 28 | typedef struct gapped_timeseries_struct |
| 29 | { |
| 30 | int first_is_good; |
| 31 | int n_data, m_data; |
| 32 | int n_gap, m_gap; |
| 33 | interval *gap_int, *data_int; |
| 34 | } gapped_timeseries; |
| 35 | |
| 36 | |
| 37 | static int maxorder(gapped_timeseries *ts, int *data_length, |
| 38 | int minpercentage) |
| 39 | { |
| 40 | int i,j; |
| 41 | int *len,idx,tmp,order; |
| 42 | // Sort gaps by length and determine the largest AR model order |
| 43 | // that includes minpercentage percent of the data to be used |
| 44 | // in the estimation of the AR coefficients. |
| 45 | len = malloc(ts->m_data*sizeof(int)); |
| 46 | memcpy(len, data_length, ts->m_data*sizeof(int)); |
| 47 | /* for (i=0; i<ts->m_data; i++) |
| 48 | printf("len[%3d] = [%d:%d] = %5d\n",i,ts->data_int[i].first, |
| 49 | ts->data_int[i].last,len[i]); */ |
| 50 | for (i=0; i<ts->m_data; i++) |
| 51 | { |
| 52 | idx = i; |
| 53 | for (j=i+1; j<ts->m_data; j++) |
| 54 | if (len[j]>len[idx]) |
| 55 | idx = j; |
| 56 | tmp = len[idx]; |
| 57 | len[idx] = len[i]; |
| 58 | len[i] = tmp; |
| 59 | } |
| 60 | #ifdef DEBUGOUT |
| 61 | { |
| 62 | for (i=0; i<ts->m_data; i++) |
| 63 | printf("%5d ",len[i]); |
| 64 | printf("\n"); |
| 65 | } |
| 66 | #endif |
| 67 | tmp = 0; i = 0; |
| 68 | while(i<ts->n_data && tmp<(minpercentage*ts->n_data)/100) |
| 69 | { |
| 70 | tmp += len[i++]; |
| 71 | } |
| 72 | order = len[i-1]-1; |
| 73 | #ifdef DEBUGOUT |
| 74 | { |
| 75 | printf("maximal initial order chosen = %d\n",order); |
| 76 | printf("using %d out of %d data points\n",tmp,ts->n_data); |
| 77 | } |
| 78 | #endif |
| 79 | free(len); |
| 80 | return order; |
| 81 | } |
| 82 | |
| 83 | |
| 84 | |
| 85 | |
| 86 | #ifdef FLOAT |
| 87 | #define TYPE FLOAT |
| 88 | #include "fahlman_ulrych_code_C99.h" |
| 89 | #undef TYPE |
| 90 | #endif |
| 91 | |
| 92 | #ifdef DOUBLE |
| 93 | #define TYPE DOUBLE |
| 94 | #include "fahlman_ulrych_code_C99.h" |
| 95 | #undef TYPE |
| 96 | #endif |
| 97 | |
| 98 | #ifdef COMPLEXFLOAT |
| 99 | #define TYPE COMPLEXFLOAT |
| 100 | #include "fahlman_ulrych_code_C99.h" |
| 101 | #undef TYPE |
| 102 | #endif |
| 103 | |
| 104 | #ifdef COMPLEXDOUBLE |
| 105 | #define TYPE COMPLEXDOUBLE |
| 106 | #include "fahlman_ulrych_code_C99.h" |
| 107 | #undef TYPE |
| 108 | #endif |
| 109 | |
| 110 |