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