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