source: git/src/datain.c @ 2a0f70e

RELEASE/1.2debug-cidebug-ci-sanitisersstereowalls-data
Last change on this file since 2a0f70e was 2a0f70e, checked in by Olly Betts <olly@…>, 9 years ago

src/datain.c: Add FIXME about warnings when the backclino has range
0-180 degrees.

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