From 586699fa813c176cd271793bcd375a6f722c56f6 Mon Sep 17 00:00:00 2001 From: Joaquim Date: Fri, 26 Jul 2024 11:49:36 +0100 Subject: [PATCH] Add new module under "seis" to read seismicity files stored with the ISF format (#8552) * Add new module under "seis" to read seismicity files stored with the ISF format * Missed calling GMT_End_IO() and that was why it was not working from Julia. * Fix bug in the date filter. * One more fix to the date filter. --- doc/rst/source/modules.rst | 4 + .../reference/supplemental-packages.rst | 1 + .../module_supplements_purpose.rst_ | 2 + doc/rst/source/supplements/seis/gmtisf.rst | 102 + src/seis/CMakeLists.txt | 2 +- src/seis/gmtisf.c | 442 ++++ src/seis/isf_head.h | 267 +++ src/seis/longopt/readisf_inc.h | 38 + src/seis/read_isf.c | 1834 +++++++++++++++++ 9 files changed, 2691 insertions(+), 1 deletion(-) create mode 100644 doc/rst/source/supplements/seis/gmtisf.rst create mode 100644 src/seis/gmtisf.c create mode 100644 src/seis/isf_head.h create mode 100644 src/seis/longopt/readisf_inc.h create mode 100644 src/seis/read_isf.c diff --git a/doc/rst/source/modules.rst b/doc/rst/source/modules.rst index 26667a537fd..69123d14b37 100644 --- a/doc/rst/source/modules.rst +++ b/doc/rst/source/modules.rst @@ -160,6 +160,7 @@ All modules are requested via a call to the :doc:`gmt` program. supplements/seis/meca supplements/seis/polar supplements/seis/sac + supplements/seis/gmtisf supplements/seis/grdshake supplements/seis/grdvs30 supplements/spotter/backtracker @@ -337,6 +338,7 @@ Supplemental Modules - :doc:`/supplements/segy/segyz` - :doc:`/supplements/seis/coupe` - :doc:`/supplements/seis/meca` + - :doc:`/supplements/seis/gmtisf` - :doc:`/supplements/seis/polar` - :doc:`/supplements/seis/sac` - :doc:`/supplements/spotter/backtracker` @@ -749,6 +751,8 @@ seis +-----------------------------------+--------------------+ | :doc:`/supplements/seis/sac` | |sac_purpose| | +-----------------------------------+--------------------+ +| :doc:`/supplements/seis/gmtisf` | |gmtisf_purpose| | ++-----------------------------------+--------------------+ | :doc:`/supplements/seis/grdshake` | |grdshake_purpose| | +-----------------------------------+--------------------+ | :doc:`/supplements/seis/grdvs30` | |grdvs30_purpose| | diff --git a/doc/rst/source/reference/supplemental-packages.rst b/doc/rst/source/reference/supplemental-packages.rst index 105655bc342..7fc926743a8 100644 --- a/doc/rst/source/reference/supplemental-packages.rst +++ b/doc/rst/source/reference/supplemental-packages.rst @@ -98,6 +98,7 @@ This package contains the programs :doc:`meca `, :doc:`polar `, :doc:`sac `, +:doc:`sac `, :doc:`grdvs30 `, and :doc:`grdshake `, which are used by seismologists for plotting focal mechanisms (including cross-sections diff --git a/doc/rst/source/supplements/module_supplements_purpose.rst_ b/doc/rst/source/supplements/module_supplements_purpose.rst_ index a8a2181346f..bf5333ea0bd 100644 --- a/doc/rst/source/supplements/module_supplements_purpose.rst_ +++ b/doc/rst/source/supplements/module_supplements_purpose.rst_ @@ -84,6 +84,8 @@ .. |pssac_purpose| replace:: Plot seismograms in SAC format +.. |gmtisf_purpose| replace:: Read seismicity data in the ISF formated file + .. |backtracker_purpose| replace:: Generate forward and backward flowlines and hotspot tracks .. |gmtpmodeler_purpose| replace:: Evaluate a plate motion model at given locations diff --git a/doc/rst/source/supplements/seis/gmtisf.rst b/doc/rst/source/supplements/seis/gmtisf.rst new file mode 100644 index 00000000000..4d2c342411b --- /dev/null +++ b/doc/rst/source/supplements/seis/gmtisf.rst @@ -0,0 +1,102 @@ +.. index:: ! gmtisf +.. include:: ../module_supplements_purpose.rst_ + +****** +gmtisf +****** + +|gmtisf_purpose| + +Synopsis +-------- + +.. include:: ../../common_SYN_OPTs.rst_ + +**gmt gmtisf** *ISFfile* +|SYN_OPT-R| +[ |-D|\ *date_start*\ [/*date_end*] ] +[ |-F|\ [**a**] ] +[ |-N| ] +[ |SYN_OPT--| ] + +Description +----------- + +Reads seismicity data from an ISC (https://www.isc.ac.uk/iscbulletin) formated *file.isc* and output [lon lat depth mag ...] +to standard output. + +Optional Arguments +------------------ + +.. _-R: + +.. |Add_-Rgeo| replace:: |Add_-R_auto_table| +.. include:: ../../explain_-Rgeo.rst_ + +.. _-D: + +**-D**\ *date_start*\ [/*date_end*] + Limit the output locations to data hose date is >= date1, or between date1 and date2. must be in ISO format, e.g, 2000-04-25" + +.. _-F: + +**-F**\ [*a*] + Select only events that have focal mechanisms. The default is output in Global CMT convention. Append 'a' for the AKI convention. + + Focal mechanisms in Global CMT convention. + + **1**,\ **2**: + longitude, latitude of event (**-:** option interchanges order) + + **3**: + depth of event in kilometers + + **4**,\ **5**,\ **6**: + strike, dip, and rake of plane 1 + + **7**,\ **8**,\ **9**: + strike, dip, and rake of plane 2 + + **10**,\ **11**: + mantissa and exponent of moment in dyne-cm + + + Focal mechanisms in Aki and Richards convention. + + **1**,\ **2**: + longitude, latitude of event (**-:** option interchanges order) + + **3**: + depth of event in kilometers + + **4**,\ **5**,\ **6**: + strike, dip and rake in degrees + + **7**: + magnitude + +.. _-N: + +**-N** + The default is to output time information [year month day hour minute] as the last 5 columns. + Use this option to skip those last 5 columns. + +Examples +-------- + +Extract the `lon lat depth mag` seismicity from file `file.isf` obtained at (https://www.isc.ac.uk/iscbulletin/) +and limiting over a geographic region and a date interval:: + + gmt isf file.isf -R-15/0/30/45 -D2001-01-01/2005-12-31 -N > seismicity.dat + +The above command can be piped directly to `psxy` to create a seismicity map. _e.,g._:: + + gmt isf file.isf -R-15/0/30/45 -D2001-01-01/2005-12-31 -N | gmt psxy -R -JM14c -Sc0.01 -Gblack -Ba -P > seis.ps + +See Also +-------- + +:doc:`psmeca`, +:doc:`pspolar`, +:doc:`pscoupe`, +:doc:`gmt `, :doc:`psbasemap `, :doc:`psxy ` diff --git a/src/seis/CMakeLists.txt b/src/seis/CMakeLists.txt index e3f9f155946..e662fdc9f5b 100644 --- a/src/seis/CMakeLists.txt +++ b/src/seis/CMakeLists.txt @@ -27,6 +27,6 @@ set (SUPPL_NAME seis) set (SUPPL_HEADERS meca.h meca_symbol.h utilmeca.h seis_defaults.h sacio.h) AUX_SOURCE_DIRECTORY (longopt SUPPL_LONG_OPT_H) -set (SUPPL_PROGS_SRCS psmeca.c pspolar.c pscoupe.c pssac.c grdshake.c grdvs30.c ${SUPPL_LONG_OPT_H}) +set (SUPPL_PROGS_SRCS psmeca.c pspolar.c pscoupe.c pssac.c grdshake.c grdvs30.c gmtisf.c ${SUPPL_LONG_OPT_H}) set (SUPPL_LIB_SRCS ${SUPPL_PROGS_SRCS} utilmeca.c sacio.c) set (SUPPL_EXAMPLE_FILES README.seis) diff --git a/src/seis/gmtisf.c b/src/seis/gmtisf.c new file mode 100644 index 00000000000..9369fd08c94 --- /dev/null +++ b/src/seis/gmtisf.c @@ -0,0 +1,442 @@ +/*-------------------------------------------------------------------- + * + * Copyright (c) 1991-2024 by the GMT Team (https://www.generic-mapping-tools.org/team.html) + * See LICENSE.TXT file for copying and redistribution conditions. + * + * This program is free software; you can redistribute it and/or modify + * it under the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation; version 3 or any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU Lesser General Public License for more details. + * + * Contact info: www.generic-mapping-tools.org + *--------------------------------------------------------------------*/ +/* + * Brief synopsis: + * + * Joaquim Luis + * Date: 24-JUL-2024 + * Version: 6 API + */ + +#include "gmt_dev.h" +#include "longopt/readisf_inc.h" + +#define THIS_MODULE_CLASSIC_NAME "gmtisf" +#define THIS_MODULE_MODERN_NAME "gmtisf" +#define THIS_MODULE_LIB "seis" +#define THIS_MODULE_PURPOSE "Read seismicity data in the ISF formated file" +#define THIS_MODULE_KEYS ">D}" +#define THIS_MODULE_NEEDS "" +#define THIS_MODULE_OPTIONS "->RV" + +#include "read_isf.c" + +/* Control structure */ + +struct READISF_CTRL { + struct READISF_In { /* Input files */ + bool active; + char *file; + } In; + struct READISF_D { /* For date control */ + bool active; + bool two_dates; + struct GMT_GCAL date1; + struct GMT_GCAL date2; + } D; + struct READISF_F { /* Fault mechanism */ + bool active; + bool aki; + } F; + struct READISF_N { /* No times */ + bool active; + } N; +}; + + +static void *New_Ctrl(struct GMT_CTRL *GMT) { /* Allocate and initialize a new control structure */ + struct READISF_CTRL *C; + + C = gmt_M_memory(GMT, NULL, 1, struct READISF_CTRL); + + return (C); +} + +static void Free_Ctrl(struct GMT_CTRL *GMT, struct READISF_CTRL *C) { /* Deallocate control structure */ + if (!C) return; + gmt_M_free(GMT, C->In.file); + gmt_M_free(GMT, C); +} + +static int usage(struct GMTAPI_CTRL *API, int level) { + /* This displays the pssac synopsis and optionally full usage information */ + + const char *name = gmt_show_name_and_purpose (API, THIS_MODULE_LIB, THIS_MODULE_CLASSIC_NAME, THIS_MODULE_PURPOSE); + if (level == GMT_MODULE_PURPOSE) return (GMT_NOERROR); + GMT_Usage(API, 0, "usage: %s isffile [-Ddate1[/date2]] [%s] [-F[a]] [-N] [%s]\n", name, GMT_Rgeoz_OPT, GMT_PAR_OPT); + + if (level == GMT_SYNOPSIS) return EXIT_FAILURE; + + GMT_Message(API, GMT_TIME_NONE, " OPTIONAL ARGUMENTS:\n"); + GMT_Option(API, "R"); + GMT_Usage(API, 1, "\n-Ddate1[/date2]"); + GMT_Usage(API, -2, "Limit the output to data >= date1, or between date1 and date2. must be in ISO format, e.g, 2000-04-25"); + GMT_Usage(API, 1, "\n-F[a]"); + GMT_Usage(API, -2, "Select only events that have focal mechanisms. The default is Global CMT convention. Append 'a' for the AKI convention."); + GMT_Usage(API, 1, "\n-N"); + GMT_Usage(API, -2, "The default is to output time information [year month day hour minute] as the last 5 columns. Use this option to skip those last 5 columns."); + + return EXIT_FAILURE; +} + +static int parse(struct GMT_CTRL *GMT, struct READISF_CTRL *Ctrl, struct GMT_OPTION *options) { + /* This parses the options provided to pssac and sets parameters in Ctrl. + * Note Ctrl has already been initialized and non-zero default values set. + * Any GMT common options will override values set previously by other commands. + * It also replaces any file names specified as input or output with the data ID + * returned when registering these sources/destinations with the API. + */ + + unsigned int n_errors = 0, n_files = 0, pos = 0; + int i, j, k; + size_t n_alloc = 0, len; + char txt_a[GMT_LEN256] = {""}, txt_b[GMT_LEN256] = {""}, p[GMT_BUFSIZ] = {""}; + char *pch = NULL; + double t; + struct GMT_OPTION *opt = NULL; + struct GMTAPI_CTRL *API = GMT->parent; + + for (opt = options; opt; opt = opt->next) { /* Process all the options given */ + switch (opt->option) { + + case '<': /* Input file (only one is accepted) */ + if (n_files++ > 0) break; + if ((Ctrl->In.active = gmt_check_filearg (GMT, '<', opt->arg, GMT_IN, GMT_IS_GRID)) != 0) + Ctrl->In.file = strdup (opt->arg); + else + n_errors++; + break; + case 'D': + n_errors += gmt_M_repeated_module_option(API, Ctrl->D.active); + if ((pch = strchr(opt->arg, '/')) != NULL) { /* Have two dates */ + gmt_scanf_arg(GMT, (&pch[1]), GMT_IS_ABSTIME, false, &t); + gmt_gcal_from_dt(GMT, t, &Ctrl->D.date2); /* Convert t to a complete calendar structure */ + pch[0] = '\0'; /* Hide the second date */ + Ctrl->D.two_dates = true; + } + gmt_scanf_arg(GMT, &opt->arg[0], GMT_IS_ABSTIME, false, &t); + gmt_gcal_from_dt(GMT, t, &Ctrl->D.date1); /* Convert t to a complete calendar structure */ + break; + case 'F': + n_errors += gmt_M_repeated_module_option(API, Ctrl->F.active); + if (opt->arg[0] == 'a') Ctrl->F.aki = true; + break; + case 'N': + n_errors += gmt_M_repeated_module_option(API, Ctrl->N.active); + break; + + /* Processes program-specific parameters */ + } + } + + /* Check that the options selected are mutually consistent */ + n_errors += gmt_M_check_condition (GMT, n_files != 1, "Syntax error: Must specify a single ISF file\n"); + + return (n_errors ? GMT_PARSE_ERROR : GMT_OK); +} + +static bool filter_date(struct READISF_CTRL *Ctrl, int yyyy, int mm, int dd) { + /* Returns true if the date is out of date range */ + bool is_out = false; + if (yyyy < Ctrl->D.date1.year) is_out = true; + else if (mm < Ctrl->D.date1.month) is_out = true; + else if (dd < Ctrl->D.date1.day_m) is_out = true; + if (!is_out && Ctrl->D.two_dates) { + if (yyyy > Ctrl->D.date2.year) is_out = true; + else if (mm > Ctrl->D.date2.month) is_out = true; + else if (dd > Ctrl->D.date2.day_m) is_out = true; + } + return is_out; +} + +#define bailout(code) {gmt_M_free_options(mode); return (code);} +#define Return(code) {Free_Ctrl(GMT, Ctrl); gmt_end_module(GMT, GMT_cpy); bailout (code);} + +EXTERN_MSC int GMT_gmtisf(void *V_API, int mode, void *args) { /* High-level function that implements the pssac task */ + + FILE *fp; + char timfix,epifix,depfix,antype,loctype,magind; + char *etype, *author, *origid, *magtype, line[ISF_LINE_LEN]; + char **mag_t, evid[ISF_LINE_LEN], region[ISF_LINE_LEN]; + char f_type[6], f_plane[6]; + bool got_event = false, event_end, tensor_end, got_region = false; + int error = GMT_NOERROR; + int i, in, mag_c = 0, event_c, idx_min_rms, np, ns, n_out_cols; + int export_aki = false, export_cmt = false, export_tensor = false; + int yyyy,mm,dd,hh,mi,ss,msec,strike,ndef,nsta,gap; + int *years, *months, *days, *hours, *minutes; + float stime,sdobs,lat,lon,depth,smaj,smin,sdepth,mindist,maxdist; + float mag, magerr; + float *rms, *lats, *lons, *depths, *mags; + float scalar_moment, fclvd, mrr, mtt, mpp, mrt, mtp, mpr; + float scalar_moment_unc, fclvd_unc, mrr_unc, mtt_unc, mpp_unc, mrt_unc, mtp_unc, mpr_unc, duration; + float strike1, dip1, rake1, strike2, dip2, rake2; + float t_val, t_azim, t_pl, b_val, b_azim, b_pl, p_val, p_azim, p_pl; + double west = 0.0, east = 0.0, south = 0.0, north = 0.0; + + /* Moment tensor variables */ + int scale_factor, nsta1, nsta2, ncomp1, ncomp2, got_momten_line1, got_momten_line2, centroid; + + double out[16]; + struct READISF_CTRL *Ctrl = NULL; + struct GMT_CTRL *GMT = NULL, *GMT_cpy = NULL; /* General GMT internal parameters */ + struct GMT_OPTION *options = NULL; + struct GMTAPI_CTRL *API = gmt_get_api_ptr (V_API); /* Cast from void to GMTAPI_CTRL pointer */ + struct GMT_RECORD *Out = NULL; + + /*----------------------- Standard module initialization and parsing ----------------------*/ + if (API == NULL) return GMT_NOT_A_SESSION; + if (mode == GMT_MODULE_PURPOSE) return usage(API, GMT_MODULE_PURPOSE); /* Return the purpose of program */ + options = GMT_Create_Options(API, mode, args); if (API->error) return (API->error); /* Set or get option list */ + + if ((error = gmt_report_usage(API, options, 0, usage)) != GMT_NOERROR) bailout(error); /* Give usage if requested */ + + /* Parse the command-line arguments; return if errors are encountered */ + + if ((GMT = gmt_init_module(API, THIS_MODULE_LIB, THIS_MODULE_CLASSIC_NAME, THIS_MODULE_KEYS, THIS_MODULE_NEEDS, module_kw, &options, &GMT_cpy)) == NULL) bailout(API->error); /* Save current state */ + if (GMT_Parse_Common(API, THIS_MODULE_OPTIONS, options)) Return(API->error); + + Ctrl = New_Ctrl(GMT); /* Allocate and initialize a new control structure */ + if ((error = parse(GMT, Ctrl, options)) != 0) Return(error); + + /*---------------------------- This is the pssac main code ----------------------------*/ + + /* I Don't remember where size 100 came from. ISC site example in fortran also uses 100 */ + etype = gmt_M_memory(GMT, NULL, ISF_ETYPE_LEN, char *); + author = gmt_M_memory(GMT, NULL, ISF_AUTHOR_LEN, char *); + origid = gmt_M_memory(GMT, NULL, ISF_ORIGID_LEN, char *); + magtype = gmt_M_memory(GMT, NULL, ISF_MAGTYPE_LEN, char *); + mags = gmt_M_memory(GMT, NULL, 100, float *); + rms = gmt_M_memory(GMT, NULL, 100, float *); + lats = gmt_M_memory(GMT, NULL, 100, float *); + lons = gmt_M_memory(GMT, NULL, 100, float *); + depths = gmt_M_memory(GMT, NULL, 100, float *); + years = gmt_M_memory(GMT, NULL, 100, int *); + months = gmt_M_memory(GMT, NULL, 100, int *); + days = gmt_M_memory(GMT, NULL, 100, int *); + hours = gmt_M_memory(GMT, NULL, 100, int *); + minutes = gmt_M_memory(GMT, NULL, 100, int *); + mag_t = gmt_M_memory(GMT, NULL, 100, char *); + for (i = 0; i < 100; i++) + mag_t[i] = gmt_M_memory(GMT, NULL, GMT_LEN16, char); + + n_out_cols = Ctrl->F.active ? 11+5 : 4+5; + if (Ctrl->F.aki) n_out_cols = 7+5; /* The Aki convention uses 7 columns */ + if (Ctrl->N.active) n_out_cols -= 5; + + gmt_set_geographic(GMT, GMT_OUT); /* Output lon/lat */ + + if (GMT_Init_IO(API, GMT_IS_DATASET, GMT_IS_PLP, GMT_OUT, GMT_ADD_DEFAULT, 0, options) != GMT_NOERROR) /* Establishes data output */ + Return(API->error); + + if (GMT_Begin_IO(API, GMT_IS_DATASET, GMT_OUT, GMT_HEADER_OFF) != GMT_NOERROR) /* Enables data output and sets access mode */ + Return(API->error); + + if ((error = GMT_Set_Columns(API, GMT_OUT, n_out_cols, GMT_COL_FIX_NO_TEXT)) != GMT_NOERROR) Return (error); + for (i = 2; i < n_out_cols; i++) GMT->current.io.col_type[GMT_OUT][i] = GMT_IS_FLOAT; + + if ((fp = fopen(Ctrl->In.file, "rt")) == NULL) { + GMT_Report(GMT->parent, GMT_MSG_ERROR, "Unable to open file %s [permission trouble?]\n", Ctrl->In.file); + return GMT_NOTSET; + } + + if (GMT->common.R.active[0]) { + got_region = true; + west = GMT->common.R.wesn[XLO]; + east = GMT->common.R.wesn[XHI]; + south = GMT->common.R.wesn[YLO]; + north = GMT->common.R.wesn[YHI]; + } + + Out = gmt_new_record(GMT, out, NULL); + + if (Ctrl->F.active) { /* Focal mechanims */ + event_end = tensor_end = false; + while (fgets (line, ISF_LINE_LEN, fp) != NULL) { + gmt_chop(line); /* Chops off any CR and/or LF */ + if (!read_event_id(line, evid, region)) continue; + if(!read_origin_head(line)) { + event_c = read_event_data(fp,line,&yyyy,&mm,&dd,&hh,&mi,&ss,&msec,&timfix,&stime,&sdobs, + &lat,&lon,&epifix,&smaj,&smin,&strike,&depth,&depfix,&sdepth,&ndef, + &nsta,&gap,&mindist,&maxdist,&antype,&loctype,etype,author,origid, + lats,lons,rms,depths,years,months,days,hours,minutes,&idx_min_rms); + if (event_c > 0) got_event = true; + if (!read_momten_head_1(line)) goto L1; /* Awfull, I know */ + } + else if (!read_origin_centroid(line)) { /* Harvard Moment Tensor event */ + i = 0; + got_momten_line1 = got_momten_line2 = false; + if (fgets(line, ISF_LINE_LEN, fp) != NULL) + if (!read_momten_head_1(line)) i++; + if (fgets(line, ISF_LINE_LEN, fp) != NULL) + if (!read_momten_head_2(line)) i++; + if (fgets(line, ISF_LINE_LEN, fp) != NULL) { + if (!read_momten_line_1(line, &scale_factor, &scalar_moment, &fclvd, &mrr, &mtt, + &mpp, &mrt, &mtp, &mpr, &nsta1, &nsta2, author)) + got_momten_line1 = true; + } + if (fgets(line, ISF_LINE_LEN, fp) != NULL) { + if (!read_momten_line_2(line,&scalar_moment_unc,&fclvd_unc,&mrr_unc,&mtt_unc, + &mpp_unc,&mrt_unc,&mtp_unc,&mpr_unc,&ncomp1,&ncomp2,&duration)) + got_momten_line2 = true; + } + centroid = true; + } + else if (!read_momten_head_1(line)) { + got_momten_line1 = false; +L1: + if (fgets(line, ISF_LINE_LEN, fp) != NULL) + read_momten_head_2(line); + if (fgets(line, ISF_LINE_LEN, fp) != NULL) { + if (!read_momten_line_1(line, &scale_factor, &scalar_moment, &fclvd, &mrr, &mtt, + &mpp, &mrt, &mtp, &mpr, &nsta1, &nsta2, author)) + got_momten_line1 = true; + } + centroid = false; + } + else if (!read_fault_plane_head(line)) { + if (fgets(line, ISF_LINE_LEN, fp) != NULL) + read_fault_plane (line, f_type, &strike1, &dip1, &rake1, &np, &ns, f_plane, author); + if (fgets(line, ISF_LINE_LEN, fp) != NULL) + read_fault_plane (line, f_type, &strike2, &dip2, &rake2, &np, &ns, f_plane, author); + } + else if (!read_axes_head(line)) { + if (fgets(line, ISF_LINE_LEN, fp) != NULL) + if (!read_axes(line, &scale_factor, &t_val, &t_azim, &t_pl, &b_val, &b_azim, + &b_pl, &p_val, &p_azim, &p_pl,author)); + tensor_end = true; + } + else if (tensor_end && !read_netmag_head(line)) { + mag_c = read_mags(fp,line,magtype,&magind,&mag,&magerr,&nsta,author,origid,mag_t,mags); + event_end = true; + } + + if (got_event && event_end && tensor_end) { + /* The reported values respect the last entry in the event */ + lon = lons[idx_min_rms]; + lat = lats[idx_min_rms]; + if (got_region && (lon < west || lon > east || lat < south || lat > north)) { + got_event = event_end = false; mag_c = 0; + continue; + } + if (Ctrl->F.aki) { + if (mag_c >= 1) /* Have multiple magnitudes. Choose one */ + mag = select_mag(mag_c, mags, mag_t); + else + mag = 0; + } + + depth = depths[idx_min_rms]; + if (depth == ISF_NULL) depth = 0; + yyyy = years[idx_min_rms]; + mm = months[idx_min_rms]; + dd = days[idx_min_rms]; + got_event = event_end = tensor_end = false; + mag_c = 0; + + /* See if user set date bounds */ + if (Ctrl->D.active) { + if (filter_date(Ctrl, yyyy, mm, dd)) { + got_event = event_end = false; mag_c = 0; + continue; + } + } + + out[GMT_X] = lon; out[GMT_Y] = lat; out[2] = depth; + if (Ctrl->F.aki) { + out[3] = strike1; out[4] = dip1; out[5] = rake1; out[6] = mag; + out[7] = yyyy; out[8] = mm; out[9] = dd; out[10] = hh; out[11] = mi; + } + else { + out[3] = strike1; out[4] = dip1; out[5] = rake1; + out[6] = strike2; out[7] = dip2; out[8] = rake2; out[9] = scalar_moment; out[10] = scale_factor; + out[11] = yyyy; out[12] = mm; out[13] = dd; out[14] = hh; out[15] = mi; + } + GMT_Put_Record(GMT->parent, GMT_WRITE_DATA, Out); + } + } + } + else { /* Just the event data (no focal mechanism) */ + event_end = true; + while (fgets (line, ISF_LINE_LEN, fp) != NULL) { + gmt_chop(line); /* Chops off any CR and/or LF */ + if (!read_event_id(line, evid, region)) continue; + if (!read_origin_head(line)) { + event_c = read_event_data(fp,line,&yyyy,&mm,&dd,&hh,&mi,&ss,&msec,&timfix,&stime,&sdobs, + &lat,&lon,&epifix,&smaj,&smin,&strike,&depth,&depfix,&sdepth,&ndef, + &nsta,&gap,&mindist,&maxdist,&antype,&loctype,etype,author,origid, + lats,lons,rms,depths,years,months,days,hours,minutes,&idx_min_rms); + if (event_c > 0) got_event = true; + } + else if (!read_netmag_head(line)) { + mag_c = read_mags(fp,line,magtype,&magind,&mag,&magerr,&nsta,author,origid,mag_t,mags); + event_end = true; + } + if (got_event && event_end) { + /* Select the event detection that has the minimum RMS */ + lon = lons[idx_min_rms]; + lat = lats[idx_min_rms]; + if (got_region && (lon < west || lon > east || lat < south || lat > north)) { + got_event = event_end = false; mag_c = 0; + continue; + } + depth = depths[idx_min_rms]; + if (depth == ISF_NULL) depth = 0; + yyyy = years[idx_min_rms]; + mm = months[idx_min_rms]; + dd = days[idx_min_rms]; + hh = hours[idx_min_rms]; + mi = minutes[idx_min_rms]; + + /* See if user set date bounds */ + if (Ctrl->D.active) { + if (filter_date(Ctrl, yyyy, mm, dd)) { + got_event = event_end = false; mag_c = 0; + continue; + } + } + + if (mag_c >= 1) /* Have multiple magnitudes. Choose one */ + mag = select_mag(mag_c, mags, mag_t); + else + mag = 0; + if (depth == ISF_NULL) depth = 0; + + got_event = event_end = false; + mag_c = 0; + out[GMT_X] = lon; out[GMT_Y] = lat; out[2] = depth; out[3] = mag; + out[4] = yyyy; out[5] = mm; out[6] = dd; out[7] = hh; out[8] = mi; + GMT_Put_Record(GMT->parent, GMT_WRITE_DATA, Out); + } + } + } + + fclose(fp); + gmt_M_free(GMT, etype); gmt_M_free(GMT, author); gmt_M_free(GMT, origid); gmt_M_free(GMT, magtype); + gmt_M_free(GMT, mags); gmt_M_free(GMT, rms); gmt_M_free(GMT, lats); gmt_M_free(GMT, lons); + gmt_M_free(GMT, depths); gmt_M_free(GMT, years); gmt_M_free(GMT, months); gmt_M_free(GMT, days); + gmt_M_free(GMT, hours); gmt_M_free(GMT, minutes); + for (i = 0; i < 100; i++) gmt_M_free(GMT, mag_t[i]); + gmt_M_free(GMT, mag_t); + gmt_M_free(GMT, Out); + + if (GMT_End_IO(API, GMT_OUT, 0) != GMT_NOERROR) { /* Disables further data output */ + Return(API->error); + } + + Return(GMT_OK); +} diff --git a/src/seis/isf_head.h b/src/seis/isf_head.h new file mode 100644 index 00000000000..70040d765d5 --- /dev/null +++ b/src/seis/isf_head.h @@ -0,0 +1,267 @@ +/*-------------------------------------------------------------------- + * Copyright (c) 2004-2012 by J. Luis + * + * This program is part of Mirone and is free software; you can redistribute + * it and/or modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * Contact info: w3.ualg.pt/~jluis/mirone + *--------------------------------------------------------------------*/ + +#define ISF_NULL 9999999 +#define ISF_LINE_LEN 140 +#define ISF_COMM_LEN 80 +#define ISF_EVID_LEN 8 +#define ISF_REGION_LEN 65 +#define ISF_ETYPE_LEN 2 +#define ISF_AUTHOR_LEN 9 +#define ISF_ORIGID_LEN 8 +#define ISF_MAGTYPE_LEN 5 +#define ISF_STA_LEN 5 +#define ISF_NET_LEN 9 +#define ISF_CHAN_LEN 3 +#define ISF_PHASE_LEN 8 +#define ISF_ARRID_LEN 8 +#define ISF_F_TYPE_LEN 3 +#define ISF_F_PLANE_LEN 5 +#define ISF_I_LOCTYPE_LEN 6 +#define ISF_COUNTRY_LEN 3 +#define ISF_POSTCODE_LEN 10 +#define ISF_I_SCALE_LEN 5 +#define ISF_AUXID_LEN 4 +#define ISF_GROUPID_LEN 8 + +#define ISF_NUM_STA 200 +#define ISF_NUM_PARAM 100 + +char isf_error[ISF_LINE_LEN*2]; +char isf_prev_line_type[ISF_LINE_LEN]; + +/* Utilities used by other functions */ + +int check_prev_line_type(char *line_type); +void print_float(FILE *fp,float x,int width,int max_prec); +int is_null(int i); +int check_int(char *s); +int check_float(char *s); +int partline(char *part,char *line,int offset,int numchars); +int check_whole(char *s); +int all_blank(char *s); + + +/* Functions to set variables to impossible values */ + +int nullify_data_type(char *data_type, char *subtype, char *data_format, char *subformat); + +int nullify_event_id(char *evid, char *region); + +int nullify_origin(int *yyyy, int *mm, int *dd, int *hh, int *mi, + int *ss, int *msec, char *timfix, float *stime, float *sdobs, + float *lat, float *lon, char *epifix, float *smaj, float *smin, + int *strike, float *depth, char *depfix, float *sdepth, + int *ndef, int *nsta, int *gap, float *mindist, + float *maxdist, char *antype, char *loctype, char *etype, + char *author, char *origid); + +int nullify_origin_param(char **param, char **value, int *numparam); + +int nullify_momten(int *scale_factor, float *scalar_moment, float *fclvd, + float *mrr, float *mtt, float *mpp, float *mrt, + float *mtp, float *mpr,int *nsta1, int *nsta2, char *author, + float *scalar_moment_unc, float *fclvd_unc, + float *mrr_unc, float *mtt_unc, float *mpp_unc, + float *mrt_unc, float *mtp_unc, float *mpr_unc, + int *ncomp1, int *ncomp2, float *duration); + +int nullify_fault_plane (char *f_type, float *strike, float *dip, float *rake, + int *np, int *ns, char *f_plane, char *author); + +int nullify_axes(int *scale_factor, float *t_val, float *t_azim, + float *t_pl, float *b_val, float *b_azim, float *b_pl, + float *p_val, float *p_azim, float *p_pl, char *author); + +int nullify_axes_err(float *t_val_unc, float *t_azim_unc, float *t_pl_unc, + float *b_val_unc, float *b_azim_unc, float *b_pl_unc, + float *p_val_unc, float *p_azim_unc, float *p_pl_unc, + float *fclvd); + +int nullify_effects(char *heard, char *felt, char *damage, + char *casualties, char *uplift, char *subsidence, + char *fault, char *tsunami, char *seiche, char *volcano, + char *acoustic, char *gravity, char *t_wave, + char *liquification, char *geyser, char *landslide, + char *sandblow, char *cracks, char *lights, char *odours, + char *loctype, float *lat, float *lon, float *dist, + float *azim, char *country, char *postcode, char *net, + char *sta, float *intensity1, char *modifier, + float *intensity2, char *scale, char* author); + +int nullify_netmag(char *magtype, char* magind, float* mag, float* magerr, int* nsta, char *author, char* origid); +int nullify_netmag_sta(char **sta, int *n); +int nullify_netmag_basis(char *param, char *value); +int nullify_phase(char *sta, float *dist, float *esaz, char *phase, + int *hh, int *mi, int *ss, + int *msec, float *timeres, float *azim, float *azimres, + float *slow, float *slowres, char *timedef, char *azimdef, + char *slowdef, float *snr, float *amp, float *per, + char *picktype, char *sp_fm, char *detchar, char *magtype, + char *magind, float *mag, char *arrid); + +int nullify_phase_measure(char **param, char **value, int *numparam); +int nullify_phase_origid(char *origid); +int nullify_phase_info(char *net, char *chan, char *filter, float *filter_min, + float *filter_max, char *phase, int *yyyy, int *mm, + int *dd, float *time_unc, float *time_weight, + float *azim_unc, float *azim_weight, float *slow_unc, + float *slow_weight, float *amp_unc, float *per_unc, + float *mag_unc, char *author, char *arrid); + +int nullify_phase_min(float *timeoffset, float *azoffset, float *slowoffset, float *ampoffset, float *peroffset, + float *magoffset); + +int nullify_phase_max(float *timeoffset, float *azoffset, float *slowoffset, float *ampoffset, float *peroffset, + float *magoffset); + +int nullify_phase_correc(float *timecorr, float *azcorr,float *slowcorr, float *ampcorr, float *percorr,float *magcorr); + +int nullify_phase_original(char *chan, char *sta, int *yyyy, int *mm, + int *dd, int *hh, int *mi, int *ss, int *msec, + float *azim, float *slow, float *amp, float *per, + float *mag); + +int nullify_comment(char *comment); + +/* Functions for reading ISF format files. */ + +int read_data_type(char *line, char *data_type, char *subtype, char *data_format,char *subformat); +int read_event_id(char *line, char *evid, char *region); +int read_origin_head(char *line); +int read_origin(char *line, int *yyyy, int *mm, int *dd, int *hh, int *mi, + int *ss, int *msec, char *timfix, float *stime, float *sdobs, + float *lat, float *lon, char *epifix, float *smaj, float *smin, + int *strike, float *depth, char *depfix, float *sdepth, + int *ndef, int *nsta, int *gap, float *mindist, + float *maxdist, char *antype, char *loctype, char *etype, + char *author, char *origid); + +int read_origin_prime(char *line); +int read_origin_centroid(char *line); +int read_origin_param(char *line, char **param, char **value, char **error, int *n); +int read_momten_head_1(char *line); +int read_momten_head_2(char *line); + +int read_momten_line_1(char *line, int *scale_factor, float *scalar_moment, + float *fclvd, float *mrr, float *mtt, float *mpp, + float *mrt, float *mtp, float *mpr, int *nsta1, + int *nsta2, char *author); + + +int read_momten_line_2(char *line, float *scalar_moment_unc, float *fclvd_unc, + float *mrr_unc, float *mtt_unc, float *mpp_unc, + float *mrt_unc, float *mtp_unc, float *mpr_unc, + int *ncomp1, int *ncomp2, float *duration); + +int read_fault_plane_head (char *line); + +int read_fault_plane (char *line, char *f_type, float *strike, float *dip, + float *rake, int *np, int *ns, char *f_plane, + char *author); + +int read_axes_head(char *line); + +int read_axes_err_head(char *line); + +int read_axes(char *line, int *scale_factor, float *t_val, float *t_azim, + float *t_pl, float *b_val, float *b_azim, float *b_pl, + float *p_val, float *p_azim, float *p_pl, char *author); + +int read_axes_err(char *line, float *t_val_unc, float *t_azim_unc, + float *t_pl_unc, float *b_val_unc, float *b_azim_unc, + float *b_pl_unc, float *p_val_unc, float *p_azim_unc, + float *p_pl_unc, float *fclvd); + +int read_netmag_head(char *line); + +int read_netmag(char *line, char *magtype, char* magind, float* mag, + float* magerr, int* nsta, char *author, char* origid); + +int read_netmag_sta(char *line, char **sta, int* n); + +int read_netmag_basis(char *line, char *param, char *value); + +int read_effects_head(char *line); + +int read_effects(char *line, char *heard, char *felt, char *damage, + char *casualties, char *uplift, char *subsidence, + char *fault, char *tsunami, char *seiche, char *volcano, + char *acoustic, char *gravity, char *t_wave, + char *liquification, char *geyser, char *landslide, + char *sandblow, char *cracks, char *lights, char *odours, + char *loctype, float *lat, float *lon, float *dist, + float *azim, char *country, char *postcode, char *net, + char *sta, float *intensity1, char *modifier, + float *intensity2, char *scale, char* author); + +int read_phase_head(char *line); + +int read_phase_origid(char *line, char *origid); + +int read_phase(char *line, char *sta, float *dist, float *esaz, + char *phase,int *hh, + int *mi, int *ss, int *msec, float *timeres, float *azim, + float *azimres, float *slow, float *slowres, char *timedef, + char *azimdef, char *slowdef, float *snr, float *amp, + float *per, char *picktype, char *sp_fm, char *detchar, + char *magtype, char *magind, float *mag, char *arrid); + +int read_phase_info_head(char *line); + +int read_phase_info(char *line, char *net, char *chan, char *filter, + float *filter_min, + float *filter_max, char *phase, int *yyyy, int *mm, + int *dd, float *time_unc, float *time_weight, + float *azim_unc, float *azim_weight, float *slow_unc, + float *slow_weight, float *amp_unc, float *per_unc, + float *mag_unc, char *author, char *arrid); + +int read_phase_measure(char *line, char **param, char **value, char **error, int *n); + +int read_phase_min(char *line, float *timeoffset, float *azoffset, + float *slowoffset, float *ampoffset, float *peroffset, + float *magoffset); + +int read_phase_max(char *line, float *timeoffset, float *azoffset, + float *slowoffset, float *ampoffset, float *peroffset, + float *magoffset); + +int read_phase_correc(char *line, float *timecorr, float *azcorr, + float *slowcorr, float *ampcorr, float *percorr, + float *magcorr ); + +int read_phase_original(char *line, char *chan, char *sta, int *yyyy, int *mm, + int *dd, int *hh, int *mi, int *ss, int *msec, + float *azim, float *slow, float *amp, float *per, + float *mag); + + +int read_comment(char *line,char *comment); + +int read_stop(char *line); + +int read_event_data(FILE *fp, char *line, int *yyyy, int *mm, int *dd, int *hh, int *mi, + int *ss, int *msec, char *timfix, float *stime, float *sdobs, + float *lat, float *lon, char *epifix, float *smaj, float *smin, + int *strike, float *depth, char *depfix, float *sdepth, + int *ndef, int *nsta, int *gap, float *mindist, float *maxdist, + char *antype, char *loctype, char *etype, char *author, char *origid, + float *lats, float *lons, float *rms, float *depths, int *years, + int *months, int *days, int *hours, int *minutes, int *idx_min_rms); + +int read_mags(FILE *fp, char *line, char *magtype, char* magind, float* mag, float* magerr, + int* nsta, char* author, char* origid, char **mag_t, float *mags); diff --git a/src/seis/longopt/readisf_inc.h b/src/seis/longopt/readisf_inc.h new file mode 100644 index 00000000000..27c7bc0d569 --- /dev/null +++ b/src/seis/longopt/readisf_inc.h @@ -0,0 +1,38 @@ +/*-------------------------------------------------------------------- + * + * Copyright (c) 1991-2024 by the GMT Team (https://www.generic-mapping-tools.org/team.html) + * See LICENSE.TXT file for copying and redistribution conditions. + * + * This program is free software; you can redistribute it and/or modify + * it under the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation; version 3 or any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU Lesser General Public License for more details. + * + * Contact info: www.generic-mapping-tools.org + *--------------------------------------------------------------------*/ + +#ifndef READISF_INC_H +#define READISF_INC_H + +/* Translation table from long to short module options, directives and modifiers */ + +static struct GMT_KEYWORD_DICTIONARY module_kw[] = { /* Local options for this module */ + /* separator, short_option, long_option, + short_directives, long_directives, + short_modifiers, long_modifiers, + transproc_mask */ + { 0, 'D', "date", "", "", "", "", + GMT_TP_STANDARD }, + { 0, 'F', "focal", "a", "akima", "", "", + GMT_TP_STANDARD }, + { 0, 'N', "notime", "", "", "", "", + GMT_TP_STANDARD }, + + { 0, '\0', "", "", "", "", "", 0} /* End of list marked with empty option and strings */ +}; + +#endif /* !READISF_INC_H */ diff --git a/src/seis/read_isf.c b/src/seis/read_isf.c new file mode 100644 index 00000000000..e718afa6189 --- /dev/null +++ b/src/seis/read_isf.c @@ -0,0 +1,1834 @@ +/* From https://www.isc.ac.uk/standards/isf/download/isf_c.tar */ + +#include "isf_head.h" + +float select_mag(int mag_c, float *mags, char **mag_t); + +int read_event_data(FILE *fp, char *line, int *yyyy, int *mm, int *dd, int *hh, int *mi, + int *ss, int *msec, char *timfix, float *stime, float *sdobs, + float *lat, float *lon, char *epifix, float *smaj, float *smin, + int *strike, float *depth, char *depfix, float *sdepth, + int *ndef, int *nsta, int *gap, float *mindist, float *maxdist, + char *antype, char *loctype, char *etype, char *author, char *origid, + float *lats, float *lons, float *rms, float *depths, int *years, + int *months, int *days, int *hours, int *minutes, int *idx_min_rms) { + + int done = false, event_c = 0; + float rms_min = 1e9, gap_min = 360; + *idx_min_rms = 0; /* I'm currently using GAP instead of RMS */ + + while (!done && (fgets (line, ISF_LINE_LEN, fp) != NULL)) { + if (!read_origin(line,yyyy,mm,dd,hh,mi,ss,msec,timfix,stime,sdobs, + lat,lon,epifix,smaj,smin,strike,depth,depfix,sdepth,ndef, + nsta,gap,mindist,maxdist,antype,loctype,etype,author,origid)) { + lats[event_c] = *lat; lons[event_c] = *lon; + rms[event_c] = *sdobs; depths[event_c] = *depth; + years[event_c] = *yyyy; months[event_c] = *mm; + days[event_c] = *dd; hours[event_c] = *hh; + minutes[event_c] = *mi; + if (*gap != ISF_NULL && *gap < gap_min) { + gap_min = *gap; + *idx_min_rms = event_c; + } + event_c++; + } + else + done = true; + } + return event_c; +} + +int read_mags(FILE *fp, char *line, char *magtype, char* magind, float* mag, float* magerr, + int* nsta, char* author, char* origid, char **mag_t, float *mags) { + + int done = false, mag_c = 0; + + while (!done && (fgets (line, ISF_LINE_LEN, fp) != NULL)) { + if(!read_netmag(line,magtype,magind,mag,magerr,nsta,author,origid)) { + sscanf (magtype, "%s",mag_t[mag_c]); + mags[mag_c] = *mag; + mag_c++; + } + else + done = true; + } + return(mag_c); +} + +float select_mag(int mag_c, float *mags, char **mag_t) { + /* Select among the following magnitudes in this orther */ + int i; + float mag; + for (i = 0; i < mag_c; i++) { + gmt_str_toupper (mag_t[i]); + if (!strncmp(mag_t[i],"MW",2)) { + mag = mags[i]; + break; + } + else if (!strncmp(mag_t[i],"MB",2)) { + mag = mags[i]; + break; + } + else if (!strncmp(mag_t[i],"MS",2)) { + mag = mags[i]; + break; + } + else if (!strncmp(mag_t[i],"MD",2)) { + mag = mags[i]; + break; + } + else if (!strncmp(mag_t[i],"ML",2)) { + mag = mags[i]; + break; + } + else /* A magnitude not in this list */ + mag = mags[0]; + } + return (mag); +} + +/* Parses a line asuming it to be an event title line. + Requires event ID to be present but allows lines with no region. + + Returns 0 if the line is a properly formatted event ID line. + Returns 20 and writes a diagnostic to isf_error on error. +*/ +int read_event_id(char *line, char *evid, char *region) { + char substr[ISF_LINE_LEN]; + + /* Chars 1-5: must be the word 'Event'. Char 6: must be a space. */ + if (strncmp(line,"Event ",6) && strncmp(line,"EVENT ",6)){ + sprintf(isf_error,"not an event title line: %s",line); + return 20; + } + + /* Chars 7-14: event ID */ + if (!partline(evid,line,6,8)){ + sprintf(isf_error,"missing evid: %s",line); + return 20; + } + if (check_whole(evid)){ + sprintf(isf_error,"bad evid: %s",line); + return 20; + } + + /* Not quite right but lots of people hit CR after evid */ + if (strlen(line) < 15) return 0; + + /* Char 15: must be a space */ + if (line[14] != ' ' ){ + sprintf(isf_error,"bad format, char 15: %s",line); + return 20; + } + + /* Chars 16-80: geographic region if there */ + partline(region,line,15,65); + + /* Check for extra characters after char 80. */ + if (partline(substr,line,80,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + + return 0; +} + + +/* Tests a line to discover if it is an origin header line. + + Returns 0 if the line is a properly formatted origin header line. + Returns 20 and writes a diagnostic to isf_error otherwise. +*/ +int read_origin_head(char *line) { + char substr[ISF_LINE_LEN]; + char head[] = " Date Time Err RMS Latitude Longitude Smaj Smin Az Depth Err Ndef Nsta Gap mdist Mdist Qual Author OrigID"; + int headlen = 136; + + if (strncmp(line,head,headlen) != 0){ + sprintf(isf_error,"not an origin header: %s",line); + return 20; + } + + /* Check for extra characters after char 136. */ + if (partline(substr,line,headlen,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + + +/* Parses a line asuming it to be an origin line. + Values are asigned to variables, the pointers to which have been sent + as arguments. If an optional parameter is not given then the + corresponding variable will have ISF_NULL assigned to it. + + Returns 0 if the line is a properly formatted origin line. + Returns 20 and writes a diagnostic to isf_error on error. +*/ +int read_origin(char *line, int *yyyy, int *mm, int *dd, int *hh, int *mi, + int *ss, int *msec, char *timfix, float *stime, float *sdobs, + float *lat, float *lon, char *epifix, float *smaj, float *smin, + int *strike, float *depth, char *depfix, float *sdepth, + int *ndef, int *nsta, int *gap, float *mindist, + float *maxdist, char *antype, char *loctype, char *etype, + char *author, char *origid) { + + char substr[ISF_LINE_LEN]; + + /* Chars 1-4: year. */ + if (!partline(substr,line,0,4)){ + sprintf(isf_error,"missing year: %s",line); + return 20; + } + if (check_int(substr)){ + sprintf(isf_error,"bad year: %s",line); + return 20; + } + *yyyy = atoi(substr); + + /* Char 5: '/' character. */ + if (line[4] != '/'){ + sprintf(isf_error,"bad date: %s",line); + return 20; + } + + /* Chars 6-7: month. */ + if (!partline(substr,line,5,2)){ + sprintf(isf_error,"missing month: %s",line); + return 20; + } + if (check_int(substr)){ + sprintf(isf_error,"bad month: %s",line); + return 20; + } + *mm = atoi(substr); + + /* Char 8: '/' character. */ + if (line[7] != '/'){ + sprintf(isf_error,"bad date: %s",line); + return 20; + } + + /* Chars 9-10: day. */ + if (!partline(substr,line,8,2)){ + sprintf(isf_error,"missing day: %s",line); + return 20; + } + if (check_int(substr)){ + sprintf(isf_error,"bad day: %s",line); + return 20; + } + *dd = atoi(substr); + + /* Char 11: space. */ + if (line[10] != ' '){ + sprintf(isf_error,"bad date: %s",line); + return 20; + } + + /* Chars 12,13: hour. */ + if (!partline(substr,line,11,2)){ + sprintf(isf_error,"missing hour: %s",line); + return 20; + } + if (check_int(substr)){ + sprintf(isf_error,"bad hour: %s",line); + return 20; + } + *hh = atoi(substr); + + /* Char 14: ':' character. */ + if (line[13] != ':'){ + sprintf(isf_error,"bad date: %s",line); + return 20; + } + + /* Chars 15,16: minute. */ + if (!partline(substr,line,14,2)){ + sprintf(isf_error,"missing minute: %s",line); + return 20; + } + if (check_int(substr)){ + sprintf(isf_error,"bad minute: %s",line); + return 20; + } + *mi = atoi(substr); + + /* Char 17: ':' character. */ + if (line[16] != ':'){ + sprintf(isf_error,"bad date: %s",line); + return 20; + } + + /* Chars 18,19: integral second. */ + if (!partline(substr,line,17,2)){ + sprintf(isf_error,"missing second: %s",line); + return 20; + } + if (check_int(substr)){ + sprintf(isf_error,"bad second: %s",line); + return 20; + } + *ss = atoi(substr); + + /* Char 20-22: msec or spaces. */ + /* Allow decimal place with no numbers after it. */ + if (partline(substr,line,20,2)){ + /* Char 20: '.' character */ + if (line[19] != '.'){ + sprintf(isf_error,"bad date: %s",line); + return 20; + } + /* Chars 21,22: 10s of msec. */ + if (!isdigit(line[20])){ + sprintf(isf_error,"bad date: %s",line); + return 20; + } + *msec = (line[20]-'0')*100; + + if (isdigit(line[21])){ + *msec += (line[21]-'0')*10; + } + else if (line[21] != ' ') { + sprintf(isf_error,"bad date: %s",line); + return 20; + } + } + else { + /* Char 20: '.' character or space */ + if (line[19] != '.' && line[19] != ' '){ + sprintf(isf_error,"bad date: %s",line); + return 20; + } + *msec = ISF_NULL; + } + + /* Char 23: timfix - either f or space. */ + if (line[22] == ' ' || line[22] == 'f'){ + *timfix = line[22]; + } + else { + sprintf(isf_error,"bad timfix: %s",line); + return 20; + } + + /* Char 24: space. */ + if (line[23] != ' '){ + sprintf(isf_error,"bad format, char 24: %s",line); + return 20; + } + + /* Chars 25-29: origin time error - float if anything. */ + if (partline(substr,line,24,5)){ + if (check_float(substr)){ + sprintf(isf_error,"bad stime: %s",line); + return 20; + } + *stime = (float)atof(substr); + } + else { + *stime = ISF_NULL; + } + + /* Char 30: space. */ + if (line[29] != ' '){ + sprintf(isf_error,"bad format, char 30: %s",line); + return 20; + } + + /* Chars 31-35: rms (sdobs) - float if anything. */ + if (partline(substr,line,30,5)){ + if (check_float(substr)){ + sprintf(isf_error,"bad sdobs: %s",line); + return 20; + } + *sdobs = (float)atof(substr); + } + else { + *sdobs = ISF_NULL; + } + + /* Char 36: space. */ + if (line[35] != ' '){ + sprintf(isf_error,"bad format, char 36: %s",line); + return 20; + } + + /* Chars 37-44: lattitude - must be there. */ + if (!partline(substr,line,36,8)){ + sprintf(isf_error,"missing lattitude: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad lattitude: %s",line); + return 20; + } + *lat = (float)atof(substr); + + /* Char 45: space. */ + if (line[44] != ' '){ + sprintf(isf_error,"bad format, char 45: %s",line); + return 20; + } + + /* Chars 46-54: longitude - must be float. */ + if (!partline(substr,line,45,9)){ + sprintf(isf_error,"missing longitude: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad longitude: %s",line); + return 20; + } + *lon = (float)atof(substr); + + /* Char 55: epifix - either f or space. */ + if (line[54] == ' ' || line[54] == 'f'){ + *epifix = line[54]; + } + else { + sprintf(isf_error,"bad epifix: %s",line); + return 20; + } + + /* Chars 56-60: semi-major axis of error ellipse - float if there. */ + /* This is departure from format but smaj < smin is daft. */ + if (partline(substr,line,55,5)){ + if (check_float(substr)){ + sprintf(isf_error,"bad smaj: %s",line); + return 20; + } + *smaj = (float)atof(substr); + } + else { + *smaj = ISF_NULL; + } + + /* Char 61: space. */ + if (line[60] != ' '){ + sprintf(isf_error,"bad format, char 61: %s",line); + return 20; + } + + /* Chars 62-66: semi-minor axis of error ellipse - float if there. */ + if (partline(substr,line,61,5)){ + if (check_float(substr)){ + sprintf(isf_error,"bad smin: %s",line); + return 20; + } + *smin = (float)atof(substr); + } + else { + *smin = ISF_NULL; + } + + /* Char 67: space. */ + if (line[66] != ' '){ + sprintf(isf_error,"bad format, char 67: %s",line); + return 20; + } + + /* Chars 68-70: strike - integer if there. */ + if (partline(substr,line,67,3)){ + if (check_int(substr)){ + sprintf(isf_error,"bad strike: %s",line); + return 20; + } + *strike = atoi(substr); + } + else { + *strike = ISF_NULL; + } + + /* Char 71: space. */ + if (line[70] != ' '){ + sprintf(isf_error,"bad format, char 71: %s",line); + return 20; + } + + /* Chars 72-76: depth - float if there. */ + if (partline(substr,line,71,5)){ + if (check_float(substr)){ + sprintf(isf_error,"bad depth: %s",line); + return 20; + } + *depth = (float)atof(substr); + } + else { + *depth = ISF_NULL; + } + + /* Char 77: depfix - either f,d, or space. */ + if (line[76] == ' ' || line[76] == 'f' || line[76] == 'd'){ + *depfix = line[76]; + } + else { + sprintf(isf_error,"bad depfix: %s",line); + return 20; + } + + /* Char 78: space. */ + if (line[77] != ' '){ + sprintf(isf_error,"bad format, char 78: %s",line); + return 20; + } + + /* Chars 79-82: depth error - float if there. */ + if (partline(substr,line,78,4)){ + if (check_float(substr)){ + sprintf(isf_error,"bad sdepth: %s",line); + return 20; + } + *sdepth = (float)atof(substr); + } + else { + *sdepth = ISF_NULL; + } + + /* Char 83: space. */ + if (line[82] != ' '){ + sprintf(isf_error,"bad format, char 83: %s",line); + return 20; + } + + /* Chars 84-87: ndef - integer if there. */ + if (partline(substr,line,83,4)){ + if (check_int(substr)){ + sprintf(isf_error,"bad ndef: %s",line); + return 20; + } + *ndef = atoi(substr); + } + else { + *ndef = ISF_NULL; + } + + /* Char 88: space. */ + if (line[87] != ' '){ + sprintf(isf_error,"bad format, char 88: %s",line); + return 20; + } + + /* Chars 89-92: nsta - integer if there. */ + if (partline(substr,line,88,4)){ + if (check_int(substr)){ + sprintf(isf_error,"bad nsta: %s",line); + return 20; + } + *nsta = atoi(substr); + } + else { + *nsta = ISF_NULL; + } + + /* Char 93: space. */ + if (line[92] != ' '){ + sprintf(isf_error,"bad format, char 93: %s",line); + return 20; + } + + /* Chars 94-96: gap - integer if there */ + if (partline(substr,line,93,3)){ + if (check_int(substr)){ + sprintf(isf_error,"bad gap: %s",line); + return 20; + } + *gap = atoi(substr); + } + else { + *gap = ISF_NULL; + } + + /* Char 97: space. */ + if (line[96] != ' '){ + sprintf(isf_error,"bad format, char 97: %s",line); + return 20; + } + + /* Chars 98-103: minimum distance - float if there. */ + if (partline(substr,line,97,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mindist: %s",line); + return 20; + } + *mindist = (float)atof(substr); + } + else { + *mindist = ISF_NULL; + } + + /* Char 104: space. */ + if (line[103] != ' '){ + sprintf(isf_error,"bad format, char 104: %s",line); + return 20; + } + + /* Chars 105-110: maximum distance - float if there. */ + if (partline(substr,line,104,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad maxdist: %s",line); + return 20; + } + *maxdist = (float)atof(substr); + } + else { + *maxdist = ISF_NULL; + } + + /* Char 111: space. */ + if (line[110] != ' '){ + sprintf(isf_error,"bad format, char 111: %s",line); + return 20; + } + + /* Char 112: analysis type - either space, a, m, or g. */ + if (line[111] == ' ' || line[111] == 'a' || line[111] == 'm' || \ + line[111] =='g' ){ + *antype = line[111]; + } + else { + sprintf(isf_error,"bad antype: %s",line); + return 20; + } + + /* Char 113: space. */ + if (line[112] != ' '){ + sprintf(isf_error,"bad format, char 113: %s",line); + return 20; + } + + /* Char 114: location method - either space, i, p, g, or o. */ + if (line[113] == ' ' || line[113] == 'i' || line[113] == 'p' || \ + line[113] =='g' || line[113] == 'o'){ + *loctype = line[113]; + } + else { + sprintf(isf_error,"bad loctype: %s",line); + return 20; + } + + /* Char 115: space. */ + if (line[114] != ' '){ + sprintf(isf_error,"bad format, char 115: %s",line); + return 20; + } + + /* Chars 116-117: event type, any characters allowed but must be there. */ + if (!partline(etype,line,115,2)){ + strcpy(etype,""); + } + else if (strlen(etype) != 2){ + sprintf(isf_error,"bad etype: %s",line); + return 20; + } + + /* Char 118: space. */ + if (line[117] != ' '){ + sprintf(isf_error,"bad format, char 118: %s",line); + return 20; + } + + /* Chars 119-127: author, any characters allowed but must be there. */ + if (!partline(author,line,118,9)){ + sprintf(isf_error,"missing author: %s",line); + return 20; + } + if (check_whole(author)){ + sprintf(isf_error,"bad author: %s",line); + return 20; + } + + /* Char 128: space */ + if (line[127] != ' '){ + sprintf(isf_error,"bad format, char 128: %s",line); + return 20; + } + + /* Chars 129-136: origin ID, any characters allowed but must be there. */ + if (!partline(origid,line,128,8)){ + sprintf(isf_error,"missing origid: %s",line); + return 20; + } + if (check_whole(origid)){ + sprintf(isf_error,"bad origid: %s",line); + return 20; + } + + /* Check for extra stuff after char 136. */ + if (partline(substr,line,136,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* To check that a line is a good magnitude block header line. + + Returns 0 if the line is a magnitude block header. + Returns 20 and writes a diagnostic to isf_error otherwise. +*/ +int read_netmag_head(char *line) { + char substr[ISF_LINE_LEN]; + + char head[] = "Magnitude Err Nsta Author OrigID"; + int headlen = 38; + + if (strncmp(line,head,headlen) != 0){ + sprintf(isf_error,"not a netmag header: %s",line); + return 20; + } + if (partline(substr,line,headlen,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Parses a line assuming that it is a magnitude sub-block data line. + Values are asigned to variables, the pointers to which have been sent + as arguments. If an optional parameter is not given then the + corresponding variable will have ISF_NULL assigned to it. + + Returns 0 if the line is a properly formatted magnitude line, + Returns 20 and writes a diagnostic to isf_error otherwise. +*/ +int read_netmag( char *line, char *magtype, char* magind, float* mag, + float* magerr, int* nsta, char* author, char* origid) +{ + char substr[ISF_LINE_LEN]; + + /* Chars 1-5: magnitude type, any characters allowed but must be there. */ + if (!partline(magtype,line,0,5)){ + sprintf(isf_error,"missing magtype: %s",line); + return 20; + } + if (check_whole(magtype)){ + sprintf(isf_error,"bad magtype: %s",line); + return 20; + } + + /* Char 6: less than or greater than indicator or space only. */ + if (line[5] == ' ' || line[5] == '<' || line[5] == '>'){ + *magind = line[5]; + } + else { + sprintf(isf_error,"bad magind: %s",line); + return 20; + } + + /* Chars 7-10: magnitude, must be float. */ + if (!partline(substr,line,6,4)){ + sprintf(isf_error,"missing magnitude: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad magnitude: %s",line); + return 20; + } + *mag = (float)atof(substr); + + /* Char 11: must be a space. */ + if (line[10] != ' ' ){ + sprintf(isf_error,"bad format, char 11: %s",line); + return 20; + } + + /* Chars 12-14: magnitude error, float if anything. */ + if (partline(substr,line,11,3)){ + if (check_float(substr)){ + sprintf(isf_error,"bad magnitude error: %s",line); + return 20; + } + *magerr = (float)atof(substr); + } + else { + *magerr = ISF_NULL; + } + + /* Char 15: must be a space. */ + if (line[14] != ' ' ){ + sprintf(isf_error,"bad format, char 15: %s",line); + return 20; + } + + /* Chars 16-19: number of stations, integer if anything. */ + if (partline(substr,line,15,4)){ + if (check_float(substr)){ + sprintf(isf_error,"bad nsta: %s",line); + return 20; + } + *nsta = atoi(substr); + } + else { + *nsta = ISF_NULL; + } + + /* Char 20: must be a space. */ + if (line[19] != ' ' ){ + sprintf(isf_error,"bad format, char 20: %s",line); + return 20; + } + + /* Chars 21-29: author, any characters allowed but must be there. */ + if (!partline(author,line,20,9)){ + sprintf(isf_error,"missing author: %s",line); + return 20; + } + if (check_whole(author)){ + sprintf(isf_error,"bad author: %s",line); + return 20; + } + + /* Char 30: must be a space. */ + if (line[29] != ' ' ){ + sprintf(isf_error,"bad format, char 30: %s",line); + return 20; + } + + /* Chars 31-38: origin ID, any characters allowed but must be there. */ + if (!partline(origid,line,30,8)){ + sprintf(isf_error,"missing origid: %s",line); + return 20; + } + if (check_whole(origid)){ + sprintf(isf_error,"bad origid: %s",line); + return 20; + } + + /* Check for extra stuff after char 38. */ + if (partline(substr,line,38,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Parses a line to test whether it is a centroid origin label. + + Returns 0 if the line is a properly formatted centroid origin line. + Returns 20 and writes a diagnostic to isf_error if not. +*/ +int read_origin_centroid(char *line) { + char substr[ISF_LINE_LEN]; + + if (strncmp(line," (#CENTROID)",12)){ + sprintf(isf_error,"not a centroid comment: %s",line); + return 20; + } + if (partline(substr,line,13,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Tests a line to discover if it is a first moment tensor header comment. + + Returns 0 if the line is a first moment tensor header line. + Returns 20 and writes a diagnostic to isf_error otherwise. +*/ +int read_momten_head_1(char *line) { + char substr[ISF_LINE_LEN]; + char head[] = " (#MOMTENS sc M0 fCLVD MRR MTT MPP MRT MTP MPR NST1 NST2 Author )"; + int headlen = 88; + + if (strncmp(line,head,headlen) != 0){ + sprintf(isf_error,"not a momten header: %s",line); + return 20; + } + if (partline(substr,line,headlen,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + + +/* Tests a line to discover if it is a second moment tensor header comment. + + Returns 0 if the line is a second moment tensor header line. + Returns 20 and writes a diagnostic to isf_error otherwise. +*/ +int read_momten_head_2(char *line) { + char substr[ISF_LINE_LEN]; + char head[] = " (# eM0 eCLVD eRR eTT ePP eRT eTP ePR NCO1 NCO2 Duration )"; + int headlen = 88; + + if (strncmp(line,head,headlen) != 0){ + sprintf(isf_error,"not a momten header2: %s",line); + return 20; + } + if (partline(substr,line,headlen,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Parses a line asuming it to be a first moment tensor data comment. + Values are asigned to variables, the pointers to which have been sent + as arguments. If an optional parameter is not given then the + corresponding variable will have ISF_NULL assigned to it. + + Returns 0 if the line is a properly formatted first moment tensor data line. + Returns 20 and writes a diagnostic to isf_error on error. +*/ +int read_momten_line_1(char *line, int *scale_factor, float *scalar_moment, + float *fclvd, float *mrr, float *mtt, float *mpp, + float *mrt, float *mtp, float *mpr, int *nsta1, + int *nsta2, char *author) +{ + char substr[ISF_LINE_LEN]; + + /* Chars 1-11: should be the string ' (# ' */ + if (strncmp(line," (# ",11) != 0){ + sprintf(isf_error,"not a moment tensor line: %s",line); + return 20; + } + + /* Chars 12,13: scale factor - integer */ + if (!partline(substr,line,11,2)){ + sprintf(isf_error,"missing scale factor: %s",line); + return 20; + } + if (check_int(substr)){ + sprintf(isf_error,"bad scale factor: %s",line); + return 20; + } + *scale_factor = atoi(substr); + + /* Char 14: must be a space */ + if (line[13] != ' ' ){ + sprintf(isf_error,"bad format, char 14: %s",line); + return 20; + } + + /* Chars 15-19: scalar seismic moment - must be float. */ + if (!partline(substr,line,14,5)){ + sprintf(isf_error,"missing moment: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad moment: %s",line); + return 20; + } + *scalar_moment = (float)atof(substr); + + /* Char 20: must be a space */ + if (line[19] != ' ' ){ + sprintf(isf_error,"bad format, char 20: %s",line); + return 20; + } + + /* Chars 21-25: fCLVD, float if anything */ + if (partline(substr,line,20,5)){ + if (check_float(substr)){ + sprintf(isf_error,"bad fclvd: %s",line); + return 20; + } + *fclvd = (float)atof(substr); + } + else { + *fclvd = ISF_NULL; + } + + /* Char 26: must be a space */ + if (line[25] != ' ' ){ + sprintf(isf_error,"bad format, char 26: %s",line); + return 20; + } + + /* Chars 27-32: radial-radial element, float if anything */ + if (partline(substr,line,26,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mrr: %s",line); + return 20; + } + *mrr = (float)atof(substr); + } + else { + *mrr = ISF_NULL; + } + + /* Char 33: must be a space */ + if (line[32] != ' ' ){ + sprintf(isf_error,"bad format, char 33: %s",line); + return 20; + } + + /* Chars 34-39: theta-theta element, float if anything */ + if (partline(substr,line,33,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mtt: %s",line); + return 20; + } + *mtt = (float)atof(substr); + } + else { + *mtt = ISF_NULL; + } + + /* Char 40: must be a space */ + if (line[39] != ' ' ){ + sprintf(isf_error,"bad format, char 40: %s",line); + return 20; + } + + /* Chars 41-46: phi-phi element, float if anything */ + if (partline(substr,line,40,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mpp: %s",line); + return 20; + } + *mpp = (float)atof(substr); + } + else { + *mpp = ISF_NULL; + } + + /* Char 47: must be a space */ + if (line[46] != ' ' ){ + sprintf(isf_error,"bad format, char 47: %s",line); + return 20; + } + + /* Chars 48-53: radial-theta element, float if anything */ + if (partline(substr,line,47,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mrt: %s",line); + return 20; + } + *mrt = (float)atof(substr); + } + else { + *mrt = ISF_NULL; + } + + /* Char 54: must be a space */ + if (line[53] != ' ' ){ + sprintf(isf_error,"bad format, char 54: %s",line); + return 20; + } + + /* Chars 55-60: theta-phi element, float if anything */ + if (partline(substr,line,54,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mtp: %s",line); + return 20; + } + *mtp = (float)atof(substr); + } + else { + *mtp = ISF_NULL; + } + + /* Char 61: must be a space */ + if (line[60] != ' ' ){ + sprintf(isf_error,"bad format, char 61: %s",line); + return 20; + } + + /* Chars 62-67: phi-radial element, float if anything */ + if (partline(substr,line,61,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mpr: %s",line); + return 20; + } + *mpr = (float)atof(substr); + } + else { + *mpr = ISF_NULL; + } + + /* Char 68: must be a space */ + if (line[67] != ' ' ){ + sprintf(isf_error,"bad format, char 68: %s",line); + return 20; + } + + /* Chars 69-72: nsta1, int if anything */ + if (partline(substr,line,68,4)){ + if (check_int(substr)){ + sprintf(isf_error,"bad nsta1: %s",line); + return 20; + } + *nsta1 = atoi(substr); + } + else { + *nsta1 = ISF_NULL; + } + + /* Char 73: must be a space */ + if (line[72] != ' ' ){ + sprintf(isf_error,"bad format, char 73: %s",line); + return 20; + } + + /* Chars 74-77: nsta2, int if anything */ + if (partline(substr,line,73,4)){ + if (check_int(substr)){ + sprintf(isf_error,"bad nsta2: %s",line); + return 20; + } + *nsta2 = atoi(substr); + } + else { + *nsta2 = ISF_NULL; + } + + /* Char 78: must be a space */ + if (line[77] != ' ' ){ + sprintf(isf_error,"bad format, char 78: %s",line); + return 20; + } + + /* Chars 79-87: author, any characters allowed but must be there. */ + if (!partline(author,line,78,9)){ + sprintf(isf_error,"missing author: %s",line); + return 20; + } + if (check_whole(author)){ + sprintf(isf_error,"bad author: %s",line); + return 20; + } + + /* Check for extra characters - could be close bracket somewhere. */ + if (partline(substr,line,87,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Parses a line asuming it to be a second moment tensor data comment. + Values are asigned to variables, the pointers to which have been sent + as arguments. If an optional parameter is not given then the + corresponding variable will have ISF_NULL assigned to it. + + Returns 0 if the line is properly formatted second moment tensor data. + Returns 20 and writes a diagnostic to isf_error on error. +*/ +int read_momten_line_2(char *line, float *scalar_moment_unc, float *fclvd_unc, + float *mrr_unc, float *mtt_unc, float *mpp_unc, + float *mrt_unc, float *mtp_unc, float *mpr_unc, + int *ncomp1, int *ncomp2, float *duration) +{ + char substr[ISF_LINE_LEN]; + + /* Chars 1-14: should be the string ' (# '. */ + if (strncmp(line," (# ",14) != 0){ + sprintf(isf_error,"not a moment tensor line: %s",line); + return 20; + } + + /* Chars 15-19: uncertainty in scalar seismic moment - float if there. */ + if (partline(substr,line,14,5)){ + if (check_float(substr)){ + sprintf(isf_error,"bad scalar_moment_unc: %s",line); + return 20; + } + *scalar_moment_unc = (float)atof(substr); + } + else { + *scalar_moment_unc = ISF_NULL; + } + + /* Char 20: must be a space. */ + if (line[19] != ' ' ){ + sprintf(isf_error,"bad format, char 20: %s",line); + return 20; + } + + /* Chars 21-25: uncertainty in fCLVD, float if anything. */ + if (partline(substr,line,20,5)){ + if (check_float(substr)){ + sprintf(isf_error,"bad fclvd_unc: %s",line); + return 20; + } + *fclvd_unc = (float)atof(substr); + } + else { + *fclvd_unc = ISF_NULL; + } + + /* Char 26: must be a space. */ + if (line[25] != ' ' ){ + sprintf(isf_error,"bad format, char 26: %s",line); + return 20; + } + + /* Chars 27-32: uncertainty in radial-radial element, float if anything. */ + if (partline(substr,line,26,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mrr_unc: %s",line); + return 20; + } + *mrr_unc = (float)atof(substr); + } + else { + *mrr_unc = ISF_NULL; + } + + /* Char 33: must be a space. */ + if (line[32] != ' ' ){ + sprintf(isf_error,"bad format, char 33: %s",line); + return 20; + } + + /* Chars 34-39: uncertainty in theta-theta element, float if anything. */ + if (partline(substr,line,33,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mtt_unc: %s",line); + return 20; + } + *mtt_unc = (float)atof(substr); + } + else { + *mtt_unc = ISF_NULL; + } + + /* Char 40: must be a space */ + if (line[39] != ' ' ){ + sprintf(isf_error,"bad format, char 40: %s",line); + return 20; + } + + /* Chars 41-46: uncertainty in phi-phi element, float if anything. */ + if (partline(substr,line,40,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mpp_unc: %s",line); + return 20; + } + *mpp_unc = (float)atof(substr); + } + else { + *mpp_unc = ISF_NULL; + } + + /* Char 47: must be a space */ + if (line[46] != ' ' ){ + sprintf(isf_error,"bad format, char 47: %s",line); + return 20; + } + + /* Chars 48-53: uncertainty in radial-theta element, float if anything. */ + if (partline(substr,line,47,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mrt_unc: %s",line); + return 20; + } + *mrt_unc = (float)atof(substr); + } + else { + *mrt_unc = ISF_NULL; + } + + /* Char 54: must be a space. */ + if (line[53] != ' ' ){ + sprintf(isf_error,"bad format: %s",line); + return 20; + } + + /* Chars 55-60: uncertainty in theta-phi element, float if anything. */ + if (partline(substr,line,54,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mtp_unc: %s",line); + return 20; + } + *mtp_unc = (float)atof(substr); + } + else { + *mtp_unc = ISF_NULL; + } + + /* Char 61: must be a space. */ + if (line[60] != ' ' ){ + sprintf(isf_error,"bad format, char 61: %s",line); + return 20; + } + + /* Chars 62-67: uncertainty in phi-radial element, float if anything. */ + if (partline(substr,line,61,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad mpr_unc: %s",line); + return 20; + } + *mpr_unc = (float)atof(substr); + } + else { + *mpr_unc = ISF_NULL; + } + + /* Char 68: must be a space. */ + if (line[67] != ' ' ){ + sprintf(isf_error,"bad format, char 68: %s",line); + return 20; + } + + /* Chars 69-72: ncomp1, int if anything. */ + if (partline(substr,line,68,4)){ + if (check_int(substr)){ + sprintf(isf_error,"bad ncomp1: %s",line); + return 20; + } + *ncomp1 = atoi(substr); + } + else { + *ncomp1 = ISF_NULL; + } + + /* Char 73: must be a space */ + if (line[72] != ' ' ){ + sprintf(isf_error,"bad format, char 73: %s",line); + return 20; + } + + /* Chars 74-77: ncomp2, int if anything. */ + if (partline(substr,line,73,4)){ + if (check_int(substr)){ + sprintf(isf_error,"bad ncomp2: %s",line); + return 20; + } + *ncomp2 = atoi(substr); + } + else { + *ncomp2 = ISF_NULL; + } + + /* Char 78: must be a space. */ + if (line[77] != ' ' ){ + sprintf(isf_error,"bad format, char 78: %s",line); + return 20; + } + + /* Chars 79-86: duration, float if anything. */ + if (partline(substr,line,78,8)){ + if (check_float(substr)){ + sprintf(isf_error,"bad duration: %s",line); + return 20; + } + *duration = (float)atof(substr); + } + else { + *duration = ISF_NULL; + } + + /* Check for extra stuff - not including close bracket. */ + if (partline(substr,line,86,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Tests a line to discover if it is a fault plane header comment. + + Returns 0 if the line is a fault plane header. + Returns 20 and writes a diagnostic to isf_error otherwise. +*/ +int read_fault_plane_head(char *line) +{ + char substr[ISF_LINE_LEN]; + char head[] = " (#FAULT_PLANE Typ Strike Dip Rake NP NS Plane Author )"; + int headlen = 64; + + if (strncmp(line,head,headlen) != 0){ + sprintf(isf_error,"not a fault plane header: %s",line); + return 20; + } + if (partline(substr,line,headlen,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Parses a line asuming it to be a fault plane data comment. + Could be first or second plane, the only difference is whether + author field is expected or not. + Values are asigned to variables, the pointers to which have been sent + as arguments. If an optional parameter is not given then the + corresponding variable will have ISF_NULL assigned to it. + + Returns 0 if the line is a properly formatted fault plane data line. + Returns 20 and writes a diagnostic to isf_error on error. +*/ +int read_fault_plane (char *line, char *f_type, float *strike, float *dip, + float *rake, int *np, int *ns, char *f_plane, char *author) +{ + char substr[ISF_LINE_LEN]; + int line_num; + + /* Chars 1-3: the strings ' (#' or ' (+', */ + /* depending on whether this is the first or second plane given. */ + /* Chars 4-15: spaces. */ + if (strncmp(line," (# ",15) == 0){ + line_num = 1; + } + else if (strncmp(line," (+ ",15) == 0){ + line_num = 2; + } + else { + sprintf(isf_error,"not a fault plane line: %s",line); + return 20; + } + + /* Chars 16-18: fault plane solution type. */ + if (partline(f_type,line,15,3)){ + if (strcmp(f_type,"FM") && strcmp(f_type,"BB") && strcmp(f_type,"BDC")){ + sprintf(isf_error,"bad f_type: %s",line); + return 20; + } + } + else { + strcpy(f_type,""); + } + + /* Char 19: must be a space */ + if (line[18] != ' ' ){ + sprintf(isf_error,"bad format, char 19: %s",line); + return 20; + } + + /* Chars 20-25: strike, must be float. */ + if (!partline(substr,line,19,6)){ + sprintf(isf_error,"missing strike: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad strike: %s",line); + return 20; + } + *strike = (float)atof(substr); + + /* Char 26: must be a space. */ + if (line[25] != ' ' ){ + sprintf(isf_error,"bad format, char 26: %s",line); + return 20; + } + + /* Chars 27-31: dip, must be float. */ + if (!partline(substr,line,26,5)){ + sprintf(isf_error,"missing dip: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad dip: %s",line); + return 20; + } + *dip = (float)atof(substr); + + /* Char 32: must be a space. */ + if (line[31] != ' ' ){ + sprintf(isf_error,"bad format, char 32: %s",line); + return 20; + } + + /* Chars 33-39: rake, float - need not be there if both planes given. */ + if (partline(substr,line,32,7)){ + if (check_float(substr)){ + sprintf(isf_error,"bad rake: %s",line); + return 20; + } + *rake = (float)atof(substr); + } + else { + *rake = ISF_NULL; + } + + /* Char 40: must be a space. */ + if (line[39] != ' ' ){ + sprintf(isf_error,"bad format, char 40: %s",line); + return 20; + } + + /* Chars 41-43: np, int if there. */ + if (partline(substr,line,40,3)){ + if (check_int(substr)){ + sprintf(isf_error,"bad np: %s",line); + return 20; + } + *np = atoi(substr); + } + else { + *np = ISF_NULL; + } + + /* Char 44: must be a space. */ + if (line[43] != ' ' ){ + sprintf(isf_error,"bad format, char 44: %s",line); + return 20; + } + + /* Chars 45-47: ns, int if there. */ + if (partline(substr,line,44,3)){ + if (check_int(substr)){ + sprintf(isf_error,"bad ns: %s",line); + return 20; + } + *ns = atoi(substr); + } + else { + *ns = ISF_NULL; + } + + /* Char 48: must be a space. */ + if (line[47] != ' ' ){ + sprintf(isf_error,"bad format, char 48: %s",line); + return 20; + } + + /* Chars 49-53: plane identification. */ + if (partline(f_plane,line,48,5)){ + if (strcmp(f_plane,"FAULT") != 0 && strcmp(f_plane,"AUXIL") != 0){ + sprintf(isf_error,"bad f_plane: %s",line); + return 20; + } + } + else { + strcpy(f_plane,""); + } + + /* Chars 54-63: First plane has author, don't read for second plane. */ + if (line_num==1){ + + /* Char 54: must be a space */ + if (line[53] != ' ' ){ + sprintf(isf_error,"bad format, char 54: %s",line); + return 20; + } + + /* Chars 55-63: author, any characters allowed but must be there */ + if (!partline(author,line,54,9)){ + sprintf(isf_error,"missing author: %s",line); + return 20; + } + if (check_whole(author)){ + sprintf(isf_error,"bad author: %s",line); + return 20; + } + } + + /* Check for extra stuff - not including close bracket */ + if (partline(substr,line,63,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Tests a line to discover if it is a principal axes header comment. + + Returns 0 if the line is a principal axes header. + Returns 20 and writes a diagnostic to isf_error otherwise. +*/ +int read_axes_head(char *line) { + char substr[ISF_LINE_LEN]; + char head[] = " (#PRINAX sc T_val T_azim T_pl B_val B_azim B_pl P_val P_azim P_pl Author )"; + int headlen = 83; + + if (strncmp(line,head,headlen) != 0){ + sprintf(isf_error,"not an axes header: %s",line); + return 20; + } + if (partline(substr,line,headlen,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Parses a line asuming it to be a principal axes data comment. + Values are asigned to variables, the pointers to which have been sent + as arguments. If an optional parameter is not given then the + corresponding variable will have ISF_NULL assigned to it. + + Returns 0 if the line is a properly formatted principal axes data line. + Returns 20 and writes a diagnostic to isf_error on error. +*/ +int read_axes(char *line, int *scale_factor, float *t_val, float *t_azim, + float *t_pl, float *b_val, float *b_azim, float *b_pl, + float *p_val, float *p_azim, float *p_pl, char *author) +{ + + char substr[ISF_LINE_LEN]; + + /* Chars 1-10: should be the string ' (# '. */ + if (strncmp(line," (# ",10) != 0){ + sprintf(isf_error,"not an axes line: %s",line); + return 20; + } + + /* Chars 11,12: scale factor - int if there. */ + if (partline(substr,line,10,2)){ + if (check_int(substr)){ + sprintf(isf_error,"bad scale factor: %s",line); + return 20; + } + *scale_factor = atoi(substr); + } + else { + *scale_factor = ISF_NULL; + } + + /* Char 13: must be a space. */ + if (line[12] != ' ' ){ + sprintf(isf_error,"bad format, char 13: %s",line); + return 20; + } + + /* Chars 14-19: t value - float if there. */ + if (partline(substr,line,13,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad t_val: %s",line); + return 20; + } + *t_val = (float)atof(substr); + } + else { + *t_val = ISF_NULL; + } + + /* Char 20: must be a space. */ + if (line[19] != ' ' ){ + sprintf(isf_error,"bad format, char 20: %s",line); + return 20; + } + + /* Chars 21-26: t azim, must be float. */ + if (!partline(substr,line,20,6)){ + sprintf(isf_error,"missing t_azim: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad t_azim: %s",line); + return 20; + } + *t_azim = (float)atof(substr); + + /* Char 27: must be a space. */ + if (line[26] != ' ' ){ + sprintf(isf_error,"bad format, char 27: %s",line); + return 20; + } + + /* Chars 28-32: t plunge, must be float. */ + if (!partline(substr,line,27,5)){ + sprintf(isf_error,"missing t_pl: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad t_pl: %s",line); + return 20; + } + *t_pl = (float)atof(substr); + + /* Char 33: must be a space. */ + if (line[32] != ' ' ){ + sprintf(isf_error,"bad format, char 33: %s",line); + return 20; + } + + /* Chars 34-39: b value - float if there. */ + if (partline(substr,line,33,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad b_val: %s",line); + return 20; + } + *b_val = (float)atof(substr); + } + else { + *b_val = ISF_NULL; + } + + /* Char 40: must be a space. */ + if (line[39] != ' ' ){ + sprintf(isf_error,"bad format, char 40: %s",line); + return 20; + } + + /* Chars 41-46: b azim, must be float. */ + if (!partline(substr,line,40,6)){ + sprintf(isf_error,"missing b_azim: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad b_azim: %s",line); + return 20; + } + *b_azim = (float)atof(substr); + + /* Char 47: must be a space. */ + if (line[46] != ' ' ){ + sprintf(isf_error,"bad format, char 47: %s",line); + return 20; + } + + /* Chars 48-52: b plunge, must be float. */ + if (!partline(substr,line,47,5)){ + sprintf(isf_error,"missing b_pl: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad b_pl: %s",line); + return 20; + } + *b_pl = (float)atof(substr); + + + /* Char 53: must be a space. */ + if (line[52] != ' ' ){ + sprintf(isf_error,"bad format, char 53: %s",line); + return 20; + } + + /* Chars 54-59: p value - float if there. */ + if (partline(substr,line,53,6)){ + if (check_float(substr)){ + sprintf(isf_error,"bad p_val: %s",line); + return 20; + } + *p_val = (float)atof(substr); + } + else { + *p_val = ISF_NULL; + } + + /* Char 60: must be a space. */ + if (line[59] != ' ' ){ + sprintf(isf_error,"bad format, char 60: %s",line); + return 20; + } + + /* Chars 61-66: p azim, must be float. */ + if (!partline(substr,line,60,6)){ + sprintf(isf_error,"missing p_azim: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad p_azim: %s",line); + return 20; + } + *p_azim = (float)atof(substr); + + /* Char 67: must be a space. */ + if (line[66] != ' ' ){ + sprintf(isf_error,"bad format, char 67: %s",line); + return 20; + } + + /* Chars 68-72: p plunge, must be float */ + if (!partline(substr,line,67,5)){ + sprintf(isf_error,"missing p_pl: %s",line); + return 20; + } + if (check_float(substr)){ + sprintf(isf_error,"bad p_pl: %s",line); + return 20; + } + *p_pl = (float)atof(substr); + + /* Char 73: must be a space. */ + if (line[72] != ' ' ){ + sprintf(isf_error,"bad format, char 73: %s",line); + return 20; + } + + /* Chars 74-82: author, any characters allowed but must be there. */ + if (!partline(author,line,73,9)){ + sprintf(isf_error,"missing author: %s",line); + return 20; + } + if (check_whole(author)){ + sprintf(isf_error,"bad author: %s",line); + return 20; + } + + /* Check for extra stuff - not including close bracket. */ + if (partline(substr,line,82,0)){ + sprintf(isf_error,"extra characters at end: %s",line); + return 20; + } + return 0; +} + +/* Get a substring, removing leading and trailing white space. + Expects a string, an offset from the start of the string, and a maximum + length for the resulting substring. If this length is 0 it will take up + to the end of the input string. + + Need to allow for ')' to come after the required field at the end of a + comment. Discard ')' at end of a string as long as there's no '(' + before it anywhere. + + Returns the length of the resulting substring. +*/ +int partline(char *part, char *line, int offset, int numchars) { + int i,j,len; + int bracket=0; + + len = strlen(line); + if (len < offset) return 0; + + if (numchars == 0) numchars = len-offset; + + for (i=offset, j=0; i < offset+numchars; i++){ + if ( j == 0 && (line[i] == ' ' || line[i] == '\t')) continue; + if ( line[i] == '\0' || line[i] == '\n') break; + part[j++] = line[i]; + if ( line[i] == '(' ) bracket=1; + } + if (!bracket){ + while (--j != -1 && (part[j]==' ' || part[j]=='\t' || part[j]==')' )); + part[++j]='\0'; + } + else if (j){ + while (part[--j] == ' ' || part[j] == '\t'); + part[++j]='\0'; + } + return j; +} + +/* To check that a string has no spaces in it. + Returns 0 if there are no spaces or 1 if there is a space. +*/ +int check_whole(char *substr) { + int i; + for (i=0; i < strlen(substr); i++){ + if (substr[i] == ' ' || substr[i] == '\t') return 1; + } + return 0; +} + +/* To check that a string contains only sign/number characters and so + is suitable for atoi - atoi itself does no checking. + Returns 0 if OK, 1 if not. +*/ +int check_int(char *substr) { + int i; + for (i=0; i < strlen(substr); i++){ + if (isdigit(substr[i])) continue; + if (i==0) + if (substr[i] == '+' || substr[i] == '-') continue; + return 1; + } + return 0; +} + + +/* To check that a string contains only sign/number/point characters and so + is suitable for atof - atof itself does no checking. + Returns 0 if OK, 1 if not. +*/ +int check_float(char *substr) { + int i; + for (i=0; i < strlen(substr); i++){ + if (isdigit(substr[i]) || substr[i] == '.') continue; + if (i==0) + if (substr[i] == '+' || substr[i] == '-') continue; + return 1; + } + return 0; +} +