src/datain.c

   1 /* datain.c
   2  * Reads in survey files, dealing with special characters, keywords & data
   3  * Copyright (C) 1991-2003,2005,2009,2010,2011,2012,2013,2014,2015,2016 Olly Betts
   4  * Copyright (C) 2004 Simeon Warner
   5  *
   6  * This program is free software; you can redistribute it and/or modify
   7  * it under the terms of the GNU General Public License as published by
   8  * the Free Software Foundation; either version 2 of the License, or
   9  * (at your option) any later version.
  10  *
  11  * This program is distributed in the hope that it will be useful,
  12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14  * GNU General Public License for more details.
  15  *
  16  * You should have received a copy of the GNU General Public License
  17  * along with this program; if not, write to the Free Software
  18  * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
  19  */
  20
  21 #ifdef HAVE_CONFIG_H
  22 #include <config.h>
  23 #endif
  24
  25 #include <limits.h>
  26 #include <stdarg.h>
  27
  28 #include "debug.h"
  29 #include "cavern.h"
  30 #include "date.h"
  31 #include "filename.h"
  32 #include "message.h"
  33 #include "filelist.h"
  34 #include "netbits.h"
  35 #include "netskel.h"
  36 #include "readval.h"
  37 #include "datain.h"
  38 #include "commands.h"
  39 #include "out.h"
  40 #include "str.h"
  41 #include "thgeomag.h"
  42
  43 #define EPSILON (REAL_EPSILON * 1000)
  44
  45 #define var(I) (pcs->Var[(I)])
  46
  47 /* true if x is not-a-number value in Compass (999.0 or -999.0)    */
  48 /* Compass uses 999.0 but understands Karst data which used -999.0 */
  49 #define is_compass_NaN(x) ( fabs(fabs(x)-999.0) <  EPSILON )
  50
  51 int ch;
  52
  53 typedef enum {
  54     CTYPE_OMIT, CTYPE_READING, CTYPE_PLUMB, CTYPE_INFERPLUMB, CTYPE_HORIZ
  55 } clino_type;
  56
  57 /* Don't explicitly initialise as we can't set the jmp_buf - this has
  58  * static scope so will be initialised like this anyway */
  59 parse file /* = { NULL, NULL, 0, fFalse, NULL } */ ;
  60
  61 bool f_export_ok;
  62
  63 static real value[Fr - 1];
  64 #define VAL(N) value[(N)-1]
  65 static real variance[Fr - 1];
  66 #define VAR(N) variance[(N)-1]
  67 static long location[Fr - 1];
  68 #define LOC(N) location[(N)-1]
  69 static int location_width[Fr - 1];
  70 #define WID(N) location_width[(N)-1]
  71
  72 /* style functions */
  73 static void data_normal(void);
  74 static void data_cartesian(void);
  75 static void data_passage(void);
  76 static void data_nosurvey(void);
  77 static void data_ignore(void);
  78
  79 void
  80 get_pos(filepos *fp)
  81 {
  82    fp->ch = ch;
  83    fp->offset = ftell(file.fh);
  84    if (fp->offset == -1)
  85       fatalerror_in_file(file.filename, 0, /*Error reading file*/18);
  86 }
  87
  88 void
  89 set_pos(const filepos *fp)
  90 {
  91    ch = fp->ch;
  92    if (fseek(file.fh, fp->offset, SEEK_SET) == -1)
  93       fatalerror_in_file(file.filename, 0, /*Error reading file*/18);
  94 }
  95
  96 static void
  97 report_parent(parse * p) {
  98     if (p->parent)
  99         report_parent(p->parent);
 100     /* Force re-report of include tree for further errors in
 101      * parent files */
 102     p->reported_where = fFalse;
 103     /* TRANSLATORS: %s is replaced by the filename of the parent file, and %u
 104      * by the line number in that file.  Your translation should also contain
 105      * %s:%u so that automatic parsing of error messages to determine the file
 106      * and line number still works. */
 107     fprintf(STDERR, msg(/*In file included from %s:%u:\n*/5), p->filename, p->line);
 108 }
 109
 110 static void
 111 error_list_parent_files(void)
 112 {
 113    if (!file.reported_where && file.parent) {
 114       report_parent(file.parent);
 115       /* Suppress reporting of full include tree for further errors
 116        * in this file */
 117       file.reported_where = fTrue;
 118    }
 119 }
 120
 121 static void
 122 show_line(int col, int width)
 123 {
 124    /* Rewind to beginning of line. */
 125    long cur_pos = ftell(file.fh);
 126    int tabs = 0;
 127    if (cur_pos < 0 || fseek(file.fh, file.lpos, SEEK_SET) == -1)
 128       fatalerror_in_file(file.filename, 0, /*Error reading file*/18);
 129
 130    /* Read the whole line and write it out. */
 131    PUTC(' ', STDERR);
 132    while (!feof(file.fh)) {
 133       int c = GETC(file.fh);
 134       if (isEol(c)) break;
 135       if (c == '\t') ++tabs;
 136       PUTC(c, STDERR);
 137    }
 138    fputnl(STDERR);
 139
 140    /* If we have a location in the line for the error, indicate it. */
 141    if (col) {
 142       PUTC(' ', STDERR);
 143       if (tabs == 0) {
 144          while (--col) PUTC(' ', STDERR);
 145       } else {
 146          /* Copy tabs from line, replacing other characters with spaces - this
 147           * means that the caret should line up correctly. */
 148          if (fseek(file.fh, file.lpos, SEEK_SET) == -1)
 149             fatalerror_in_file(file.filename, 0, /*Error reading file*/18);
 150          while (--col) {
 151             int c = GETC(file.fh);
 152             if (c != '\t') c = ' ';
 153             PUTC(c, STDERR);
 154          }
 155       }
 156       PUTC('^', STDERR);
 157       while (width > 1) {
 158          PUTC('~', STDERR);
 159          --width;
 160       }
 161       fputnl(STDERR);
 162    }
 163
 164    /* Revert to where we were. */
 165    if (fseek(file.fh, cur_pos, SEEK_SET) == -1)
 166       fatalerror_in_file(file.filename, 0, /*Error reading file*/18);
 167 }
 168
 169 static int caret_width = 0;
 170
 171 static void
 172 compile_v_report_fpos(int severity, long fpos, int en, va_list ap)
 173 {
 174    int col = 0;
 175    error_list_parent_files();
 176    if (fpos >= file.lpos)
 177       col = fpos - file.lpos - caret_width;
 178    v_report(severity, file.filename, file.line, col, en, ap);
 179    if (file.fh) show_line(col, caret_width);
 180 }
 181
 182 static void
 183 compile_v_report(int diag_flags, int en, va_list ap)
 184 {
 185    int severity = (diag_flags & DIAG_SEVERITY_MASK);
 186    if (diag_flags & (DIAG_COL|DIAG_BUF)) {
 187       if (file.fh) {
 188          if (diag_flags & DIAG_BUF) caret_width = strlen(buffer);
 189          compile_v_report_fpos(severity, ftell(file.fh), en, ap);
 190          if (diag_flags & DIAG_BUF) caret_width = 0;
 191          if (diag_flags & DIAG_SKIP) skipline();
 192          return;
 193       }
 194    }
 195    error_list_parent_files();
 196    v_report(severity, file.filename, file.line, 0, en, ap);
 197    if (file.fh) {
 198       if (diag_flags & DIAG_BUF) {
 199          show_line(0, strlen(buffer));
 200       } else {
 201          show_line(0, caret_width);
 202       }
 203    }
 204    if (diag_flags & DIAG_SKIP) skipline();
 205 }
 206
 207 void
 208 compile_diagnostic(int diag_flags, int en, ...)
 209 {
 210    va_list ap;
 211    va_start(ap, en);
 212    if (diag_flags & (DIAG_TOKEN|DIAG_UINT|DIAG_DATE|DIAG_NUM)) {
 213       char *p = NULL;
 214       int len = 0;
 215       skipblanks();
 216       if (diag_flags & DIAG_TOKEN) {
 217          while (!isBlank(ch) && !isEol(ch)) {
 218             s_catchar(&p, &len, (char)ch);
 219             nextch();
 220          }
 221       } else if (diag_flags & DIAG_UINT) {
 222          while (isdigit(ch)) {
 223             s_catchar(&p, &len, (char)ch);
 224             nextch();
 225          }
 226       } else if (diag_flags & DIAG_DATE) {
 227          while (isdigit(ch) || ch == '.') {
 228             s_catchar(&p, &len, (char)ch);
 229             nextch();
 230          }
 231       } else {
 232          if (isMinus(ch) || isPlus(ch)) {
 233             s_catchar(&p, &len, (char)ch);
 234             nextch();
 235          }
 236          while (isdigit(ch)) {
 237             s_catchar(&p, &len, (char)ch);
 238             nextch();
 239          }
 240          if (isDecimal(ch)) {
 241             s_catchar(&p, &len, (char)ch);
 242             nextch();
 243          }
 244          while (isdigit(ch)) {
 245             s_catchar(&p, &len, (char)ch);
 246             nextch();
 247          }
 248       }
 249       if (p) {
 250          caret_width = strlen(p);
 251          osfree(p);
 252       }
 253       compile_v_report(diag_flags|DIAG_COL, en, ap);
 254       caret_width = 0;
 255    } else {
 256       compile_v_report(diag_flags, en, ap);
 257    }
 258    va_end(ap);
 259 }
 260
 261 static void
 262 compile_diagnostic_reading(int diag_flags, reading r, int en, ...)
 263 {
 264    va_list ap;
 265    int severity = (diag_flags & DIAG_SEVERITY_MASK);
 266    va_start(ap, en);
 267    caret_width = WID(r);
 268    compile_v_report_fpos(severity, LOC(r) + caret_width, en, ap);
 269    caret_width = 0;
 270    va_end(ap);
 271 }
 272
 273 static void
 274 compile_error_reading_skip(reading r, int en, ...)
 275 {
 276    va_list ap;
 277    va_start(ap, en);
 278    caret_width = WID(r);
 279    compile_v_report_fpos(1, LOC(r) + caret_width, en, ap);
 280    caret_width = 0;
 281    va_end(ap);
 282    skipline();
 283 }
 284
 285 void
 286 compile_diagnostic_at(int diag_flags, const char * filename, unsigned line, int en, ...)
 287 {
 288    va_list ap;
 289    int severity = (diag_flags & DIAG_SEVERITY_MASK);
 290    va_start(ap, en);
 291    v_report(severity, filename, line, 0, en, ap);
 292    va_end(ap);
 293 }
 294
 295 void
 296 compile_diagnostic_pfx(int diag_flags, const prefix * pfx, int en, ...)
 297 {
 298    va_list ap;
 299    int severity = (diag_flags & DIAG_SEVERITY_MASK);
 300    va_start(ap, en);
 301    v_report(severity, pfx->filename, pfx->line, 0, en, ap);
 302    va_end(ap);
 303 }
 304
 305 void
 306 compile_diagnostic_token_show(int diag_flags, int en)
 307 {
 308    char *p = NULL;
 309    int len = 0;
 310    skipblanks();
 311    while (!isBlank(ch) && !isEol(ch)) {
 312       s_catchar(&p, &len, (char)ch);
 313       nextch();
 314    }
 315    if (p) {
 316       caret_width = strlen(p);
 317       compile_diagnostic(diag_flags|DIAG_COL, en, p);
 318       caret_width = 0;
 319       osfree(p);
 320    } else {
 321       compile_diagnostic(DIAG_ERR|DIAG_COL, en, "");
 322    }
 323 }
 324
 325 static void
 326 compile_error_string(const char * s, int en, ...)
 327 {
 328    va_list ap;
 329    va_start(ap, en);
 330    caret_width = strlen(s);
 331    compile_v_report(DIAG_ERR|DIAG_COL, en, ap);
 332    va_end(ap);
 333    caret_width = 0;
 334 }
 335
 336 /* This function makes a note where to put output files */
 337 static void
 338 using_data_file(const char *fnm)
 339 {
 340    if (!fnm_output_base) {
 341       /* was: fnm_output_base = base_from_fnm(fnm); */
 342       fnm_output_base = baseleaf_from_fnm(fnm);
 343    } else if (fnm_output_base_is_dir) {
 344       /* --output pointed to directory so use the leaf basename in that dir */
 345       char *lf, *p;
 346       lf = baseleaf_from_fnm(fnm);
 347       p = use_path(fnm_output_base, lf);
 348       osfree(lf);
 349       osfree(fnm_output_base);
 350       fnm_output_base = p;
 351       fnm_output_base_is_dir = 0;
 352    }
 353 }
 354
 355 static void
 356 skipword(void)
 357 {
 358    while (!isBlank(ch) && !isEol(ch)) nextch();
 359 }
 360
 361 extern void
 362 skipblanks(void)
 363 {
 364    while (isBlank(ch)) nextch();
 365 }
 366
 367 extern void
 368 skipline(void)
 369 {
 370    while (!isEol(ch)) nextch();
 371 }
 372
 373 static void
 374 process_eol(void)
 375 {
 376    int eolchar;
 377
 378    skipblanks();
 379
 380    if (!isEol(ch)) {
 381       if (!isComm(ch))
 382          compile_diagnostic(DIAG_ERR|DIAG_COL, /*End of line not blank*/15);
 383       skipline();
 384    }
 385
 386    eolchar = ch;
 387    file.line++;
 388    /* skip any different eol characters so we get line counts correct on
 389     * DOS text files and similar, but don't count several adjacent blank
 390     * lines as one */
 391    while (ch != EOF) {
 392       nextch();
 393       if (ch == eolchar || !isEol(ch)) {
 394          break;
 395       }
 396       if (ch == '\n') eolchar = ch;
 397    }
 398    file.lpos = ftell(file.fh) - 1;
 399 }
 400
 401 static bool
 402 process_non_data_line(void)
 403 {
 404    skipblanks();
 405
 406    if (isData(ch)) return fFalse;
 407
 408    if (isKeywd(ch)) {
 409       nextch();
 410       handle_command();
 411    }
 412
 413    process_eol();
 414
 415    return fTrue;
 416 }
 417
 418 static void
 419 read_reading(reading r, bool f_optional)
 420 {
 421    int n_readings;
 422    q_quantity q;
 423    switch (r) {
 424       case Tape: q = Q_LENGTH; break;
 425       case BackTape: q = Q_BACKLENGTH; break;
 426       case Comp: q = Q_BEARING; break;
 427       case BackComp: q = Q_BACKBEARING; break;
 428       case Clino: q = Q_GRADIENT; break;
 429       case BackClino: q = Q_BACKGRADIENT; break;
 430       case FrDepth: case ToDepth: q = Q_DEPTH; break;
 431       case Dx: q = Q_DX; break;
 432       case Dy: q = Q_DY; break;
 433       case Dz: q = Q_DZ; break;
 434       case FrCount: case ToCount: q = Q_COUNT; break;
 435       case Left: q = Q_LEFT; break;
 436       case Right: q = Q_RIGHT; break;
 437       case Up: q = Q_UP; break;
 438       case Down: q = Q_DOWN; break;
 439       default:
 440         q = Q_NULL; /* Suppress compiler warning */;
 441         BUG("Unexpected case");
 442    }
 443    LOC(r) = ftell(file.fh);
 444    VAL(r) = read_numeric_multi(f_optional, &n_readings);
 445    WID(r) = ftell(file.fh) - LOC(r);
 446    VAR(r) = var(q);
 447    if (n_readings > 1) VAR(r) /= sqrt(n_readings);
 448 }
 449
 450 static void
 451 read_bearing_or_omit(reading r)
 452 {
 453    int n_readings;
 454    q_quantity q = Q_NULL;
 455    LOC(r) = ftell(file.fh);
 456    VAL(r) = read_numeric_multi_or_omit(&n_readings);
 457    WID(r) = ftell(file.fh) - LOC(r);
 458    switch (r) {
 459       case Comp: q = Q_BEARING; break;
 460       case BackComp: q = Q_BACKBEARING; break;
 461       default:
 462         q = Q_NULL; /* Suppress compiler warning */;
 463         BUG("Unexpected case");
 464    }
 465    VAR(r) = var(q);
 466    if (n_readings > 1) VAR(r) /= sqrt(n_readings);
 467 }
 468
 469 /* For reading Compass MAK files which have a freeform syntax */
 470 static void
 471 nextch_handling_eol(void)
 472 {
 473    nextch();
 474    while (ch != EOF && isEol(ch)) {
 475       process_eol();
 476    }
 477 }
 478
 479 #define LITLEN(S) (sizeof(S"") - 1)
 480 #define has_ext(F,L,E) ((L) > LITLEN(E) + 1 &&\
 481                         (F)[(L) - LITLEN(E) - 1] == FNM_SEP_EXT &&\
 482                         strcasecmp((F) + (L) - LITLEN(E), E) == 0)
 483 extern void
 484 data_file(const char *pth, const char *fnm)
 485 {
 486    int begin_lineno_store;
 487    parse file_store;
 488    volatile enum {FMT_SVX, FMT_DAT, FMT_MAK} fmt = FMT_SVX;
 489
 490    {
 491       char *filename;
 492       FILE *fh;
 493       size_t len;
 494
 495       if (!pth) {
 496          /* file specified on command line - don't do special translation */
 497          fh = fopenWithPthAndExt(pth, fnm, EXT_SVX_DATA, "rb", &filename);
 498       } else {
 499          fh = fopen_portable(pth, fnm, EXT_SVX_DATA, "rb", &filename);
 500       }
 501
 502       if (fh == NULL) {
 503          compile_error_string(fnm, /*Couldn’t open file “%s”*/24, fnm);
 504          return;
 505       }
 506
 507       len = strlen(filename);
 508       if (has_ext(filename, len, "dat")) {
 509          fmt = FMT_DAT;
 510       } else if (has_ext(filename, len, "mak")) {
 511          fmt = FMT_MAK;
 512       }
 513
 514       file_store = file;
 515       if (file.fh) file.parent = &file_store;
 516       file.fh = fh;
 517       file.filename = filename;
 518       file.line = 1;
 519       file.lpos = 0;
 520       file.reported_where = fFalse;
 521       nextch();
 522    }
 523
 524    using_data_file(file.filename);
 525
 526    begin_lineno_store = pcs->begin_lineno;
 527    pcs->begin_lineno = 0;
 528
 529    if (fmt == FMT_DAT) {
 530       short *t;
 531       int i;
 532       settings *pcsNew;
 533
 534       pcsNew = osnew(settings);
 535       *pcsNew = *pcs; /* copy contents */
 536       pcsNew->begin_lineno = 0;
 537       pcsNew->next = pcs;
 538       pcs = pcsNew;
 539       default_units(pcs);
 540       default_calib(pcs);
 541
 542       pcs->style = STYLE_NORMAL;
 543       pcs->units[Q_LENGTH] = METRES_PER_FOOT;
 544       t = ((short*)osmalloc(ossizeof(short) * 257)) + 1;
 545
 546       t[EOF] = SPECIAL_EOL;
 547       memset(t, 0, sizeof(short) * 33);
 548       for (i = 33; i < 127; i++) t[i] = SPECIAL_NAMES;
 549       t[127] = 0;
 550       for (i = 128; i < 256; i++) t[i] = SPECIAL_NAMES;
 551       t['\t'] |= SPECIAL_BLANK;
 552       t[' '] |= SPECIAL_BLANK;
 553       t['\032'] |= SPECIAL_EOL; /* Ctrl-Z, so olde DOS text files are handled ok */
 554       t['\n'] |= SPECIAL_EOL;
 555       t['\r'] |= SPECIAL_EOL;
 556       t['.'] |= SPECIAL_DECIMAL;
 557       t['-'] |= SPECIAL_MINUS;
 558       t['+'] |= SPECIAL_PLUS;
 559       pcs->Translate = t;
 560       pcs->Case = OFF;
 561       pcs->Truncate = INT_MAX;
 562       pcs->infer = BIT(INFER_EQUATES)|BIT(INFER_EXPORTS)|BIT(INFER_PLUMBS);
 563    } else if (fmt == FMT_MAK) {
 564       short *t;
 565       int i;
 566       settings *pcsNew;
 567
 568       pcsNew = osnew(settings);
 569       *pcsNew = *pcs; /* copy contents */
 570       pcsNew->begin_lineno = 0;
 571       pcsNew->next = pcs;
 572       pcs = pcsNew;
 573
 574       t = ((short*)osmalloc(ossizeof(short) * 257)) + 1;
 575
 576       t[EOF] = SPECIAL_EOL;
 577       memset(t, 0, sizeof(short) * 33);
 578       for (i = 33; i < 127; i++) t[i] = SPECIAL_NAMES;
 579       t[127] = 0;
 580       for (i = 128; i < 256; i++) t[i] = SPECIAL_NAMES;
 581       t['['] = t[','] = t[';'] = 0;
 582       t['\t'] |= SPECIAL_BLANK;
 583       t[' '] |= SPECIAL_BLANK;
 584       t['\032'] |= SPECIAL_EOL; /* Ctrl-Z, so olde DOS text files are handled ok */
 585       t['\n'] |= SPECIAL_EOL;
 586       t['\r'] |= SPECIAL_EOL;
 587       t['.'] |= SPECIAL_DECIMAL;
 588       t['-'] |= SPECIAL_MINUS;
 589       t['+'] |= SPECIAL_PLUS;
 590       pcs->Translate = t;
 591       pcs->Case = OFF;
 592       pcs->Truncate = INT_MAX;
 593    }
 594
 595 #ifdef HAVE_SETJMP_H
 596    /* errors in nested functions can longjmp here */
 597    if (setjmp(file.jbSkipLine)) {
 598       skipline();
 599       process_eol();
 600    }
 601 #endif
 602
 603    if (fmt == FMT_DAT) {
 604       while (!feof(file.fh) && !ferror(file.fh)) {
 605          static reading compass_order[] = {
 606             Fr, To, Tape, CompassDATComp, CompassDATClino,
 607             CompassDATLeft, CompassDATRight, CompassDATUp, CompassDATDown,
 608             CompassDATFlags, IgnoreAll
 609          };
 610          static reading compass_order_backsights[] = {
 611             Fr, To, Tape, CompassDATComp, CompassDATClino,
 612             CompassDATLeft, CompassDATRight, CompassDATUp, CompassDATDown,
 613             CompassDATBackComp, CompassDATBackClino,
 614             CompassDATFlags, IgnoreAll
 615          };
 616          /* <Cave name> */
 617          skipline();
 618          process_eol();
 619          /* SURVEY NAME: <Short name> */
 620          get_token();
 621          get_token();
 622          /* if (ch != ':') ... */
 623          nextch();
 624          get_token();
 625          skipline();
 626          process_eol();
 627          /* SURVEY DATE: 7 10 79  COMMENT:<Long name> */
 628          get_token();
 629          get_token();
 630          copy_on_write_meta(pcs);
 631          if (ch == ':') {
 632              int year, month, day;
 633
 634              nextch();
 635
 636              /* NB order is *month* *day* year */
 637              month = read_uint();
 638              day = read_uint();
 639              year = read_uint();
 640              /* Note: Larry says a 2 digit year is always 19XX */
 641              if (year < 100) year += 1900;
 642
 643              pcs->meta->days1 = pcs->meta->days2 = days_since_1900(year, month, day);
 644          } else {
 645              pcs->meta->days1 = pcs->meta->days2 = -1;
 646          }
 647          pcs->declination = HUGE_REAL;
 648          skipline();
 649          process_eol();
 650          /* SURVEY TEAM: */
 651          get_token();
 652          get_token();
 653          skipline();
 654          process_eol();
 655          /* <Survey team> */
 656          nextch();
 657          skipline();
 658          process_eol();
 659          /* DECLINATION: 1.00  FORMAT: DDDDLUDRADLN  CORRECTIONS: 2.00 3.00 4.00 */
 660          get_token();
 661          nextch(); /* : */
 662          skipblanks();
 663          pcs->z[Q_DECLINATION] = -read_numeric(fFalse);
 664          pcs->z[Q_DECLINATION] *= pcs->units[Q_DECLINATION];
 665          get_token();
 666          pcs->ordering = compass_order;
 667          if (strcmp(buffer, "FORMAT") == 0) {
 668             nextch(); /* : */
 669             get_token();
 670             if (strlen(buffer) >= 12 && buffer[11] == 'B') {
 671                /* We have backsights for compass and clino */
 672                pcs->ordering = compass_order_backsights;
 673             }
 674             get_token();
 675          }
 676          if (strcmp(buffer, "CORRECTIONS") == 0) {
 677             nextch(); /* : */
 678             pcs->z[Q_BEARING] = -rad(read_numeric(fFalse));
 679             pcs->z[Q_GRADIENT] = -rad(read_numeric(fFalse));
 680             pcs->z[Q_LENGTH] = -read_numeric(fFalse);
 681          } else {
 682             pcs->z[Q_BEARING] = 0;
 683             pcs->z[Q_GRADIENT] = 0;
 684             pcs->z[Q_LENGTH] = 0;
 685          }
 686          skipline();
 687          process_eol();
 688          /* BLANK LINE */
 689          skipline();
 690          process_eol();
 691          /* heading line */
 692          skipline();
 693          process_eol();
 694          /* BLANK LINE */
 695          skipline();
 696          process_eol();
 697          while (!feof(file.fh)) {
 698             if (ch == '\x0c') {
 699                nextch();
 700                process_eol();
 701                break;
 702             }
 703             data_normal();
 704          }
 705          clear_last_leg();
 706       }
 707       {
 708          settings *pcsParent = pcs->next;
 709          SVX_ASSERT(pcsParent);
 710          pcs->ordering = NULL;
 711          free_settings(pcs);
 712          pcs = pcsParent;
 713       }
 714    } else if (fmt == FMT_MAK) {
 715       nextch_handling_eol();
 716       while (!feof(file.fh) && !ferror(file.fh)) {
 717          if (ch == '#') {
 718             /* include a file */
 719             int ch_store;
 720             char *dat_pth = path_from_fnm(file.filename);
 721             char *dat_fnm = NULL;
 722             int dat_fnm_len;
 723             nextch_handling_eol();
 724             while (ch != ',' && ch != ';' && ch != EOF) {
 725                while (isEol(ch)) process_eol();
 726                s_catchar(&dat_fnm, &dat_fnm_len, (char)ch);
 727                nextch_handling_eol();
 728             }
 729             while (ch != ';' && ch != EOF) {
 730                prefix *name;
 731                nextch_handling_eol();
 732                name = read_prefix(PFX_STATION|PFX_OPT);
 733                if (name) {
 734                   skipblanks();
 735                   if (ch == '[') {
 736                      /* fixed pt */
 737                      node *stn;
 738                      real x, y, z;
 739                      name->sflags |= BIT(SFLAGS_FIXED);
 740                      nextch_handling_eol();
 741                      while (!isdigit(ch) && ch != '+' && ch != '-' &&
 742                             ch != '.' && ch != ']' && ch != EOF) {
 743                         nextch_handling_eol();
 744                      }
 745                      x = read_numeric(fFalse);
 746                      while (!isdigit(ch) && ch != '+' && ch != '-' &&
 747                             ch != '.' && ch != ']' && ch != EOF) {
 748                         nextch_handling_eol();
 749                      }
 750                      y = read_numeric(fFalse);
 751                      while (!isdigit(ch) && ch != '+' && ch != '-' &&
 752                             ch != '.' && ch != ']' && ch != EOF) {
 753                         nextch_handling_eol();
 754                      }
 755                      z = read_numeric(fFalse);
 756                      stn = StnFromPfx(name);
 757                      if (!fixed(stn)) {
 758                         POS(stn, 0) = x;
 759                         POS(stn, 1) = y;
 760                         POS(stn, 2) = z;
 761                         fix(stn);
 762                      } else {
 763                         if (x != POS(stn, 0) || y != POS(stn, 1) ||
 764                             z != POS(stn, 2)) {
 765                            compile_diagnostic(DIAG_ERR, /*Station already fixed or equated to a fixed point*/46);
 766                         } else {
 767                            compile_diagnostic(DIAG_WARN, /*Station already fixed at the same coordinates*/55);
 768                         }
 769                      }
 770                      while (ch != ']' && ch != EOF) nextch_handling_eol();
 771                      if (ch == ']') {
 772                         nextch_handling_eol();
 773                         skipblanks();
 774                      }
 775                   } else {
 776                      /* FIXME: link station - ignore for now */
 777                      /* FIXME: perhaps issue warning? */
 778                   }
 779                   while (ch != ',' && ch != ';' && ch != EOF)
 780                      nextch_handling_eol();
 781                }
 782             }
 783             if (dat_fnm) {
 784                ch_store = ch;
 785                data_file(dat_pth, dat_fnm);
 786                ch = ch_store;
 787                osfree(dat_fnm);
 788             }
 789          } else {
 790             /* FIXME: also check for % and $ later */
 791             nextch_handling_eol();
 792          }
 793       }
 794       {
 795          settings *pcsParent = pcs->next;
 796          SVX_ASSERT(pcsParent);
 797          free_settings(pcs);
 798          pcs = pcsParent;
 799       }
 800    } else {
 801       while (!feof(file.fh) && !ferror(file.fh)) {
 802          if (!process_non_data_line()) {
 803             f_export_ok = fFalse;
 804             switch (pcs->style) {
 805              case STYLE_NORMAL:
 806              case STYLE_DIVING:
 807              case STYLE_CYLPOLAR:
 808                data_normal();
 809                break;
 810              case STYLE_CARTESIAN:
 811                data_cartesian();
 812                break;
 813              case STYLE_PASSAGE:
 814                data_passage();
 815                break;
 816              case STYLE_NOSURVEY:
 817                data_nosurvey();
 818                break;
 819              case STYLE_IGNORE:
 820                data_ignore();
 821                break;
 822              default:
 823                BUG("bad style");
 824             }
 825          }
 826       }
 827       clear_last_leg();
 828    }
 829
 830    /* don't allow *BEGIN at the end of a file, then *EXPORT in the
 831     * including file */
 832    f_export_ok = fFalse;
 833
 834    if (pcs->begin_lineno) {
 835       error_in_file(file.filename, pcs->begin_lineno,
 836                     /*BEGIN with no matching END in this file*/23);
 837       /* Implicitly close any unclosed BEGINs from this file */
 838       do {
 839          settings *pcsParent = pcs->next;
 840          SVX_ASSERT(pcsParent);
 841          free_settings(pcs);
 842          pcs = pcsParent;
 843       } while (pcs->begin_lineno);
 844    }
 845
 846    pcs->begin_lineno = begin_lineno_store;
 847
 848    if (ferror(file.fh))
 849       fatalerror_in_file(file.filename, 0, /*Error reading file*/18);
 850
 851    (void)fclose(file.fh);
 852
 853    file = file_store;
 854
 855    /* don't free this - it may be pointed to by prefix.file */
 856    /* osfree(file.filename); */
 857 }
 858
 859 static real
 860 mod2pi(real a)
 861 {
 862    return a - floor(a / (2 * M_PI)) * (2 * M_PI);
 863 }
 864
 865 static real
 866 handle_plumb(clino_type *p_ctype)
 867 {
 868    typedef enum {
 869       CLINO_NULL=-1, CLINO_UP, CLINO_DOWN, CLINO_LEVEL
 870    } clino_tok;
 871    static sztok clino_tab[] = {
 872       {"D",     CLINO_DOWN},
 873       {"DOWN",  CLINO_DOWN},
 874       {"H",     CLINO_LEVEL},
 875       {"LEVEL", CLINO_LEVEL},
 876       {"U",     CLINO_UP},
 877       {"UP",    CLINO_UP},
 878       {NULL,    CLINO_NULL}
 879    };
 880    static real clinos[] = {(real)M_PI_2, (real)(-M_PI_2), (real)0.0};
 881    clino_tok tok;
 882
 883    skipblanks();
 884    if (isalpha(ch)) {
 885       filepos fp;
 886       get_pos(&fp);
 887       get_token();
 888       tok = match_tok(clino_tab, TABSIZE(clino_tab));
 889       if (tok != CLINO_NULL) {
 890          *p_ctype = (tok == CLINO_LEVEL ? CTYPE_HORIZ : CTYPE_PLUMB);
 891          return clinos[tok];
 892       }
 893       set_pos(&fp);
 894    } else if (isSign(ch)) {
 895       int chOld = ch;
 896       nextch();
 897       if (toupper(ch) == 'V') {
 898          nextch();
 899          *p_ctype = CTYPE_PLUMB;
 900          return (!isMinus(chOld) ? M_PI_2 : -M_PI_2);
 901       }
 902
 903       if (isOmit(chOld)) {
 904          *p_ctype = CTYPE_OMIT;
 905          /* no clino reading, so assume 0 with large sd */
 906          return (real)0.0;
 907       }
 908    } else if (isOmit(ch)) {
 909       /* OMIT char may not be a SIGN char too so we need to check here as
 910        * well as above... */
 911       nextch();
 912       *p_ctype = CTYPE_OMIT;
 913       /* no clino reading, so assume 0 with large sd */
 914       return (real)0.0;
 915    }
 916    return HUGE_REAL;
 917 }
 918
 919 static void
 920 warn_readings_differ(int msgno, real diff, int units)
 921 {
 922    char buf[64];
 923    char *p;
 924    diff /= get_units_factor(units);
 925    sprintf(buf, "%.2f", fabs(diff));
 926    for (p = buf; *p; ++p) {
 927       if (*p == '.') {
 928          char *z = p;
 929          while (*++p) {
 930             if (*p != '0') z = p + 1;
 931          }
 932          p = z;
 933          break;
 934       }
 935    }
 936    strcpy(p, get_units_string(units));
 937    compile_diagnostic(DIAG_WARN, msgno, buf);
 938 }
 939
 940 static bool
 941 handle_comp_units(void)
 942 {
 943    bool fNoComp = fTrue;
 944    if (VAL(Comp) != HUGE_REAL) {
 945       fNoComp = fFalse;
 946       VAL(Comp) *= pcs->units[Q_BEARING];
 947       if (VAL(Comp) < (real)0.0 || VAL(Comp) - M_PI * 2.0 > EPSILON) {
 948          /* TRANSLATORS: Suspicious means something like 410 degrees or -20
 949           * degrees */
 950          compile_diagnostic_reading(DIAG_WARN, Comp, /*Suspicious compass reading*/59);
 951          VAL(Comp) = mod2pi(VAL(Comp));
 952       }
 953    }
 954    if (VAL(BackComp) != HUGE_REAL) {
 955       fNoComp = fFalse;
 956       VAL(BackComp) *= pcs->units[Q_BACKBEARING];
 957       if (VAL(BackComp) < (real)0.0 || VAL(BackComp) - M_PI * 2.0 > EPSILON) {
 958          /* FIXME: different message for BackComp? */
 959          compile_diagnostic_reading(DIAG_WARN, BackComp, /*Suspicious compass reading*/59);
 960          VAL(BackComp) = mod2pi(VAL(BackComp));
 961       }
 962    }
 963    return fNoComp;
 964 }
 965
 966 static real
 967 handle_compass(real *p_var)
 968 {
 969    real compvar = VAR(Comp);
 970    real comp = VAL(Comp);
 971    real backcomp = VAL(BackComp);
 972    real declination;
 973    if (pcs->z[Q_DECLINATION] != HUGE_REAL) {
 974       declination = -pcs->z[Q_DECLINATION];
 975    } else if (pcs->declination != HUGE_REAL) {
 976       /* Cached value calculated for a previous compass reading taken on the
 977        * same date (by the 'else' just below).
 978        */
 979       declination = pcs->declination;
 980    } else {
 981       if (!pcs->meta || pcs->meta->days1 == -1) {
 982           compile_diagnostic(DIAG_WARN, /*No survey date specified - using 0 for magnetic declination*/304);
 983           declination = 0;
 984       } else {
 985           int avg_days = (pcs->meta->days1 + pcs->meta->days2) / 2;
 986           double dat = julian_date_from_days_since_1900(avg_days);
 987           /* thgeomag() takes (lat, lon, h, dat) - i.e. (y, x, z, date). */
 988           declination = thgeomag(pcs->dec_y, pcs->dec_x, pcs->dec_z, dat);
 989       }
 990       declination -= pcs->convergence;
 991       /* We cache the calculated declination as the calculation is relatively
 992        * expensive.  We also cache an "assumed 0" answer so that we only
 993        * warn once per such survey rather than for every line with a compass
 994        * reading. */
 995       pcs->declination = declination;
 996    }
 997    if (comp != HUGE_REAL) {
 998       comp = (comp - pcs->z[Q_BEARING]) * pcs->sc[Q_BEARING];
 999       comp += declination;
1000    }
1001    if (backcomp != HUGE_REAL) {
1002       backcomp = (backcomp - pcs->z[Q_BACKBEARING])
1003               * pcs->sc[Q_BACKBEARING];
1004       backcomp += declination;
1005       backcomp -= M_PI;
1006       if (comp != HUGE_REAL) {
1007          real diff = comp - backcomp;
1008          real adj = fabs(diff) > M_PI ? M_PI : 0;
1009          diff -= floor((diff + M_PI) / (2 * M_PI)) * 2 * M_PI;
1010          if (sqrd(diff / 3.0) > compvar + VAR(BackComp)) {
1011             /* fore and back readings differ by more than 3 sds */
1012             /* TRANSLATORS: %s is replaced by the amount the readings disagree
1013              * by, e.g. "2.5°" or "3ᵍ". */
1014             warn_readings_differ(/*COMPASS reading and BACKCOMPASS reading disagree by %s*/98,
1015                                  diff, get_angle_units(Q_BEARING));
1016          }
1017          comp = (comp / compvar + backcomp / VAR(BackComp));
1018          compvar = (compvar + VAR(BackComp)) / 4;
1019          comp *= compvar;
1020          comp += adj;
1021       } else {
1022          comp = backcomp;
1023          compvar = VAR(BackComp);
1024       }
1025    }
1026    *p_var = compvar;
1027    return comp;
1028 }
1029
1030 static real
1031 handle_clino(q_quantity q, reading r, real val, bool percent, clino_type *p_ctype)
1032 {
1033    bool range_0_180;
1034    real z;
1035    real diff_from_abs90;
1036    val *= pcs->units[q];
1037    /* percentage scale */
1038    if (percent) val = atan(val);
1039    /* We want to warn if there's a reading which it would be impossible
1040     * to have read from the instrument (e.g. on a -90 to 90 degree scale
1041     * you can't read "96" (it's probably a typo for "69").  However, the
1042     * gradient reading from a topofil is typically in the range 0 to 180,
1043     * with 90 being horizontal.
1044     *
1045     * Really we should allow the valid range to be specified, but for now
1046     * we infer it from the zero error - if this is within 45 degrees of
1047     * 90 then we assume the range is 0 to 180.
1048     */
1049    z = pcs->z[q];
1050    range_0_180 = (z > M_PI_4 && z < 3*M_PI_4);
1051    diff_from_abs90 = fabs(val) - M_PI_2;
1052    if (diff_from_abs90 > EPSILON) {
1053       if (!range_0_180) {
1054          int clino_units = get_angle_units(q);
1055          const char * units = get_units_string(clino_units);
1056          real right_angle = M_PI_2 / get_units_factor(clino_units);
1057          /* FIXME: different message for BackClino? */
1058          /* TRANSLATORS: %.f%s will be replaced with a right angle in the
1059           * units currently in use, e.g. "90°" or "100ᵍ".  And "absolute
1060           * value" means the reading ignoring the sign (so it might be
1061           * < -90° or > 90°. */
1062          compile_diagnostic_reading(DIAG_WARN, r, /*Clino reading over %.f%s (absolute value)*/51,
1063                                     right_angle, units);
1064       }
1065    } else if (TSTBIT(pcs->infer, INFER_PLUMBS) &&
1066               diff_from_abs90 >= -EPSILON) {
1067       *p_ctype = CTYPE_INFERPLUMB;
1068    }
1069    if (range_0_180 && *p_ctype != CTYPE_INFERPLUMB) {
1070       /* FIXME: Warning message not ideal... */
1071       if (val < 0.0 || val - M_PI > EPSILON) {
1072          int clino_units = get_angle_units(q);
1073          const char * units = get_units_string(clino_units);
1074          real right_angle = M_PI_2 / get_units_factor(clino_units);
1075          compile_diagnostic_reading(DIAG_WARN, r, /*Clino reading over %.f%s (absolute value)*/51,
1076                                     right_angle, units);
1077       }
1078    }
1079    return val;
1080 }
1081
1082 static int
1083 process_normal(prefix *fr, prefix *to, bool fToFirst,
1084                clino_type ctype, clino_type backctype)
1085 {
1086    real tape = VAL(Tape);
1087    real clin = VAL(Clino);
1088    real backclin = VAL(BackClino);
1089
1090    real dx, dy, dz;
1091    real vx, vy, vz;
1092 #ifndef NO_COVARIANCES
1093    real cxy, cyz, czx;
1094 #endif
1095
1096    bool fNoComp;
1097
1098    /* adjusted tape is negative -- probably the calibration is wrong */
1099    if (tape < (real)0.0) {
1100       /* TRANSLATE different message for topofil? */
1101       compile_diagnostic_reading(DIAG_WARN, Tape, /*Negative adjusted tape reading*/79);
1102    }
1103
1104    fNoComp = handle_comp_units();
1105
1106    if (ctype == CTYPE_READING) {
1107       clin = handle_clino(Q_GRADIENT, Clino, clin,
1108                           pcs->f_clino_percent, &ctype);
1109    }
1110
1111    if (backctype == CTYPE_READING) {
1112       backclin = handle_clino(Q_BACKGRADIENT, BackClino, backclin,
1113                               pcs->f_backclino_percent, &backctype);
1114    }
1115
1116    /* un-infer the plumb if the backsight was just a reading */
1117    if (ctype == CTYPE_INFERPLUMB && backctype == CTYPE_READING) {
1118        ctype = CTYPE_READING;
1119    }
1120
1121    if (ctype != CTYPE_OMIT && backctype != CTYPE_OMIT && ctype != backctype) {
1122       /* TRANSLATORS: In data with backsights, the user has tried to give a
1123        * plumb for the foresight and a clino reading for the backsight, or
1124        * something similar. */
1125       compile_error_reading_skip(Clino, /*CLINO and BACKCLINO readings must be of the same type*/84);
1126       return 0;
1127    }
1128
1129    if (ctype == CTYPE_PLUMB || ctype == CTYPE_INFERPLUMB ||
1130        backctype == CTYPE_PLUMB || backctype == CTYPE_INFERPLUMB) {
1131       /* plumbed */
1132       if (!fNoComp) {
1133          if (ctype == CTYPE_PLUMB ||
1134              (ctype == CTYPE_INFERPLUMB && VAL(Comp) != 0.0) ||
1135              backctype == CTYPE_PLUMB ||
1136              (backctype == CTYPE_INFERPLUMB && VAL(BackComp) != 0.0)) {
1137             /* FIXME: Different message for BackComp? */
1138             /* TRANSLATORS: A "plumbed leg" is one measured using a plumbline
1139              * (a weight on a string).  So the problem here is that the leg is
1140              * vertical, so a compass reading has no meaning! */
1141             compile_diagnostic(DIAG_WARN, /*Compass reading given on plumbed leg*/21);
1142          }
1143       }
1144
1145       dx = dy = (real)0.0;
1146       if (ctype != CTYPE_OMIT) {
1147          if (backctype != CTYPE_OMIT && (clin > 0) == (backclin > 0)) {
1148             /* TRANSLATORS: We've been told the foresight and backsight are
1149              * both "UP", or that they're both "DOWN". */
1150             compile_error_reading_skip(Clino, /*Plumbed CLINO and BACKCLINO readings can't be in the same direction*/92);
1151             return 0;
1152          }
1153          dz = (clin > (real)0.0) ? tape : -tape;
1154       } else {
1155          dz = (backclin < (real)0.0) ? tape : -tape;
1156       }
1157       vx = vy = var(Q_POS) / 3.0 + dz * dz * var(Q_PLUMB);
1158       vz = var(Q_POS) / 3.0 + VAR(Tape);
1159 #ifndef NO_COVARIANCES
1160       /* Correct values - no covariances in this case! */
1161       cxy = cyz = czx = (real)0.0;
1162 #endif
1163    } else {
1164       /* Each of ctype and backctype are either CTYPE_READING/CTYPE_HORIZ
1165        * or CTYPE_OMIT */
1166       /* clino */
1167       real L2, cosG, LcosG, cosG2, sinB, cosB, dx2, dy2, dz2, v, V;
1168       if (fNoComp) {
1169          /* TRANSLATORS: Here "legs" are survey legs, i.e. measurements between
1170           * survey stations. */
1171          compile_error_reading_skip(Comp, /*Compass reading may not be omitted except on plumbed legs*/14);
1172          return 0;
1173       }
1174       if (tape == (real)0.0) {
1175          dx = dy = dz = (real)0.0;
1176          vx = vy = vz = (real)(var(Q_POS) / 3.0); /* Position error only */
1177 #ifndef NO_COVARIANCES
1178          cxy = cyz = czx = (real)0.0;
1179 #endif
1180 #if DEBUG_DATAIN_1
1181          printf("Zero length leg: vx = %f, vy = %f, vz = %f\n", vx, vy, vz);
1182 #endif
1183       } else {
1184          real sinGcosG;
1185          /* take into account variance in LEVEL case */
1186          real var_clin = var(Q_LEVEL);
1187          real var_comp;
1188          real comp = handle_compass(&var_comp);
1189          /* ctype != CTYPE_READING is LEVEL case */
1190          if (ctype == CTYPE_READING) {
1191             clin = (clin - pcs->z[Q_GRADIENT]) * pcs->sc[Q_GRADIENT];
1192             var_clin = VAR(Clino);
1193          }
1194          if (backctype == CTYPE_READING) {
1195             backclin = (backclin - pcs->z[Q_BACKGRADIENT])
1196                * pcs->sc[Q_BACKGRADIENT];
1197             if (ctype == CTYPE_READING) {
1198                if (sqrd((clin + backclin) / 3.0) > var_clin + VAR(BackClino)) {
1199                   /* fore and back readings differ by more than 3 sds */
1200                   /* TRANSLATORS: %s is replaced by the amount the readings disagree
1201                    * by, e.g. "2.5°" or "3ᵍ". */
1202                   warn_readings_differ(/*CLINO reading and BACKCLINO reading disagree by %s*/99,
1203                                        clin + backclin, get_angle_units(Q_GRADIENT));
1204                }
1205                clin = (clin / var_clin - backclin / VAR(BackClino));
1206                var_clin = (var_clin + VAR(BackClino)) / 4;
1207                clin *= var_clin;
1208             } else {
1209                clin = -backclin;
1210                var_clin = VAR(BackClino);
1211             }
1212          }
1213
1214 #if DEBUG_DATAIN
1215          printf("    %4.2f %4.2f %4.2f\n", tape, comp, clin);
1216 #endif
1217          cosG = cos(clin);
1218          LcosG = tape * cosG;
1219          sinB = sin(comp);
1220          cosB = cos(comp);
1221 #if DEBUG_DATAIN_1
1222          printf("sinB = %f, cosG = %f, LcosG = %f\n", sinB, cosG, LcosG);
1223 #endif
1224          dx = LcosG * sinB;
1225          dy = LcosG * cosB;
1226          dz = tape * sin(clin);
1227 /*      printf("%.2f\n",clin); */
1228 #if DEBUG_DATAIN_1
1229          printf("dx = %f\ndy = %f\ndz = %f\n", dx, dy, dz);
1230 #endif
1231          dx2 = dx * dx;
1232          L2 = tape * tape;
1233          V = VAR(Tape) / L2;
1234          dy2 = dy * dy;
1235          cosG2 = cosG * cosG;
1236          sinGcosG = sin(clin) * cosG;
1237          dz2 = dz * dz;
1238          v = dz2 * var_clin;
1239 #ifdef NO_COVARIANCES
1240          vx = (var(Q_POS) / 3.0 + dx2 * V + dy2 * var_comp +
1241                (.5 + sinB * sinB * cosG2) * v);
1242          vy = (var(Q_POS) / 3.0 + dy2 * V + dx2 * var_comp +
1243                (.5 + cosB * cosB * cosG2) * v);
1244          if (ctype == CTYPE_OMIT && backctype == CTYPE_OMIT) {
1245             /* if no clino, assume sd=tape/sqrt(10) so 3sds = .95*tape */
1246             vz = var(Q_POS) / 3.0 + L2 * (real)0.1;
1247          } else {
1248             vz = var(Q_POS) / 3.0 + dz2 * V + L2 * cosG2 * var_clin;
1249          }
1250          /* for Surveyor87 errors: vx=vy=vz=var(Q_POS)/3.0; */
1251 #else
1252          vx = var(Q_POS) / 3.0 + dx2 * V + dy2 * var_comp +
1253             (sinB * sinB * v);
1254          vy = var(Q_POS) / 3.0 + dy2 * V + dx2 * var_comp +
1255             (cosB * cosB * v);
1256          if (ctype == CTYPE_OMIT && backctype == CTYPE_OMIT) {
1257             /* if no clino, assume sd=tape/sqrt(10) so 3sds = .95*tape */
1258             vz = var(Q_POS) / 3.0 + L2 * (real)0.1;
1259          } else {
1260             vz = var(Q_POS) / 3.0 + dz2 * V + L2 * cosG2 * var_clin;
1261          }
1262          /* usual covariance formulae are fine in no clino case since
1263           * dz = 0 so value of var_clin is ignored */
1264          cxy = sinB * cosB * (VAR(Tape) * cosG2 + var_clin * dz2)
1265                - var_comp * dx * dy;
1266          czx = VAR(Tape) * sinB * sinGcosG - var_clin * dx * dz;
1267          cyz = VAR(Tape) * cosB * sinGcosG - var_clin * dy * dz;
1268 #if 0
1269          printf("vx = %6.3f, vy = %6.3f, vz = %6.3f\n", vx, vy, vz);
1270          printf("cxy = %6.3f, cyz = %6.3f, czx = %6.3f\n", cxy, cyz, czx);
1271 #endif
1272 #endif
1273 #if DEBUG_DATAIN_1
1274          printf("In DATAIN.C, vx = %f, vy = %f, vz = %f\n", vx, vy, vz);
1275 #endif
1276       }
1277    }
1278 #if DEBUG_DATAIN_1
1279    printf("Just before addleg, vx = %f\n", vx);
1280 #endif
1281    /*printf("dx,dy,dz = %.2f %.2f %.2f\n\n", dx, dy, dz);*/
1282    addlegbyname(fr, to, fToFirst, dx, dy, dz, vx, vy, vz
1283 #ifndef NO_COVARIANCES
1284                 , cyz, czx, cxy
1285 #endif
1286                 );
1287    return 1;
1288 }
1289
1290 static int
1291 process_diving(prefix *fr, prefix *to, bool fToFirst, bool fDepthChange)
1292 {
1293    real tape = VAL(Tape);
1294
1295    real dx, dy, dz;
1296    real vx, vy, vz;
1297 #ifndef NO_COVARIANCES
1298    real cxy = 0, cyz = 0, czx = 0;
1299 #endif
1300
1301    handle_comp_units();
1302
1303    /* depth gauge readings increase upwards with default calibration */
1304    if (fDepthChange) {
1305       SVX_ASSERT(VAL(FrDepth) == 0.0);
1306       dz = VAL(ToDepth) * pcs->units[Q_DEPTH] - pcs->z[Q_DEPTH];
1307       dz *= pcs->sc[Q_DEPTH];
1308    } else {
1309       dz = VAL(ToDepth) - VAL(FrDepth);
1310       dz *= pcs->units[Q_DEPTH] * pcs->sc[Q_DEPTH];
1311    }
1312
1313    /* adjusted tape is negative -- probably the calibration is wrong */
1314    if (tape < (real)0.0) {
1315       compile_diagnostic(DIAG_WARN, /*Negative adjusted tape reading*/79);
1316    }
1317
1318    /* check if tape is less than depth change */
1319    if (tape < fabs(dz)) {
1320       /* FIXME: allow margin of error based on variances? */
1321       /* TRANSLATORS: This means that the data fed in said this.
1322        *
1323        * It could be a gross error (e.g. the decimal point is missing from the
1324        * depth gauge reading) or it could just be due to random error on a near
1325        * vertical leg */
1326       compile_diagnostic(DIAG_WARN, /*Tape reading is less than change in depth*/62);
1327    }
1328
1329    if (tape == (real)0.0 && dz == 0.0) {
1330       dx = dy = dz = (real)0.0;
1331       vx = vy = vz = (real)(var(Q_POS) / 3.0); /* Position error only */
1332    } else if (VAL(Comp) == HUGE_REAL &&
1333               VAL(BackComp) == HUGE_REAL) {
1334       /* plumb */
1335       dx = dy = (real)0.0;
1336       if (dz < 0) tape = -tape;
1337       /* FIXME: Should use FrDepth sometimes... */
1338       dz = (dz * VAR(Tape) + tape * 2 * VAR(ToDepth))
1339          / (VAR(Tape) * 2 * VAR(ToDepth));
1340       vx = vy = var(Q_POS) / 3.0 + dz * dz * var(Q_PLUMB);
1341       /* FIXME: Should use FrDepth sometimes... */
1342       vz = var(Q_POS) / 3.0 + VAR(Tape) * 2 * VAR(ToDepth)
1343                               / (VAR(Tape) + VAR(ToDepth));
1344    } else {
1345       real L2, sinB, cosB, dz2, D2;
1346       real var_comp;
1347       real comp = handle_compass(&var_comp);
1348       sinB = sin(comp);
1349       cosB = cos(comp);
1350       L2 = tape * tape;
1351       dz2 = dz * dz;
1352       D2 = L2 - dz2;
1353       if (D2 <= (real)0.0) {
1354          /* FIXME: Should use FrDepth sometimes... */
1355          real vsum = VAR(Tape) + 2 * VAR(ToDepth);
1356          dx = dy = (real)0.0;
1357          vx = vy = var(Q_POS) / 3.0;
1358          /* FIXME: Should use FrDepth sometimes... */
1359          vz = var(Q_POS) / 3.0 + VAR(Tape) * 2 * VAR(ToDepth) / vsum;
1360          if (dz > 0) {
1361             /* FIXME: Should use FrDepth sometimes... */
1362             dz = (dz * VAR(Tape) + tape * 2 * VAR(ToDepth)) / vsum;
1363          } else {
1364             dz = (dz * VAR(Tape) - tape * 2 * VAR(ToDepth)) / vsum;
1365          }
1366       } else {
1367          real D = sqrt(D2);
1368          /* FIXME: Should use FrDepth sometimes... */
1369          real F = VAR(Tape) * L2 + 2 * VAR(ToDepth) * D2;
1370          dx = D * sinB;
1371          dy = D * cosB;
1372
1373          vx = var(Q_POS) / 3.0 +
1374             sinB * sinB * F / D2 + var_comp * dy * dy;
1375          vy = var(Q_POS) / 3.0 +
1376             cosB * cosB * F / D2 + var_comp * dx * dx;
1377          /* FIXME: Should use FrDepth sometimes... */
1378          vz = var(Q_POS) / 3.0 + 2 * VAR(ToDepth);
1379
1380 #ifndef NO_COVARIANCES
1381          cxy = sinB * cosB * (F / D2 + var_comp * D2);
1382          /* FIXME: Should use FrDepth sometimes... */
1383          cyz = -2 * VAR(ToDepth) * dy / D;
1384          czx = -2 * VAR(ToDepth) * dx / D;
1385 #endif
1386       }
1387       /* FIXME: If there's a clino reading, check it against the depth reading,
1388        * and average.
1389        * if (VAL(Clino) != HUGE_REAL || VAL(BackClino) != HUGE_REAL) { ... }
1390        */
1391    }
1392    addlegbyname(fr, to, fToFirst, dx, dy, dz, vx, vy, vz
1393 #ifndef NO_COVARIANCES
1394                 , cxy, cyz, czx
1395 #endif
1396                 );
1397    return 1;
1398 }
1399
1400 static int
1401 process_cartesian(prefix *fr, prefix *to, bool fToFirst)
1402 {
1403    real dx = (VAL(Dx) * pcs->units[Q_DX] - pcs->z[Q_DX]) * pcs->sc[Q_DX];
1404    real dy = (VAL(Dy) * pcs->units[Q_DY] - pcs->z[Q_DY]) * pcs->sc[Q_DY];
1405    real dz = (VAL(Dz) * pcs->units[Q_DZ] - pcs->z[Q_DZ]) * pcs->sc[Q_DZ];
1406
1407    addlegbyname(fr, to, fToFirst, dx, dy, dz, VAR(Dx), VAR(Dy), VAR(Dz)
1408 #ifndef NO_COVARIANCES
1409                 , 0, 0, 0
1410 #endif
1411                 );
1412    return 1;
1413 }
1414
1415 static void
1416 data_cartesian(void)
1417 {
1418    prefix *fr = NULL, *to = NULL;
1419
1420    bool fMulti = fFalse;
1421
1422    reading first_stn = End;
1423
1424    reading *ordering;
1425
1426    again:
1427
1428    for (ordering = pcs->ordering ; ; ordering++) {
1429       skipblanks();
1430       switch (*ordering) {
1431        case Fr:
1432          fr = read_prefix(PFX_STATION|PFX_ALLOW_ROOT);
1433          if (first_stn == End) first_stn = Fr;
1434          break;
1435        case To:
1436          to = read_prefix(PFX_STATION|PFX_ALLOW_ROOT);
1437          if (first_stn == End) first_stn = To;
1438          break;
1439        case Station:
1440          fr = to;
1441          to = read_prefix(PFX_STATION);
1442          first_stn = To;
1443          break;
1444        case Dx: case Dy: case Dz:
1445          read_reading(*ordering, fFalse);
1446          break;
1447        case Ignore:
1448          skipword(); break;
1449        case IgnoreAllAndNewLine:
1450          skipline();
1451          /* fall through */
1452        case Newline:
1453          if (fr != NULL) {
1454             if (!process_cartesian(fr, to, first_stn == To))
1455                skipline();
1456          }
1457          fMulti = fTrue;
1458          while (1) {
1459             process_eol();
1460             skipblanks();
1461             if (isData(ch)) break;
1462             if (!isComm(ch)) {
1463                return;
1464             }
1465          }
1466          break;
1467        case IgnoreAll:
1468          skipline();
1469          /* fall through */
1470        case End:
1471          if (!fMulti) {
1472             process_cartesian(fr, to, first_stn == To);
1473             process_eol();
1474             return;
1475          }
1476          do {
1477             process_eol();
1478             skipblanks();
1479          } while (isComm(ch));
1480          goto again;
1481        default: BUG("Unknown reading in ordering");
1482       }
1483    }
1484 }
1485
1486 static int
1487 process_cylpolar(prefix *fr, prefix *to, bool fToFirst, bool fDepthChange)
1488 {
1489    real tape = VAL(Tape);
1490
1491    real dx, dy, dz;
1492    real vx, vy, vz;
1493 #ifndef NO_COVARIANCES
1494    real cxy = 0;
1495 #endif
1496
1497    handle_comp_units();
1498
1499    /* depth gauge readings increase upwards with default calibration */
1500    if (fDepthChange) {
1501       SVX_ASSERT(VAL(FrDepth) == 0.0);
1502       dz = VAL(ToDepth) * pcs->units[Q_DEPTH] - pcs->z[Q_DEPTH];
1503       dz *= pcs->sc[Q_DEPTH];
1504    } else {
1505       dz = VAL(ToDepth) - VAL(FrDepth);
1506       dz *= pcs->units[Q_DEPTH] * pcs->sc[Q_DEPTH];
1507    }
1508
1509    /* adjusted tape is negative -- probably the calibration is wrong */
1510    if (tape < (real)0.0) {
1511       compile_diagnostic(DIAG_WARN, /*Negative adjusted tape reading*/79);
1512    }
1513
1514    if (VAL(Comp) == HUGE_REAL && VAL(BackComp) == HUGE_REAL) {
1515       /* plumb */
1516       dx = dy = (real)0.0;
1517       vx = vy = var(Q_POS) / 3.0 + dz * dz * var(Q_PLUMB);
1518       /* FIXME: Should use FrDepth sometimes... */
1519       vz = var(Q_POS) / 3.0 + 2 * VAR(ToDepth);
1520    } else {
1521       real sinB, cosB;
1522       real var_comp;
1523       real comp = handle_compass(&var_comp);
1524       sinB = sin(comp);
1525       cosB = cos(comp);
1526
1527       dx = tape * sinB;
1528       dy = tape * cosB;
1529
1530       vx = var(Q_POS) / 3.0 +
1531          VAR(Tape) * sinB * sinB + var_comp * dy * dy;
1532       vy = var(Q_POS) / 3.0 +
1533          VAR(Tape) * cosB * cosB + var_comp * dx * dx;
1534       /* FIXME: Should use FrDepth sometimes... */
1535       vz = var(Q_POS) / 3.0 + 2 * VAR(ToDepth);
1536
1537 #ifndef NO_COVARIANCES
1538       cxy = (VAR(Tape) - var_comp * tape * tape) * sinB * cosB;
1539 #endif
1540    }
1541    addlegbyname(fr, to, fToFirst, dx, dy, dz, vx, vy, vz
1542 #ifndef NO_COVARIANCES
1543                 , cxy, 0, 0
1544 #endif
1545                 );
1546    return 1;
1547 }
1548
1549 /* Process tape/compass/clino, diving, and cylpolar styles of survey data
1550  * Also handles topofil (fromcount/tocount or count) in place of tape */
1551 static void
1552 data_normal(void)
1553 {
1554    prefix *fr = NULL, *to = NULL;
1555    reading first_stn = End;
1556
1557    bool fTopofil = fFalse, fMulti = fFalse;
1558    bool fRev;
1559    clino_type ctype, backctype;
1560    bool fDepthChange;
1561    unsigned long compass_dat_flags = 0;
1562
1563    reading *ordering;
1564
1565    VAL(Tape) = VAL(BackTape) = HUGE_REAL;
1566    VAL(Comp) = VAL(BackComp) = HUGE_REAL;
1567    VAL(FrCount) = VAL(ToCount) = 0;
1568    VAL(FrDepth) = VAL(ToDepth) = 0;
1569    VAL(Left) = VAL(Right) = VAL(Up) = VAL(Down) = HUGE_REAL;
1570
1571    fRev = fFalse;
1572    ctype = backctype = CTYPE_OMIT;
1573    fDepthChange = fFalse;
1574
1575    /* ordering may omit clino reading, so set up default here */
1576    /* this is also used if clino reading is the omit character */
1577    VAL(Clino) = VAL(BackClino) = 0;
1578
1579    again:
1580
1581    /* We clear these flags in the normal course of events, but if there's an
1582     * error in a reading, we might not, so make sure it has been cleared here.
1583     */
1584    pcs->flags &= ~(BIT(FLAGS_ANON_ONE_END) | BIT(FLAGS_IMPLICIT_SPLAY));
1585    for (ordering = pcs->ordering; ; ordering++) {
1586       skipblanks();
1587       switch (*ordering) {
1588        case Fr:
1589           fr = read_prefix(PFX_STATION|PFX_ALLOW_ROOT|PFX_ANON);
1590           if (first_stn == End) first_stn = Fr;
1591           break;
1592        case To:
1593           to = read_prefix(PFX_STATION|PFX_ALLOW_ROOT|PFX_ANON);
1594           if (first_stn == End) first_stn = To;
1595           break;
1596        case Station:
1597           fr = to;
1598           to = read_prefix(PFX_STATION);
1599           first_stn = To;
1600           break;
1601        case Dir: {
1602           typedef enum {
1603              DIR_NULL=-1, DIR_FORE, DIR_BACK
1604           } dir_tok;
1605           static sztok dir_tab[] = {
1606              {"B",     DIR_BACK},
1607              {"F",     DIR_FORE},
1608           };
1609           dir_tok tok;
1610           get_token();
1611           tok = match_tok(dir_tab, TABSIZE(dir_tab));
1612           switch (tok) {
1613            case DIR_FORE:
1614              break;
1615            case DIR_BACK:
1616              fRev = fTrue;
1617              break;
1618            default:
1619              compile_diagnostic(DIAG_ERR|DIAG_BUF|DIAG_SKIP, /*Found “%s”, expecting “F” or “B”*/131, buffer);
1620              process_eol();
1621              return;
1622           }
1623           break;
1624        }
1625        case Tape: case BackTape: {
1626           reading r = *ordering;
1627           read_reading(r, fTrue);
1628           if (VAL(r) == HUGE_REAL) {
1629              if (!isOmit(ch)) {
1630                 compile_diagnostic_token_show(DIAG_ERR, /*Expecting numeric field, found “%s”*/9);
1631                 /* Avoid also warning about omitted tape reading. */
1632                 VAL(r) = 0;
1633              } else {
1634                 nextch();
1635              }
1636           } else if (VAL(r) < (real)0.0) {
1637              compile_diagnostic_reading(DIAG_WARN, r, /*Negative tape reading*/60);
1638           }
1639           break;
1640        }
1641        case Count:
1642           VAL(FrCount) = VAL(ToCount);
1643           LOC(FrCount) = LOC(ToCount);
1644           WID(FrCount) = WID(ToCount);
1645           read_reading(ToCount, fFalse);
1646           fTopofil = fTrue;
1647           break;
1648        case FrCount:
1649           read_reading(FrCount, fFalse);
1650           break;
1651        case ToCount:
1652           read_reading(ToCount, fFalse);
1653           fTopofil = fTrue;
1654           break;
1655        case Comp: case BackComp:
1656           read_bearing_or_omit(*ordering);
1657           break;
1658        case Clino: case BackClino: {
1659           reading r = *ordering;
1660           clino_type * p_ctype = (r == Clino ? &ctype : &backctype);
1661           read_reading(r, fTrue);
1662           if (VAL(r) == HUGE_REAL) {
1663              VAL(r) = handle_plumb(p_ctype);
1664              if (VAL(r) != HUGE_REAL) break;
1665              compile_diagnostic_token_show(DIAG_ERR, /*Expecting numeric field, found “%s”*/9);
1666              skipline();
1667              process_eol();
1668              return;
1669           }
1670           *p_ctype = CTYPE_READING;
1671           break;
1672        }
1673        case FrDepth: case ToDepth:
1674           read_reading(*ordering, fFalse);
1675           break;
1676        case Depth:
1677           VAL(FrDepth) = VAL(ToDepth);
1678           LOC(FrDepth) = LOC(ToDepth);
1679           WID(FrDepth) = WID(ToDepth);
1680           read_reading(ToDepth, fFalse);
1681           break;
1682        case DepthChange:
1683           fDepthChange = fTrue;
1684           VAL(FrDepth) = 0;
1685           read_reading(ToDepth, fFalse);
1686           break;
1687        case CompassDATComp:
1688           read_bearing_or_omit(Comp);
1689           if (is_compass_NaN(VAL(Comp))) VAL(Comp) = HUGE_REAL;
1690           break;
1691        case CompassDATBackComp:
1692           read_bearing_or_omit(BackComp);
1693           if (is_compass_NaN(VAL(BackComp))) VAL(BackComp) = HUGE_REAL;
1694           break;
1695        case CompassDATClino: case CompassDATBackClino: {
1696           reading r;
1697           clino_type * p_ctype;
1698           if (*ordering == CompassDATClino) {
1699              r = Clino;
1700              p_ctype = &ctype;
1701           } else {
1702              r = BackClino;
1703              p_ctype = &backctype;
1704           }
1705           read_reading(r, fFalse);
1706           if (is_compass_NaN(VAL(r))) {
1707              VAL(r) = HUGE_REAL;
1708              *p_ctype = CTYPE_OMIT;
1709           } else {
1710              *p_ctype = CTYPE_READING;
1711           }
1712           break;
1713        }
1714        case CompassDATLeft: case CompassDATRight:
1715        case CompassDATUp: case CompassDATDown: {
1716           /* FIXME: need to actually make use of these entries! */
1717           reading actual = Left + (*ordering - CompassDATLeft);
1718           read_reading(actual, fFalse);
1719           if (VAL(actual) < 0) VAL(actual) = HUGE_REAL;
1720           break;
1721        }
1722        case CompassDATFlags:
1723           if (ch == '#') {
1724              filepos fp;
1725              get_pos(&fp);
1726              nextch();
1727              if (ch == '|') {
1728                 nextch();
1729                 while (ch >= 'A' && ch <= 'Z') {
1730                    compass_dat_flags |= BIT(ch - 'A');
1731                    /* We currently understand:
1732                     *   L (exclude from length)
1733                     *   X (exclude data)
1734                     * FIXME: but should also handle at least some of:
1735                     *   C (no adjustment) (set all (co)variances to 0?)
1736                     *   P (no plot) (new flag in 3d for "hidden by default"?)
1737                     */
1738                    nextch();
1739                 }
1740                 if (ch == '#') {
1741                    nextch();
1742                 } else {
1743                    compass_dat_flags = 0;
1744                    set_pos(&fp);
1745                 }
1746              } else {
1747                 set_pos(&fp);
1748              }
1749           }
1750           break;
1751        case Ignore:
1752           skipword(); break;
1753        case IgnoreAllAndNewLine:
1754           skipline();
1755           /* fall through */
1756        case Newline:
1757           if (fr != NULL) {
1758              int r;
1759              int save_flags;
1760              int implicit_splay;
1761              if (fTopofil) {
1762                 VAL(Tape) = VAL(ToCount) - VAL(FrCount);
1763                 LOC(Tape) = LOC(ToCount);
1764                 WID(Tape) = WID(ToCount);
1765              }
1766              /* Note: frdepth == todepth test works regardless of fDepthChange
1767               * (frdepth always zero, todepth is change of depth) and also
1768               * works for STYLE_NORMAL (both remain 0) */
1769              if (TSTBIT(pcs->infer, INFER_EQUATES) &&
1770                  (VAL(Tape) == (real)0.0 || VAL(Tape) == HUGE_REAL) &&
1771                  (VAL(BackTape) == (real)0.0 || VAL(BackTape) == HUGE_REAL) &&
1772                  VAL(FrDepth) == VAL(ToDepth)) {
1773                 process_equate(fr, to);
1774                 goto inferred_equate;
1775              }
1776              if (fRev) {
1777                 prefix *t = fr;
1778                 fr = to;
1779                 to = t;
1780              }
1781              if (fTopofil) {
1782                 VAL(Tape) *= pcs->units[Q_COUNT] * pcs->sc[Q_COUNT];
1783              } else if (VAL(Tape) != HUGE_REAL) {
1784                 VAL(Tape) *= pcs->units[Q_LENGTH];
1785                 VAL(Tape) -= pcs->z[Q_LENGTH];
1786                 VAL(Tape) *= pcs->sc[Q_LENGTH];
1787              }
1788              if (VAL(BackTape) != HUGE_REAL) {
1789                 VAL(BackTape) *= pcs->units[Q_BACKLENGTH];
1790                 VAL(BackTape) -= pcs->z[Q_BACKLENGTH];
1791                 VAL(BackTape) *= pcs->sc[Q_BACKLENGTH];
1792                 if (VAL(Tape) != HUGE_REAL) {
1793                    real diff = VAL(Tape) - VAL(BackTape);
1794                    if (sqrd(diff / 3.0) > VAR(Tape) + VAR(BackTape)) {
1795                       /* fore and back readings differ by more than 3 sds */
1796                       /* TRANSLATORS: %s is replaced by the amount the readings disagree
1797                        * by, e.g. "0.12m" or "0.2ft". */
1798                       warn_readings_differ(/*TAPE reading and BACKTAPE reading disagree by %s*/97,
1799                                            diff, get_length_units(Q_LENGTH));
1800                    }
1801                    VAL(Tape) = VAL(Tape) / VAR(Tape) + VAL(BackTape) / VAR(BackTape);
1802                    VAR(Tape) = (VAR(Tape) + VAR(BackTape)) / 4;
1803                    VAL(Tape) *= VAR(Tape);
1804                 } else {
1805                    VAL(Tape) = VAL(BackTape);
1806                    VAR(Tape) = VAR(BackTape);
1807                 }
1808              } else if (VAL(Tape) == HUGE_REAL) {
1809                 compile_diagnostic_reading(DIAG_ERR, Tape, /*Tape reading may not be omitted*/94);
1810                 goto inferred_equate;
1811              }
1812              implicit_splay = TSTBIT(pcs->flags, FLAGS_IMPLICIT_SPLAY);
1813              pcs->flags &= ~(BIT(FLAGS_ANON_ONE_END) | BIT(FLAGS_IMPLICIT_SPLAY));
1814              save_flags = pcs->flags;
1815              if (implicit_splay) {
1816                 pcs->flags |= BIT(FLAGS_SPLAY);
1817              }
1818              switch (pcs->style) {
1819               case STYLE_NORMAL:
1820                 r = process_normal(fr, to, (first_stn == To) ^ fRev,
1821                                    ctype, backctype);
1822                 break;
1823               case STYLE_DIVING:
1824                 /* FIXME: Handle any clino readings */
1825                 r = process_diving(fr, to, (first_stn == To) ^ fRev,
1826                                    fDepthChange);
1827                 break;
1828               case STYLE_CYLPOLAR:
1829                 r = process_cylpolar(fr, to, (first_stn == To) ^ fRev,
1830                                      fDepthChange);
1831                 break;
1832               default:
1833                 r = 0; /* avoid warning */
1834                 BUG("bad style");
1835              }
1836              pcs->flags = save_flags;
1837              if (!r) skipline();
1838
1839              /* Swap fr and to back to how they were for next line */
1840              if (fRev) {
1841                 prefix *t = fr;
1842                 fr = to;
1843                 to = t;
1844              }
1845           }
1846
1847           fRev = fFalse;
1848           ctype = backctype = CTYPE_OMIT;
1849           fDepthChange = fFalse;
1850
1851           /* ordering may omit clino reading, so set up default here */
1852           /* this is also used if clino reading is the omit character */
1853           VAL(Clino) = VAL(BackClino) = 0;
1854           LOC(Clino) = LOC(BackClino) = -1;
1855           WID(Clino) = WID(BackClino) = 0;
1856
1857           inferred_equate:
1858
1859           fMulti = fTrue;
1860           while (1) {
1861               process_eol();
1862               skipblanks();
1863               if (isData(ch)) break;
1864               if (!isComm(ch)) {
1865                  return;
1866               }
1867           }
1868           break;
1869        case IgnoreAll:
1870           skipline();
1871           /* fall through */
1872        case End:
1873           if (!fMulti) {
1874              int save_flags;
1875              int implicit_splay;
1876              /* Compass ignore flag is 'X' */
1877              if ((compass_dat_flags & BIT('X' - 'A'))) {
1878                 process_eol();
1879                 return;
1880              }
1881              if (fRev) {
1882                 prefix *t = fr;
1883                 fr = to;
1884                 to = t;
1885              }
1886              if (fTopofil) {
1887                 VAL(Tape) = VAL(ToCount) - VAL(FrCount);
1888                 LOC(Tape) = LOC(ToCount);
1889                 WID(Tape) = WID(ToCount);
1890              }
1891              /* Note: frdepth == todepth test works regardless of fDepthChange
1892               * (frdepth always zero, todepth is change of depth) and also
1893               * works for STYLE_NORMAL (both remain 0) */
1894              if (TSTBIT(pcs->infer, INFER_EQUATES) &&
1895                  (VAL(Tape) == (real)0.0 || VAL(Tape) == HUGE_REAL) &&
1896                  (VAL(BackTape) == (real)0.0 || VAL(BackTape) == HUGE_REAL) &&
1897                  VAL(FrDepth) == VAL(ToDepth)) {
1898                 process_equate(fr, to);
1899                 process_eol();
1900                 return;
1901              }
1902              if (fTopofil) {
1903                 VAL(Tape) *= pcs->units[Q_COUNT] * pcs->sc[Q_COUNT];
1904              } else if (VAL(Tape) != HUGE_REAL) {
1905                 VAL(Tape) *= pcs->units[Q_LENGTH];
1906                 VAL(Tape) -= pcs->z[Q_LENGTH];
1907                 VAL(Tape) *= pcs->sc[Q_LENGTH];
1908              }
1909              if (VAL(BackTape) != HUGE_REAL) {
1910                 VAL(BackTape) *= pcs->units[Q_BACKLENGTH];
1911                 VAL(BackTape) -= pcs->z[Q_BACKLENGTH];
1912                 VAL(BackTape) *= pcs->sc[Q_BACKLENGTH];
1913                 if (VAL(Tape) != HUGE_REAL) {
1914                    real diff = VAL(Tape) - VAL(BackTape);
1915                    if (sqrd(diff / 3.0) > VAR(Tape) + VAR(BackTape)) {
1916                       /* fore and back readings differ by more than 3 sds */
1917                       /* TRANSLATORS: %s is replaced by the amount the readings disagree
1918                        * by, e.g. "0.12m" or "0.2ft". */
1919                       warn_readings_differ(/*TAPE reading and BACKTAPE reading disagree by %s*/97,
1920                                            diff, get_length_units(Q_LENGTH));
1921                    }
1922                    VAL(Tape) = VAL(Tape) / VAR(Tape) + VAL(BackTape) / VAR(BackTape);
1923                    VAR(Tape) = (VAR(Tape) + VAR(BackTape)) / 4;
1924                    VAL(Tape) *= VAR(Tape);
1925                 } else {
1926                    VAL(Tape) = VAL(BackTape);
1927                    VAR(Tape) = VAR(BackTape);
1928                 }
1929              } else if (VAL(Tape) == HUGE_REAL) {
1930                 compile_diagnostic_reading(DIAG_ERR, Tape, /*Tape reading may not be omitted*/94);
1931                 process_eol();
1932                 return;
1933              }
1934              implicit_splay = TSTBIT(pcs->flags, FLAGS_IMPLICIT_SPLAY);
1935              pcs->flags &= ~(BIT(FLAGS_ANON_ONE_END) | BIT(FLAGS_IMPLICIT_SPLAY));
1936              save_flags = pcs->flags;
1937              if (implicit_splay) {
1938                 pcs->flags |= BIT(FLAGS_SPLAY);
1939              }
1940              if ((compass_dat_flags & BIT('L' - 'A'))) {
1941                 /* 'L' means "exclude from length" - map this to Survex's
1942                  * FLAGS_DUPLICATE. */
1943                 pcs->flags |= BIT(FLAGS_DUPLICATE);
1944              }
1945              switch (pcs->style) {
1946               case STYLE_NORMAL:
1947                 process_normal(fr, to, (first_stn == To) ^ fRev,
1948                                ctype, backctype);
1949                 break;
1950               case STYLE_DIVING:
1951                 /* FIXME: Handle any clino readings */
1952                 process_diving(fr, to, (first_stn == To) ^ fRev,
1953                                fDepthChange);
1954                 break;
1955               case STYLE_CYLPOLAR:
1956                 process_cylpolar(fr, to, (first_stn == To) ^ fRev,
1957                                  fDepthChange);
1958                 break;
1959               default:
1960                 BUG("bad style");
1961              }
1962              pcs->flags = save_flags;
1963
1964              process_eol();
1965              return;
1966           }
1967           do {
1968              process_eol();
1969              skipblanks();
1970           } while (isComm(ch));
1971           goto again;
1972        default:
1973           BUG("Unknown reading in ordering");
1974       }
1975    }
1976 }
1977
1978 static int
1979 process_lrud(prefix *stn)
1980 {
1981    SVX_ASSERT(next_lrud);
1982    lrud * xsect = osnew(lrud);
1983    xsect->stn = stn;
1984    xsect->l = (VAL(Left) * pcs->units[Q_LEFT] - pcs->z[Q_LEFT]) * pcs->sc[Q_LEFT];
1985    xsect->r = (VAL(Right) * pcs->units[Q_RIGHT] - pcs->z[Q_RIGHT]) * pcs->sc[Q_RIGHT];
1986    xsect->u = (VAL(Up) * pcs->units[Q_UP] - pcs->z[Q_UP]) * pcs->sc[Q_UP];
1987    xsect->d = (VAL(Down) * pcs->units[Q_DOWN] - pcs->z[Q_DOWN]) * pcs->sc[Q_DOWN];
1988    xsect->meta = pcs->meta;
1989    if (pcs->meta) ++pcs->meta->ref_count;
1990    xsect->next = NULL;
1991    *next_lrud = xsect;
1992    next_lrud = &(xsect->next);
1993
1994    return 1;
1995 }
1996
1997 static void
1998 data_passage(void)
1999 {
2000    prefix *stn = NULL;
2001    reading *ordering;
2002
2003    for (ordering = pcs->ordering ; ; ordering++) {
2004       skipblanks();
2005       switch (*ordering) {
2006        case Station:
2007          stn = read_prefix(PFX_STATION);
2008          break;
2009        case Left: case Right: case Up: case Down: {
2010          reading r = *ordering;
2011          read_reading(r, fTrue);
2012          if (VAL(r) == HUGE_REAL) {
2013             if (!isOmit(ch)) {
2014                compile_diagnostic_token_show(DIAG_ERR, /*Expecting numeric field, found “%s”*/9);
2015             } else {
2016                nextch();
2017             }
2018             VAL(r) = -1;
2019          }
2020          break;
2021        }
2022        case Ignore:
2023          skipword(); break;
2024        case IgnoreAll:
2025          skipline();
2026          /* fall through */
2027        case End: {
2028          process_lrud(stn);
2029          process_eol();
2030          return;
2031        }
2032        default: BUG("Unknown reading in ordering");
2033       }
2034    }
2035 }
2036
2037 static int
2038 process_nosurvey(prefix *fr, prefix *to, bool fToFirst)
2039 {
2040    nosurveylink *link;
2041
2042    /* Suppress "unused fixed point" warnings for these stations */
2043    fr->sflags |= BIT(SFLAGS_USED);
2044    to->sflags |= BIT(SFLAGS_USED);
2045
2046    /* add to linked list which is dealt with after network is solved */
2047    link = osnew(nosurveylink);
2048    if (fToFirst) {
2049       link->to = StnFromPfx(to);
2050       link->fr = StnFromPfx(fr);
2051    } else {
2052       link->fr = StnFromPfx(fr);
2053       link->to = StnFromPfx(to);
2054    }
2055    link->flags = pcs->flags | (STYLE_NOSURVEY << FLAGS_STYLE_BIT0);
2056    link->meta = pcs->meta;
2057    if (pcs->meta) ++pcs->meta->ref_count;
2058    link->next = nosurveyhead;
2059    nosurveyhead = link;
2060    return 1;
2061 }
2062
2063 static void
2064 data_nosurvey(void)
2065 {
2066    prefix *fr = NULL, *to = NULL;
2067
2068    bool fMulti = fFalse;
2069
2070    reading first_stn = End;
2071
2072    reading *ordering;
2073
2074    again:
2075
2076    for (ordering = pcs->ordering ; ; ordering++) {
2077       skipblanks();
2078       switch (*ordering) {
2079        case Fr:
2080           fr = read_prefix(PFX_STATION|PFX_ALLOW_ROOT);
2081           if (first_stn == End) first_stn = Fr;
2082           break;
2083        case To:
2084           to = read_prefix(PFX_STATION|PFX_ALLOW_ROOT);
2085           if (first_stn == End) first_stn = To;
2086           break;
2087        case Station:
2088           fr = to;
2089           to = read_prefix(PFX_STATION);
2090           first_stn = To;
2091           break;
2092        case Ignore:
2093          skipword(); break;
2094        case IgnoreAllAndNewLine:
2095          skipline();
2096          /* fall through */
2097        case Newline:
2098          if (fr != NULL) {
2099             if (!process_nosurvey(fr, to, first_stn == To))
2100                skipline();
2101          }
2102          if (ordering[1] == End) {
2103             do {
2104                process_eol();
2105                skipblanks();
2106             } while (isComm(ch));
2107             if (!isData(ch)) {
2108                return;
2109             }
2110             goto again;
2111          }
2112          fMulti = fTrue;
2113          while (1) {
2114             process_eol();
2115             skipblanks();
2116             if (isData(ch)) break;
2117             if (!isComm(ch)) {
2118                return;
2119             }
2120          }
2121          break;
2122        case IgnoreAll:
2123          skipline();
2124          /* fall through */
2125        case End:
2126          if (!fMulti) {
2127             (void)process_nosurvey(fr, to, first_stn == To);
2128             process_eol();
2129             return;
2130          }
2131          do {
2132             process_eol();
2133             skipblanks();
2134          } while (isComm(ch));
2135          goto again;
2136        default: BUG("Unknown reading in ordering");
2137       }
2138    }
2139 }
2140
2141 /* totally ignore a line of survey data */
2142 static void
2143 data_ignore(void)
2144 {
2145    skipline();
2146    process_eol();
2147 }