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