1 /*================================================================== 2 * evresp.h 3 *=================================================================*/ 4 5 /* 6 8/28/2001 -- [ET] Added 'WIN32' directives for Windows compiler 7 compatibility; added 'extern' to variable 8 declarations at end of file. 9 8/2/2001 -- [ET] Version 3.2.28: Modified to allow command-line 10 parameters for frequency values to be missing 11 (default values used). 12 10/21/2005 -- [ET] Version 3.2.29: Implemented interpolation of 13 amplitude/phase values from responses containing 14 List blockettes (55); modified message shown when 15 List blockette encountered; modified so as not to 16 require characters after 'units' specifiers like 17 "M" and "COUNTS"; modified to handle case where 18 file contains "B052F03 Location:" and nothing 19 after it on line; added warnings for unrecognized 20 parameters; fixed issue where program would crash 21 under Linux/UNIX ("segmentation fault") if input 22 response file contained Windows-style "CR/LF" 23 line ends; renamed 'strncasecmp()' function 24 to 'strnicmp()' when compiling under Windows. 25 11/3/2005 -- [ET] Version 3.2.30: Modified to unwrap phase values 26 before interpolation and re-wrap after (if needed); 27 modified to detect when all requested interpolated 28 frequencies are out of range. 29 1/18/2006 -- [ET] Version 3.2.31: Renamed 'regexp' functions to 30 prevent name clashes with other libraries. 31 2/13/2006 -- [ET] Version 3.2.32: Moved 'use_delay()' function from 32 'evalresp.c' to 'evresp.c'; modified to close input 33 file when a single response file is specified. 34 3/27/2006 -- [ID] Version 3.2.33: Added include_HEADERS target 35 "evr_spline.h" to "Makefile.am"; upgraded "missing" 36 script. 37 3/28/2006 -- [ET] Version 3.2.34: Added 'free()' calls to "evresp.c" 38 to fix memory leaks when program functions used 39 in external applications. 40 4/4/2006 -- [ET] Version 3.2.35: Modified to support channel-IDs with 41 location codes equal to an empty string. 42 8/21/2006 -- [IGD] Version 3.2.36: Added support for TESLA units. 43 10/16/2006 -- [ET] Version 3.2.37: Added two more 'free()' calls to 44 "evresp.c" to fix memory leaks when program functions 45 are used in external applications. 46 8/11/2006 -- [IGD] Libtoolized evresp library 47 4/03/2007 -- [IGD] Added myLabel global variable which is used to add NSLC 48 labels in evalresp logging if --enable-log-label config 49 option is used 50 5/14/2010 -- [ET] Version 3.3.3: Added "#define strcasecmp stricmp" 51 if Windows. 52 */ 53 54 #ifndef EVRESP_H 55 #define EVRESP_H 56 57 #include <stdio.h> 58 #include <math.h> 59 #include <ctype.h> 60 #include <setjmp.h> 61 #include <stdarg.h> 62 63 /* IGD 10/16/04 This is for Windows which does not use Makefile.am */ 64 #ifdef VERSION 65 #define REVNUM VERSION 66 #else 67 #define REVNUM "3.3.3" 68 #endif 69 70 #define TRUE 1 71 #define FALSE 0 72 #define QUERY_DELAY -1 73 #define STALEN 64 74 #define NETLEN 64 75 #define LOCIDLEN 64 76 #define CHALEN 64 77 #define OUTPUTLEN 256 78 #define TMPSTRLEN 64 79 #define UNITS_STR_LEN 16 80 #define DATIMLEN 23 81 #define UNITSLEN 20 82 #define BLKTSTRLEN 4 83 #define FLDSTRLEN 3 84 #define MAXFLDLEN 50 85 #define MAXLINELEN 256 86 #define FIR_NORM_TOL 0.02 87 88 89 /* enumeration representing the types of units encountered (Note: if default, 90 then the response is just given in input units to output units, no 91 interpretation is made of the units used) */ 92 93 /* IGD 02/03/01 New unit pressure added */ 94 /* IGD 08/21/06 New units TESLA added */ 95 enum units { UNDEF_UNITS, DIS, VEL, ACC, COUNTS, VOLTS, DEFAULT, PRESSURE, TESLA 96 }; 97 98 /* enumeration representing the types of filters encountered */ 99 100 enum filt_types { UNDEF_FILT, LAPLACE_PZ, ANALOG_PZ, IIR_PZ, 101 FIR_SYM_1, FIR_SYM_2, FIR_ASYM, LIST, GENERIC, DECIMATION, 102 GAIN, REFERENCE, FIR_COEFFS, IIR_COEFFS 103 }; 104 105 /* enumeration representing the types of stages that are recognized */ 106 /* IGD 05/15/02 Added GENERIC_TYPE */ 107 enum stage_types { UNDEF_STAGE, PZ_TYPE, IIR_TYPE, FIR_TYPE, 108 GAIN_TYPE, LIST_TYPE, IIR_COEFFS_TYPE, GENERIC_TYPE 109 }; 110 111 /* enumeration representing the types of error codes possible */ 112 113 enum error_codes { NON_EXIST_FLD = -2, ILLEGAL_RESP_FORMAT = -5, 114 PARSE_ERROR = -4, UNDEF_PREFIX = -3, UNDEF_SEPSTR = -6, 115 OUT_OF_MEMORY = -1, UNRECOG_FILTYPE = -7, 116 UNEXPECTED_EOF = -8, ARRAY_BOUNDS_EXCEEDED = -9, 117 OPEN_FILE_ERROR = 2, RE_COMP_FAILED = 3, 118 MERGE_ERROR = 4, SWAP_FAILED = 5, USAGE_ERROR = 6, 119 BAD_OUT_UNITS = 7, IMPROP_DATA_TYPE = -10, 120 UNSUPPORT_FILTYPE = -11, ILLEGAL_FILT_SPEC = -12, 121 NO_STAGE_MATCHED = -13, UNRECOG_UNITS = -14 122 }; 123 124 /* define structures for the compound data types used in evalesp */ 125 126 /* if Windows compiler then redefine 'complex' to */ 127 /* differentiate it from the existing struct, */ 128 /* and rename 'strcasecmp' functions: */ 129 #ifdef WIN32 130 #ifdef complex 131 #undef complex 132 #endif 133 #define complex evr_complex 134 #define strcasecmp stricmp 135 #define strncasecmp strnicmp 136 #endif 137 138 struct complex { 139 double real; 140 double imag; 141 }; 142 143 struct string_array { 144 int nstrings; 145 char **strings; 146 }; 147 148 struct scn { 149 char *station; 150 char *network; 151 char *locid; 152 char *channel; 153 int found; 154 }; 155 156 struct response { 157 char station[STALEN]; 158 char network[NETLEN]; 159 char locid[LOCIDLEN]; 160 char channel[CHALEN]; 161 struct complex *rvec; 162 int nfreqs; /*Add by I.Dricker IGD to support blockette 55 */ 163 double *freqs; /*Add by I.Dricker IGD to support blockette 55 */ 164 struct response *next; 165 }; 166 167 struct file_list { 168 char *name; 169 struct file_list *next_file; 170 }; 171 172 struct matched_files { 173 int nfiles; 174 struct file_list *first_list; 175 struct matched_files *ptr_next; 176 }; 177 178 struct scn_list { 179 int nscn; 180 struct scn **scn_vec; 181 }; 182 183 /* define structures for the various types of filters defined in seed */ 184 185 struct pole_zeroType { /* a Response (Poles & Zeros) blockette */ 186 int nzeros; /* (blockettes [43] or [53]) */ 187 int npoles; 188 double a0; 189 double a0_freq; 190 struct complex *zeros; 191 struct complex *poles; 192 }; 193 194 struct coeffType { /* a Response (Coefficients) blockette */ 195 int nnumer; /* (blockettes [44] or [54]) */ 196 int ndenom; 197 double *numer; 198 double *denom; 199 double h0; /*IGD this field is new v 3.2.17 */ 200 }; 201 202 struct firType { /* a FIR Response blockette */ 203 int ncoeffs; /* (blockettes [41] or [61])*/ 204 double *coeffs; 205 double h0; 206 }; 207 208 struct listType{ /* a Response (List) blockette */ 209 int nresp; /* (blockettes [45] or [55]) */ 210 double *freq; 211 double *amp; 212 double *phase; 213 }; 214 215 struct genericType { /* a Generic Response blockette */ 216 int ncorners; /* (blockettes [46] or [56]) */ 217 double *corner_freq; 218 double *corner_slope; 219 }; 220 221 struct decimationType { /* a Decimation blockette */ 222 double sample_int; /* (blockettes [47] or [57]) */ 223 int deci_fact; 224 int deci_offset; 225 double estim_delay; 226 double applied_corr; 227 }; 228 229 struct gainType { /* a Channel Sensitivity/Gain blockette */ 230 double gain; /* (blockettes [48] or [58]) */ 231 double gain_freq; 232 }; 233 234 struct referType { /* a Response Reference blockette */ 235 int num_stages; 236 int stage_num; 237 int num_responses; 238 }; 239 240 /* define a blkt as a stucture containing the blockette type, a union 241 (blkt_info) containing the blockette info, and a pointer to the next 242 blockette in the filter sequence. The structures will be assembled to 243 form a linked list of blockettes that make up a filter, with the last 244 blockette containing a '(struct blkt *)NULL' pointer in the 'next_blkt' 245 position */ 246 247 struct blkt { 248 int type; 249 union { 250 struct pole_zeroType pole_zero; 251 struct coeffType coeff; 252 struct firType fir; 253 struct listType list; 254 struct genericType generic; 255 struct decimationType decimation; 256 struct gainType gain; 257 struct referType reference; 258 } blkt_info; 259 struct blkt *next_blkt; 260 }; 261 262 /* define a stage as a structure that contains the sequence number, the 263 input and output units, a pointer to the first blockette of the filter, 264 and a pointer to the next stage in the response. Again, the last 265 stage in the response will be indicated by a '(struct stage *)NULL pointer 266 in the 'next_stage' position */ 267 268 struct stage { 269 int sequence_no; 270 int input_units; 271 int output_units; 272 struct blkt *first_blkt; 273 struct stage *next_stage; 274 }; 275 276 /* and define a channel as a stucture containing a pointer to the head of a 277 linked list of stages. Will access the pieces one stages at a time in the 278 same order that they were read from the input file, so a linked list is the 279 easiest way to do this (since we don't have to worry about the time penalty 280 inherent in following the chain). As an example, if the first stage read 281 was a pole-zero stage, then the parts of that stage contained in a channel 282 structure called "ch" would be accessed as: 283 284 struct stage *stage_ptr; 285 struct blkt *blkt_ptr; 286 287 stage_ptr = ch.first_stage; 288 blkt_ptr = stage_ptr->first_blkt; 289 if(blkt_ptr->type == LAPLACE_PZ || blkt_ptr->type == ANALOG_PZ || 290 blkt_ptr->type == IIR_PZ){ 291 nzeros = blkt_ptr->blkt_info.poles_zeros.nzeros; 292 ........ 293 } 294 295 */ 296 297 struct channel { 298 char staname[STALEN]; 299 char network[NETLEN]; 300 char locid[LOCIDLEN]; 301 char chaname[CHALEN]; 302 char beg_t[DATIMLEN]; 303 char end_t[DATIMLEN]; 304 char first_units[MAXLINELEN]; 305 char last_units[MAXLINELEN]; 306 double sensit; 307 double sensfreq; 308 double calc_sensit; 309 double calc_delay; 310 double estim_delay; 311 double applied_corr; 312 double sint; 313 int nstages; 314 struct stage *first_stage; 315 }; 316 317 /* structure used for time comparisons */ 318 319 struct dateTime { 320 int year; 321 int jday; 322 int hour; 323 int min; 324 float sec; 325 }; 326 327 /* IGD 2007/02/27 */ 328 extern char myLabel[20]; 329 //char myLabel[20] = "aa"; 330 331 /* utility routines that are used to parse the input file line by line and 332 convert the input to what the user wants */ 333 334 struct string_array *ev_parse_line(char *); 335 struct string_array *parse_delim_line(char *, char *); 336 int get_field(FILE *, char *, int, int, char *, int); 337 int test_field(FILE *, char *, int *, int *, char *, int); 338 int get_line(FILE *, char *, int, int, char *); /* checks blkt & fld nos */ 339 int next_line(FILE *, char *, int *, int *, char *); /* returns blkt & fld nos */ 340 int count_fields(char *); 341 int count_delim_fields(char *, char *); 342 int parse_field(char *, int, char *); 343 int parse_delim_field(char *, int, char *, char *); 344 int check_line(FILE *, int *, int *, char *); 345 int get_int(char *); 346 double get_double(char *); 347 int check_units(char *); 348 int string_match(const char *, char *, char *); 349 int is_int(const char *); 350 int is_real(const char *); 351 int is_IIR_coeffs (FILE *, int); /*IGD */ 352 353 /* routines used to load a channel's response information into a linked 354 list of filter stages, each containing a linked list of blockettes */ 355 356 int find_resp(FILE *, struct scn_list *, char *, struct channel *); 357 int get_resp(FILE *, struct scn *, char *, struct channel *); 358 int get_channel(FILE *, struct channel *); 359 int next_resp(FILE *); 360 361 /* routines used to create a list of files matching the users request */ 362 363 struct matched_files *find_files(char *, struct scn_list *, int *); 364 int get_names(char *, struct matched_files *); 365 int start_child(char *, FILE **, FILE **, FILE **); 366 367 /* routines used to allocate vectors of the basic data types used in the 368 filter stages */ 369 370 struct complex *alloc_complex(int); 371 struct response *alloc_response(int); 372 struct string_array *alloc_string_array(int); 373 struct scn *alloc_scn(void); 374 struct scn_list *alloc_scn_list(int); 375 struct file_list *alloc_file_list(void); 376 struct matched_files *alloc_matched_files(void); 377 double *alloc_double(int); 378 char *alloc_char(int); 379 char **alloc_char_ptr(int); 380 381 /* allocation routines for the various types of filters */ 382 383 struct blkt *alloc_pz(void); 384 struct blkt *alloc_coeff(void); 385 struct blkt *alloc_fir(void); 386 struct blkt *alloc_ref(void); 387 struct blkt *alloc_gain(void); 388 struct blkt *alloc_list(void); 389 struct blkt *alloc_generic(void); 390 struct blkt *alloc_deci(void); 391 struct stage *alloc_stage(void); 392 393 /* routines to free up space associated with dynamically allocated 394 structure members */ 395 396 void free_string_array(struct string_array *); 397 void free_scn(struct scn *); 398 void free_scn_list(struct scn_list *); 399 void free_matched_files(struct matched_files *); 400 void free_file_list(struct file_list *); 401 void free_pz(struct blkt *); 402 void free_coeff(struct blkt *); 403 void free_fir(struct blkt *); 404 void free_list(struct blkt *); 405 void free_generic(struct blkt *); 406 void free_gain(struct blkt *); 407 void free_deci(struct blkt *); 408 void free_ref(struct blkt *); 409 void free_stages(struct stage *); 410 void free_channel(struct channel *); 411 void free_response(struct response *); 412 413 /* simple error handling routines to standardize the output error values and 414 allow for control to return to 'evresp' if a recoverable error occurs */ 415 416 void error_exit(int, char *, ...); 417 void error_return(int, char *, ...); 418 419 /* a simple routine that parses the station information from the input file */ 420 421 int parse_channel(FILE *, struct channel *); 422 423 /* parsing routines for various types of filters */ 424 425 int parse_pref(int *, int *, char *); 426 void parse_pz(FILE *, struct blkt *, struct stage *); /* pole-zero */ 427 void parse_coeff(FILE *, struct blkt *, struct stage *); /* fir */ 428 void parse_iir_coeff(FILE *, struct blkt *, struct stage *); /*I.Dricker IGD for 2.3.17 iir */ 429 void parse_list(FILE *, struct blkt *, struct stage *); /* list */ 430 void parse_generic(FILE *, struct blkt *, struct stage *); /* generic */ 431 int parse_deci(FILE *, struct blkt *); /* decimation */ 432 int parse_gain(FILE *, struct blkt *); /* gain/sensitivity */ 433 void parse_fir(FILE *, struct blkt *, struct stage *); /* fir */ 434 void parse_ref(FILE *, struct blkt *, struct stage *); /* response reference */ 435 436 /* remove trailing white space from (if last arg is 'a') and add null character to 437 end of (if last arg is 'a' or 'e') input (FORTRAN) string */ 438 439 int add_null(char *, int, char); 440 441 /* run a sanity check on the channel's filter sequence */ 442 443 void merge_coeffs(struct blkt *, struct blkt **); 444 void merge_lists(struct blkt *, struct blkt **); /* Added by I.Dricker IGD for v 3.2.17 of evalresp */ 445 void check_channel(struct channel *); 446 void check_sym(struct blkt *, struct channel *); 447 448 /* routines used to calculate the instrument responses */ 449 450 /*void calc_resp(struct channel *, double *, int, struct complex *,char *, int, int);*/ 451 void calc_resp(struct channel *chan, double *freq, int nfreqs, struct complex *output, 452 char *out_units, int start_stage, int stop_stage, int useTotalSensitivityFlag); 453 void convert_to_units(int, char *, struct complex *, double); 454 void analog_trans(struct blkt *, double, struct complex *); 455 void fir_sym_trans(struct blkt *, double, struct complex *); 456 void fir_asym_trans(struct blkt *, double, struct complex *); 457 void iir_pz_trans(struct blkt *, double, struct complex *); 458 void calc_time_shift(double, double, struct complex *); 459 void zmul(struct complex *, struct complex *); 460 void norm_resp(struct channel *, int, int, int); 461 void calc_list(struct blkt *, int , struct complex *); /*IGD i.dricker@isti.com for version 3.2.17 */ 462 void iir_trans(struct blkt *, double , struct complex *); /* IGD for version 3.2.17 */ 463 int is_time (const char * ); 464 465 /* compare two times and determine if the first is greater, equal to, or less than the second */ 466 467 int timecmp(struct dateTime *dt1, struct dateTime *dt2); 468 469 /* print the channel info, followed by the list of filters */ 470 471 void print_chan(struct channel *, int, int, int, int, int, int); 472 473 /* print the response information to the output files */ 474 475 void print_resp(double *, int, struct response *, char *, int); 476 void print_resp_itp(double *, int, struct response *, char *, int, 477 int, double, int); 478 479 /* evaluate responses for user requested station/channel/network tuple at the 480 frequencies requested by the user */ 481 482 struct response *evresp(char *, char *, char *, char *, char *, char *, char *, 483 double *, int, char *, char *, int, int, int, int); 484 struct response *evresp_itp(char *, char *, char *, char *, char *, char *, 485 char *, double *, int, char *, char *, int, 486 int, int, int, int, double, int); 487 488 /* Interpolates amplitude and phase values from the set of frequencies 489 in the List blockette to the requested set of frequencies. */ 490 void interpolate_list_blockette(double **freq, double **amp, double **phase, 491 int *p_number_points, double *req_freq_arr, 492 int req_num_freqs, double tension); 493 494 495 extern char SEEDUNITS[][UNITS_STR_LEN]; 496 497 /* define the "first line" and "first field" arguments that are used by the 498 parsing routines (they are used to compensate for the fact that in 499 parsing the RESP files, one extra line is always read because the 500 "end" of a filter sequence cannot be determined until the first line 501 of the "next" filter sequence or station-channel pair is read from the 502 file */ 503 504 extern char FirstLine[]; 505 extern int FirstField; 506 507 /* values for use in calculating responses (defined in 'evresp_') */ 508 509 extern double Pi; 510 extern double twoPi; 511 512 /* define a global flag to use if using "default" units */ 513 514 extern int def_units_flag; 515 516 /* define a pointer to a channel structure to use in determining the input and 517 output units if using "default" units and for use in error output*/ 518 519 extern struct channel *GblChanPtr; 520 extern float unitScaleFact; 521 522 /* define global variables for use in printing error messages */ 523 524 extern char *curr_file; 525 extern int curr_seq_no; 526 527 /* and set a global variable to contain the environment for the setjmp/longjmp 528 combination for error handling */ 529 530 extern jmp_buf jump_buffer; 531 532 double unwrap_phase(double phase, double prev_phase, double range, 533 double *added_value); 534 double wrap_phase(double phase, double range, double *added_value); 535 int use_delay(int flag); 536 #endif 537 538