source: TRACY3/branches/tracy3-3.10.1b/tracy/tracy/src/t2lat.cc @ 23

Last change on this file since 23 was 23, checked in by zhangj, 10 years ago

Clean version of Tracy: SoleilVersion at the end of 2011.Use this clean version to find the correct dipole fringe field to have the correct FMAP and FMAPDP of ThomX. Modified files: tpsa_lin.cc, soleillib.cc, prtmfile.cc, rdmfile.cc, read_script.cc, physlib.cc, tracy.cc, t2lat.cc, t2elem.cc, naffutils.cc in /tracy/tracy/src folder; naffutils.h, tracy_global.h, physlib.h, tracy.h, read_script.h, solielilib.h, t2elem.h in /tracy/tracy.inc folder; soltracy.cc in tracy/tools folder
  • Property svn:executable set to *
File size: 121.5 KB
Line 
1/* Tracy-2
2
3J. Bengtsson, CBP, LBL      1990 - 1994   Pascal version
4              SLS, PSI      1995 - 1997
5M. Boege      SLS, PSI      1998          C translation
6L. Nadolski   SOLEIL        2002          Link to NAFF, Radia field maps
7J. Bengtsson  NSLS-II, BNL  2004 -       
8
9*/
10/* Current revision $Revision: 1.20 $
11 On branch $Name:  $
12 Latest change $Date: 2012/01/27 16:47:39 $ by $Author: zhang $
13*/
14
15
16#define NoBmax       200          // maximum number of blocks (NoB)
17#define NoBEmax      9000         /* maximum number of elements in
18                                     a block (Elem_NFam) */
19#define UDImax       100          // max number of user defined constants
20#define LatLLng      (132+1)
21#define Lat_nkw_max  200          // no. of key words
22
23// tables
24
25#define emax         322
26#define emin         (-292)
27#define kmax         15           // max no. of significant digits
28
29#define nmax         LONG_MAX
30
31#define nn           3                    // added by nsrl-ii
32#define tmax         100                  // added by nsrl-ii
33
34#define smax            600   /*size   of string table*/
35#define xmax            LONG_MAX   /*2**8 - 1*/
36
37typedef char Latlinetype[LatLLng];
38
39typedef enum
40{
41  bndsym, defsym, dfcsym, drfsym, elmsym, fcssym, horsym, versym,corkicksym, monsym,
42  qdsym, sexsym,  plus_, minus_, lparent, rparent, eql, comma, lbrack,
43  rbrack, neq, andsy, semicolon, times, rdiv, intcon, realcon, becomes, colon,
44  leq, pwrsym, lss, geo, gtr, period_, charcon, stringcon, ident, geq, lsym,
45  bobrhosym, bobrhovsym, bobrhohsym, kxvsym, kxhsym, phisym, ksym,
46  tsym, t1sym, t2sym,
47  gapsym, thksym, invsym, thnsym,
48  endsym, tsksym, bemsym, corsym, prnsym, tblsym, possym, prmsym,
49  udisym, squote, linsym, mthsym, celsym, matsym, cavsym, symsym, chmsym,
50  cctsym, usesym, andsym, dspsym, kicksym, wglsym, nsym, mrksym,
51  nbdsym, frgsym, latsym, mpsym, dbnsym, kssym, homsym, lmdsym, dtsym, xytsym,
52  vrfsym, harnumsym, frqsym, gstsym, typsym, rollsym, idsym,
53  fnamesym1, fnamesym2, scalingsym1, scalingsym2, fmsym, harmsym, sprsym, recsym, solsym,
54  ff1sym, ff2sym, ffscalingsym, tiltsym
55} Lat_symbol;   /*\*/
56// idsym fnamesym1 fnamesym2 scalingsym added for insertion
57// ring sym added
58// ff1sym ff2sym for quadrupole entrance and exit fringe field
59// ffscalingsym scaling factor for entrance and exit fringe field.  /*J.Zhang
60// tilt added for compatibility with Tracy II
61/* p2c: t2lat.pas, line 52:
62 * Note: Line breaker spent 0.0 seconds, 5000 tries on line 603 [251] */
63
64
65//#define SETBITS  (8*(unsigned)sizeof(long int))
66#define SETBITS  32  /* Length of set is assumed to be 32,
67                        see e.g. DealWithDefns             */
68
69const int  max_set = (solsym-bndsym+1)/SETBITS + 2;
70
71// array for set
72typedef long int  symset[max_set];
73
74typedef  alfa_       Lat_keytype[Lat_nkw_max];
75typedef  Lat_symbol  Lat_ksytype[Lat_nkw_max];
76typedef  Lat_symbol  Lat_spstype[256];
77
78typedef struct _REC_BlockStype
79{
80  partsName Bname;   /* name of a beam line */
81  long  BSTART, BOWARI;
82} _REC_BlockStype;
83
84typedef _REC_BlockStype BlockStype[NoBmax];
85
86typedef long index_;
87
88typedef enum {
89  notyp, ints, reals, bools, chars, arrays, records
90} types;
91
92typedef long typset;
93
94typedef struct _REC_UDItable {
95  partsName  Uname;
96  double     Uvalue;
97} _REC_UDItable;
98
99
100//* Local variables for Lattice_Read: */
101struct LOC_Lattice_Read
102{
103  FILE     **fi, **fo;
104  jmp_buf  _JL9999;
105  long     Symmetry;
106  bool     Ring;    /* true is CELL is a ring */
107
108  long  NoB;        /* Number of defined Blocks */
109  BlockStype BlockS;
110
111  long  Bstack[NoBEmax];
112  long  Bpointer;
113
114  long  UDIC;       /* Number of user defined constants */
115  _REC_UDItable UDItable[UDImax];
116
117  long         nkw;    /* number of key word */
118  Lat_symbol   sym;    /* last symbol read by GetSym*/
119  alfa_        id;     /* identifier from GetSym*/
120  long         inum;   /* integer from GetSym*/
121  double       rnum;   /* double number from GetSym*/
122  char         chin;   /* last character read from source program*/
123  long         cc;     /* character counter*/
124  long         lc;     /* program location counter*/
125  long         ll;     /* length of current line*/
126  long         errpos;
127  Latlinetype  line;
128
129  Lat_keytype  key;
130  Lat_ksytype  ksy;
131  Lat_spstype  sps;
132
133  symset       defbegsys, elmbegsys;
134  bool         skipflag, rsvwd;
135};
136
137
138// Set operations
139
140// val IN s
141int P_inset(register unsigned val, register long *s)
142{
143  register unsigned  bit;
144
145  bit = val % SETBITS; val /= SETBITS;
146  if (val < (unsigned)*s++ && ((1L<<bit) & (unsigned)s[val]))
147    return 1;
148  return 0;
149}
150
151
152// s := s + [val]
153long *P_addset(register long *s, register unsigned val)
154{
155  register long      *sbase = s;
156  register unsigned  bit, size;
157
158  bit = val % SETBITS; val /= SETBITS; size = *s;
159  if (++val > size) {
160    s += size;
161    while (val > size)
162      *++s = 0, size++;
163    *sbase = size;
164  } else
165    s += val;
166  *s |= 1L<<bit;
167  if (size+1 > (unsigned)max_set) {
168    cout << "P_addset: size+1 > max_set " << size+1
169         << "(" << max_set << ")" << endl;
170    exit_(1);
171  }
172  return sbase;
173}
174
175
176// d := s
177long *P_expset(register long *d, register long s)
178{
179  if (s) {
180    d[1] = s;
181    *d = 1;
182  } else
183    *d = 0;
184  return d;
185}
186
187
188// d := s1 + s2
189long *P_setunion(register long *d, register long *s1, register long *s2)
190{
191  long          *dbase = d++;
192  register int  sz1 = *s1++, sz2 = *s2++;
193 
194 while (sz1 > 0 && sz2 > 0) {
195    *d++ = *s1++ | *s2++;
196    sz1--, sz2--;
197  }
198  while (--sz1 >= 0)
199    *d++ = *s1++;
200  while (--sz2 >= 0)
201    *d++ = *s2++;
202  *dbase = d - dbase - 1;
203  return dbase;
204}
205
206
207static long CheckElementtable(const char *name, struct LOC_Lattice_Read *LINK)
208{
209  /* globval.Elem_nFam = Number of parts in a Element */
210  long  i, j, FORLIM;
211
212  j = 0;
213  if (globval.Elem_nFam > Elem_nFamMax) {
214    printf("Elem_nFamMax exceeded: %ld(%d)\n",
215           globval.Elem_nFam, Elem_nFamMax);
216    exit_(1);
217  }
218
219  //  if (strstr(LINK->line,"insertion") != NULL) return 0;
220
221  FORLIM = globval.Elem_nFam;
222  for (i = 1; i <= FORLIM; i++) {
223    if (!strncmp(ElemFam[i-1].ElemF.PName, name, sizeof(partsName)))
224      j = i;
225  }
226  return j;
227}
228
229
230static long CheckBLOCKStable(const char *name, struct LOC_Lattice_Read *LINK)
231{
232  /* NoB = Number of Block defs */
233  long  i, j, FORLIM;
234
235  j = 0;
236  if (LINK->NoB > NoBmax) {
237    printf("** NoBmax exhausted: %ld(%ld)\n", LINK->NoB, (long)NoBmax);
238    return j;
239  }
240  //  if (strstr(LINK->line,"insertion") != NULL) return 0;
241
242  FORLIM = LINK->NoB;
243  for (i = 1; i <= FORLIM; i++) {
244    if (!strncmp(LINK->BlockS[i - 1].Bname, name, sizeof(partsName)))
245      j = i;
246  }
247  return j;
248}
249
250
251//static void InitUDItable(struct LOC_Lattice_Read *LINK)
252//{
253//  LINK->UDIC = 0;
254//}
255
256
257static long CheckUDItable(const char *name, struct LOC_Lattice_Read *LINK)
258{
259  long  i, j, FORLIM;
260
261  j = 0;
262  if (LINK->UDIC > UDImax) {
263    printf("** UDImax exhausted: %ld(%d)\n", LINK->UDIC, UDImax);
264    exit_(1);
265    return j;
266  }
267  FORLIM = LINK->UDIC;
268  for (i = 1L; i <= FORLIM; i++) {
269    if (!strncmp(LINK->UDItable[i - 1L].Uname, name, sizeof(partsName)))
270      j = i;
271  }
272  return j;
273}
274
275
276static void EnterUDItable(const char *name, double X, struct LOC_Lattice_Read *LINK)
277{
278  _REC_UDItable  *WITH;
279
280  LINK->UDIC++;
281  if (LINK->UDIC > UDImax) {
282    printf("** UDImax exhausted: %ld(%d)\n", LINK->UDIC, UDImax);
283    exit_(1);
284    return;
285  }
286  WITH = &LINK->UDItable[LINK->UDIC - 1];
287  memcpy(WITH->Uname, name, sizeof(partsName));
288  WITH->Uvalue = X;
289}
290
291
292static void ModUDItable(long N, double X, struct LOC_Lattice_Read *LINK)
293{
294  _REC_UDItable *WITH;
295
296  WITH = &LINK->UDItable[N - 1];
297  WITH->Uvalue = X;
298}
299
300
301static void RefUDItable(const char *name, double *X, struct LOC_Lattice_Read *LINK)
302{
303  long k;
304
305  k = CheckUDItable(name, LINK);
306  *X = LINK->UDItable[k - 1].Uvalue;
307}
308
309
310//static void PrintUpname(char *name, struct LOC_Lattice_Read *LINK)
311//{
312//  /*(var name:partsname)*/
313//  long i;
314//  char ch;
315//
316//  for (i = 0; i < NameLength; i++) {
317//    ch = name[i];
318//    if ('a' <= ch && ch <= 'z')
319//      ch = _toupper(ch);
320//    putc(ch, *LINK->fo);
321//  }
322//}
323
324
325//static void PrintUpname1(char *name, long *pos,
326//                         struct LOC_Lattice_Read *LINK)
327//{  /*1*/
328//  /*var name : partsname; var pos : integer*/
329//  long i;
330//  char ch;
331//
332//  *pos = 0;
333//  for (i = 0; i < NameLength; i++) {  /*2*/
334//    ch = name[i];
335//    if ('a' <= ch && ch <= 'z')
336//      ch = _toupper(ch);
337//    if (ch != ' ') {
338//      putc(ch, *LINK->fo);
339//      (*pos)++;
340//    }
341//  }  /*2*/
342//}  /*1*/
343
344//void PrintUpname2(char *name, struct LOC_Lattice_Read *LINK)
345//{
346//  /*(var name:partsname)*/
347//  long i;
348//  char ch;
349//
350//  for (i = 0; i < NameLength; i++) {
351//    ch = name[i];
352//    if ('a' <= ch && ch <= 'z')
353//      ch = _toupper(ch);
354//    putc(ch, *LINK->fo);
355//  }
356//}
357
358
359static void abort_(struct LOC_Lattice_Read *LINK)
360{
361  long i;
362
363  for (i = 1; i <= nn; i++)
364    putchar('\007');
365  printf("\n>>>>> error detected in the lattice file <<<<<<\n\n");
366  ErrFlag = true;
367  /*goto 9999*/
368//  printf("% .5E\n", sqrt(-1.0));
369  exit_(1);
370}
371
372
373static void ENDskip(FILE **fo, long *errpos, long *cc, bool *skipflag,
374                    struct LOC_Lattice_Read *LINK)
375{
376  /*underline skips part of input*/
377  while (*errpos < *cc) {
378    putc('-', *fo);
379    (*errpos)++;
380  }
381  *skipflag = false;
382}
383
384
385static void Lat_Error(long n, FILE **fo, long *cc, long *errpos,
386                      struct LOC_Lattice_Read *LINK)
387{
388  if (*errpos != 0L)   /*write(fo, ' ****')*/
389    return;
390  if (*cc > *errpos) {
391    fprintf(*fo, "%*c^%2ld", (int)(*cc - *errpos), ' ', n);
392    *errpos = *cc + 3;
393  }
394}
395
396
397static void Lat_Nextch(FILE **fi, FILE **fo, long *cc, long *ll, long *errpos,
398                       long *lc, char *chin, bool *skipflag, char *line,
399                       struct LOC_Lattice_Read *LINK)
400{
401  if (*cc == *ll) {
402    if (P_eof(*fi)) {
403      fprintf(*fo, "\nprogram incomplete\n");
404      /*errormsg;*/
405      abort_(LINK);
406    }
407
408    if (*errpos != 0) {
409      if (*skipflag)
410        ENDskip(fo, errpos, cc, skipflag, LINK);
411      putc('\n', *fo);
412      *errpos = 0;
413    }
414    /* write(fo, */
415    /*lc: 5, */
416    /* '  ');*/
417    (*lc)++;
418
419    *ll = 0;
420    *cc = 0;
421
422    while (!P_eoln(*fi)) {
423      (*ll)++;
424      if ((*ll) > LatLLng) {
425        printf("Lat_Nextch: LatLLng exceeded %ld (%d)\n", (*ll)-1, LatLLng-1);
426        exit_(1);
427      }
428      *chin = getc(*fi);
429      if (*chin == '\n')
430        *chin = ' ';
431      putc(*chin, *fo);
432      line[*ll - 1] = *chin;
433    }
434    (*ll)++;
435    fscanf(*fi, "%*[^\n]");
436
437    getc(*fi);
438    line[*ll - 1] = ' ';
439    /*read(fi, line[ll]);*/
440    putc('\n', *fo);
441  }
442  (*cc)++;
443  if ((*cc) > LatLLng) {
444    printf("Lat_Nextch: LatLLng exceeded %ld (%d)\n", (*cc), LatLLng);
445    exit_(1);
446  }
447  *chin = line[*cc - 1];
448  /* upper case to lower case */
449  if (isupper(*chin))
450    *chin = _tolower(*chin);
451  /* tab */
452  if (*chin == '\t')
453    *chin = ' ';
454}  /* Lat_Nextch */
455
456
457static void Lat_errorm(const char *cmnt, FILE **fi, FILE **fo, long *cc,
458                       long *ll, long *errpos, long *lc, char *chin,
459                       bool *skipflag, char *line,
460                       struct LOC_Lattice_Read *LINK)
461{
462  /*write(fo, ' ****')*/
463  if (*cc > *errpos) {
464    fprintf(*fo, "%*c^%.80s", (int)(*cc - *errpos), ' ', cmnt);
465    *errpos = *cc + 3;
466  }
467  while (!P_eof(*fi))
468    Lat_Nextch(fi, fo, cc, ll, errpos, lc, chin, skipflag, line, LINK);
469  ErrFlag = true;
470  abort_(LINK);
471}
472
473/* Local variables for Lat_GetSym: */
474struct LOC_Lat_GetSym
475{
476  struct LOC_Lattice_Read  *LINK;
477  FILE    **fi, **fo;
478  long    *cc, *ll, *errpos, *lc, emax_, emin_;
479  char    *chin;
480  double  *rnum;
481  bool    *skipflag;
482  char    *line;
483  long    k, e;
484};
485
486
487static void NextCh(struct LOC_Lat_GetSym *LINK)
488{
489  Lat_Nextch(LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos, LINK->lc,
490             LINK->chin, LINK->skipflag, LINK->line, LINK->LINK);
491}
492
493
494static void readscale(struct LOC_Lat_GetSym *LINK)
495{
496  long s, sign;
497
498  /* readscale  */
499  NextCh(LINK);
500  while (*LINK->chin == ' ')
501    NextCh(LINK);
502  sign = 1;
503  s = 0;
504  if (*LINK->chin == '+')
505    NextCh(LINK);
506  else if (*LINK->chin == '-') {
507    NextCh(LINK);
508    sign = -1;
509  }
510  if (!isdigit(*LINK->chin))
511    Lat_Error(40, LINK->fo, LINK->cc, LINK->errpos, LINK->LINK);
512  else {
513    do {
514      s = s * 10 + *LINK->chin - '0';
515      NextCh(LINK);
516    }
517    while (isdigit(*LINK->chin));
518  }
519  LINK->e += s * sign;
520}
521
522
523static void adjustscale(struct LOC_Lat_GetSym *LINK)
524{
525  long s;
526  double d, t;
527
528  /*  adjustscale  */
529  if (LINK->k + LINK->e > LINK->emax_) {
530    Lat_Error(21, LINK->fo, LINK->cc, LINK->errpos, LINK->LINK);
531    return;
532  }
533  if (LINK->k + LINK->e < LINK->emin_) {
534    *LINK->rnum = 0.0;
535    return;
536  }
537  s = abs(LINK->e);
538  t = 1.0;
539  d = 10.0;
540  do {
541    while (!(s & 1)) {
542      s /= 2;
543      d *= d;
544    }
545    s--;
546    t = d * t;
547  }
548  while (s != 0);
549
550  if (LINK->e >= 0)
551    *LINK->rnum *= t;
552  else
553    *LINK->rnum /= t;
554}
555
556
557static void Lat_GetSym(FILE **fi_, FILE **fo_, long *cc_, long *ll_,
558                       long *errpos_, long *lc_, long *nkw, long *inum,
559                       long emax__, long emin__, long kmax_, long nmax_,
560                       char *chin_, char *id, double *rnum_, bool *skipflag_,
561                       bool *rsvwd, char *line_, Lat_symbol *sym, alfa_ *key,
562                       Lat_symbol *ksy, Lat_symbol *sps,
563                       struct LOC_Lattice_Read *LINK)
564{  /*GetSym*/
565  struct LOC_Lat_GetSym V;
566
567  long i, j, mysign;
568  bool parsename;
569
570
571  V.LINK = LINK;
572  /*GetSym*/
573  V.fi = fi_;
574  V.fo = fo_;
575  V.cc = cc_;
576  V.ll = ll_;
577  V.errpos = errpos_;
578  V.lc = lc_;
579  V.emax_ = emax__;
580  V.emin_ = emin__;
581  V.chin = chin_;
582  V.rnum = rnum_;
583  V.skipflag = skipflag_;
584  V.line = line_;
585  *rsvwd = false;
586  mysign = 1;
587  parsename = false;
588 _L1:
589  while (*V.chin == ' ')
590
591    NextCh(&V);
592
593  switch (*V.chin) {
594
595  case 'a':
596  case 'b':
597  case 'c':
598  case 'd':
599  case 'e':
600  case 'f':
601  case 'g':
602  case 'h':
603  case 'i':
604  case 'j':
605  case 'k':
606  case 'l':
607  case 'm':
608  case 'n':
609  case 'o':
610  case 'p':
611  case 'q':
612  case 'r':
613  case 's':
614  case 't':
615  case 'u':
616  case 'v':
617  case 'w':
618  case 'x':
619  case 'y':
620  case 'z':
621  case '"':  /*identifier or wordsymbol*/
622    V.k = 0;
623    memcpy(id, "               ", sizeof(alfa_));
624    do {
625      if (*V.chin == '"')
626        parsename = !parsename;
627      if (V.k < NameLength) {
628        V.k++; id[V.k - 1] = *V.chin;
629      } else {
630        printf("In Lat_GetSym, symbol: %s too long, max value is %d\n",
631                id, NameLength);
632        exit_(1);
633      }
634      NextCh(&V);
635    } while (parsename || *V.chin == '_' || islower(*V.chin) ||
636             isdigit(*V.chin));
637
638    /*writeln(fo, 'GetSym detected: id=', id);*/
639
640    i = 1;
641    j = *nkw;   /*binary search*/
642    do {
643      V.k = (i + j) / 2;
644      if (strncmp(id, key[V.k - 1], sizeof(alfa_)) <= 0)
645        j = V.k - 1;
646      if (strncmp(id, key[V.k - 1], sizeof(alfa_)) >= 0)
647        i = V.k + 1;
648      /*  writeln(fo, '  bunary: id=', id, '  key[', k:3, ']=', key[k],
649          'i=', i:4, ' j=', j:4, ' k=', k:4, ' i-1-j=', (i-1-j):4);*/
650    }
651    while (i <= j);
652    if (i - 1 > j) {
653      *sym = ksy[V.k - 1]; *rsvwd = true;
654      /*  writeln(fo, 'GetSym detected reserved word: id=', id,
655          '  k=', k:4, '  key[', k:4, ']=', key[k]);*/
656    } else {
657      if (!strncmp(id, "t              ", sizeof(alfa_)))
658        *sym = tsym;
659      else if (!strncmp(id, "gap            ", sizeof(alfa_)))
660        *sym = gapsym;
661      else if (!strncmp(id, "l              ", sizeof(alfa_)))
662        *sym = lsym;
663      else if (!strncmp(id, "n              ", sizeof(alfa_)))
664        *sym = nsym;
665      else if (!strncmp(id, "bobrho         ", sizeof(alfa_)))
666        *sym = bobrhosym;
667      else if (!strncmp(id, "bobrhov        ", sizeof(alfa_)))
668        *sym = bobrhovsym;
669      else if (!strncmp(id, "bobrhoh        ", sizeof(alfa_)))
670        *sym = bobrhohsym;
671      else if (!strncmp(id, "kxv            ", sizeof(alfa_)))
672        *sym = kxvsym;
673      else if (!strncmp(id, "kxh            ", sizeof(alfa_)))
674        *sym = kxhsym;
675      else if (!strncmp(id, "phi            ", sizeof(alfa_)))
676        *sym = phisym;
677      else if (!strncmp(id, "k              ", sizeof(alfa_)))
678        *sym = ksym;
679      else if (!strncmp(id, "harnum         ", sizeof(alfa_)))
680        *sym = harnumsym;
681      else
682        *sym = ident;
683    }
684    break;
685
686  case '0':
687  case '1':
688  case '2':
689  case '3':
690  case '4':
691  case '5':
692  case '6':
693  case '7':
694  case '8':
695  case '9':  /*number*/
696    V.k = 0;
697    *inum = 0;
698    *sym = intcon;
699    do {
700      *inum = *inum * 10 + *V.chin - '0';
701      V.k++;
702      NextCh(&V);
703    }
704    while (isdigit(*V.chin));
705    if (V.k > kmax_ || *inum > nmax_) {
706      Lat_Error(21, V.fo, V.cc, V.errpos, LINK);
707      *inum = 0;
708      V.k = 0;
709    }
710    if (*V.chin == '.') {
711      NextCh(&V);
712      if (*V.chin == '.')
713        *V.chin = ':';
714      else {
715        *sym = realcon;
716        *V.rnum = *inum;
717        V.e = 0;
718        while (isdigit(*V.chin)) {
719          V.e--;
720          *V.rnum = 10.0 * *V.rnum + *V.chin - '0';
721          NextCh(&V);
722        }
723        while (*V.chin == ' ')
724          NextCh(&V);
725        if (V.e == 0)
726          Lat_Error(40, V.fo, V.cc, V.errpos, LINK);
727        if (*V.chin == 'd' || *V.chin == 'D' || *V.chin == 'e' ||
728            *V.chin == 'E')
729          readscale(&V);
730        if (V.e != 0)
731          adjustscale(&V);
732      }
733    } else {
734      if (*V.chin == 'd' || *V.chin == 'D' || *V.chin == 'e' ||
735          *V.chin == 'E') {
736        *sym = realcon;
737        *V.rnum = *inum;
738        V.e = 0;
739        readscale(&V);
740        if (V.e != 0)
741          adjustscale(&V);
742      }
743    }
744    if (*sym == intcon)
745      *inum *= mysign;
746    else {
747      if (*sym == realcon)
748        *V.rnum = mysign * *V.rnum;
749    }
750    break;
751
752  case ':':   /*, col*/
753    NextCh(&V);
754    if (*V.chin == '=') {
755      *sym = becomes;
756      NextCh(&V);
757    }
758    else
759      *sym = colon;
760    break;
761
762  case '<':
763    NextCh(&V);
764    if (*V.chin == '=')
765      {
766        *sym = leq;
767        NextCh(&V);
768      }
769    else {
770      if (*V.chin == '>') {
771        *sym = neq;
772        NextCh(&V);
773      }
774      else
775        *sym = lss;
776    }
777    break;
778
779  case '>':
780    NextCh(&V);
781    if (*V.chin == '=') {
782      *sym = geq;
783      NextCh(&V);
784    }
785    else
786      *sym = gtr;
787    break;
788
789  case '.':
790    NextCh(&V);
791    *sym = period_;
792    break;
793
794  case '*':
795    NextCh(&V);
796    if (*V.chin == '*') {
797      *sym = pwrsym;
798      NextCh(&V);
799    }
800    else
801      *sym = times;
802    break;
803
804  case '{':
805    do {
806      NextCh(&V);
807    }
808    while (*V.chin != '}');
809    NextCh(&V);
810    goto _L1;
811    break;
812
813  case '+':
814  case '-':
815  case '/':
816  case '(':
817  case ')':
818  case '=':
819  case ',':
820  case ';':
821  case '[':
822  case ']':
823  case '\'':
824    *sym = sps[*V.chin];
825    /* IF chin='+' THEN BEGIN nextch; goto 1 END ELSE
826       IF chin='-' THEN BEGIN nextch; mysign:=-1; goto 1 END ELSE*/
827    NextCh(&V);
828    break;
829
830  case '$':
831  case '!':
832  case '@':
833  case '?':
834  case '_':
835  case '&':
836  case '\\':
837  case '^':
838    Lat_Error(24L, V.fo, V.cc, V.errpos, LINK);
839    NextCh(&V);
840    goto _L1;
841    break;
842
843  default:
844    Lat_Error(24L, V.fo, V.cc, V.errpos, LINK);
845    NextCh(&V);
846    goto _L1;
847    break;
848  }
849}
850
851/* Local variables for Lat_EVAL: */
852struct LOC_Lat_EVAL
853{
854  struct LOC_Lattice_Read *LINK;
855  FILE **fi, **fo;
856  long *cc, *ll, *errpos, *lc, *nkw, *inum, emax_, emin_, kmax_, nmax_;
857  char *chin;
858  char *id;
859  double *rnum;
860  bool *skipflag, *rsvwd;
861  char *line;
862  Lat_symbol *sym;
863  alfa_ *key;
864  Lat_symbol *ksy;
865  Lat_symbol *sps;
866  jmp_buf _JL999;
867
868  double S[tmax + 1];
869  long t;
870  symset facbegsys;
871};
872
873static void Expression(struct LOC_Lat_EVAL *LINK);
874
875
876static void GetSym(struct LOC_Lat_EVAL *LINK)
877{
878  /* reads next symbol  */
879  Lat_GetSym(LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos, LINK->lc,
880             LINK->nkw, LINK->inum, LINK->emax_, LINK->emin_, LINK->kmax_,
881             LINK->nmax_, LINK->chin, LINK->id, LINK->rnum, LINK->skipflag,
882             LINK->rsvwd, LINK->line, LINK->sym, LINK->key, LINK->ksy,
883             LINK->sps, LINK->LINK);
884}
885
886
887static void errorm(const char *cmnt, struct LOC_Lat_EVAL *LINK)
888{
889  Lat_errorm(cmnt, LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos,
890             LINK->lc, LINK->chin, LINK->skipflag, LINK->line, LINK->LINK);
891}
892
893static void test(long *s1, const char *cmnt, struct LOC_Lat_EVAL *LINK)
894{
895  if (!P_inset(*LINK->sym, s1))
896    errorm(cmnt, LINK);
897}
898
899
900static void getest(long *s1, const char *cmnt, struct LOC_Lat_EVAL *LINK)
901{
902  GetSym(LINK);
903  if (!P_inset(*LINK->sym, s1))
904    errorm(cmnt, LINK);
905}
906
907
908static double ArcSin(double x, struct LOC_Lat_EVAL *LINK)
909{
910  if (fabs(x) > 1e0)
911    longjmp(LINK->_JL999, 1);
912  if (x == 1e0)
913    return (M_PI / 2e0);
914  else if (x == -1e0)
915    return (M_PI / -2e0);
916  else
917    return atan(x / sqrt(1e0 - x * x));
918}
919
920
921static double ArcCos(double x, struct LOC_Lat_EVAL *LINK)
922{
923  if (fabs(x) > 1e0)
924    longjmp(LINK->_JL999, 1);
925  if (x == 1e0)
926    return 0e0;
927  else if (x == -1e0)
928    return M_PI;
929  else
930    return atan(sqrt(1e0 - x * x) / x);
931}
932
933
934static void writes(struct LOC_Lat_EVAL *LINK)
935{
936  /*writeln('PUSH:  s[', t:3, ']=', s[t]);*/
937}
938
939
940static void PUSH(double x, struct LOC_Lat_EVAL *LINK)
941{
942  LINK->t++;
943  if (LINK->t == tmax)
944    {
945      printf("** Lat_Eval: stack overflow\n");
946      longjmp(LINK->_JL999, 1);
947    }
948  LINK->S[LINK->t] = x;
949  writes(LINK);
950}
951
952
953static double BlockLength(long ii, struct LOC_Lat_EVAL *LINK)
954{
955  long    k1, k2, k3;
956  double  S;
957
958  S = 0.0;
959  if (ii == 0)
960    return S;
961  k2 = LINK->LINK->BlockS[ii - 1].BSTART;
962  k3 = LINK->LINK->BlockS[ii - 1].BOWARI;
963  for (k1 = k2 - 1; k1 < k3; k1++)
964    S += ElemFam[LINK->LINK->Bstack[k1] - 1].ElemF.PL;
965  return S;
966}
967
968/* Local variables for Expression: */
969struct LOC_Expression
970{
971  struct LOC_Lat_EVAL *LINK;
972};
973
974/* Local variables for Term: */
975struct LOC_Term
976{
977  struct LOC_Expression *LINK;
978};
979
980/* Local variables for Factor: */
981struct LOC_Factor
982{
983  struct LOC_Term *LINK;
984  long i;
985};
986
987
988static double GetKparm(long direction, struct LOC_Factor *LINK)
989{
990  double Result;
991  symset SET;
992
993  getest(P_expset(SET, 1 << ((long)lbrack)), "<[> expected",
994         LINK->LINK->LINK->LINK);
995  GetSym(LINK->LINK->LINK->LINK);
996  Expression(LINK->LINK->LINK->LINK);
997  test(P_expset(SET, 1 << ((long)rbrack)), "<]> expected",
998       LINK->LINK->LINK->LINK);
999  if (direction == 1)
1000    Result = ElemFam[LINK->i-1].ElemF.M->PBpar[(long)((long)LINK->
1001             LINK->LINK->LINK->S[LINK->LINK->LINK->LINK->t])+HOMmax];
1002  else
1003    Result = ElemFam[LINK->i-1].ElemF.M->PBpar[HOMmax-(long)LINK->
1004            LINK->LINK->LINK->S[LINK->LINK->LINK->LINK->t]];
1005  LINK->LINK->LINK->LINK->t--;
1006  /* GetSym;*/
1007  return Result;
1008}
1009
1010
1011static void Factor(struct LOC_Term *LINK)
1012{
1013  struct LOC_Factor V;
1014  double x = 0.0;
1015  partsName fname;
1016  long SET[(long)period_ / 32 + 2];
1017  elemtype *WITH;
1018  MpoleType *WITH1;
1019  symset SET1;
1020  long SET2[(long)lsym / 32 + 2];
1021
1022  V.LINK = LINK;
1023  /* factor */
1024  if (!P_inset(*LINK->LINK->LINK->sym, LINK->LINK->LINK->facbegsys))
1025    longjmp(LINK->LINK->LINK->_JL999, 1);
1026  /*while sym in facbegsys do*/
1027  if (*LINK->LINK->LINK->sym == ident) {  /*1:  of ident */
1028    V.i = CheckUDItable(LINK->LINK->LINK->id, LINK->LINK->LINK->LINK);
1029    if (V.i != 0) {   /* UDI */
1030      x = LINK->LINK->LINK->LINK->UDItable[V.i - 1].Uvalue;
1031      PUSH(x, LINK->LINK->LINK);
1032      GetSym(LINK->LINK->LINK);
1033    } else {
1034      V.i = CheckElementtable(LINK->LINK->LINK->id, LINK->LINK->LINK->LINK);
1035      if (V.i != 0) {
1036        getest(P_addset(P_expset(SET, 0), (long)period_), "<.> expected",
1037               LINK->LINK->LINK);
1038        /*--> new */
1039        GetSym(LINK->LINK->LINK);
1040        memcpy(fname, LINK->LINK->LINK->id, sizeof(alfa_));
1041        memset(fname + sizeof(alfa_), ' ',
1042               sizeof(partsName) - sizeof(alfa_));
1043        WITH = &ElemFam[V.i - 1].ElemF;
1044        WITH1 = WITH->M;
1045        if (!strncmp(fname, "l              ", sizeof(partsName)))
1046          x = WITH->PL;
1047        else if (!strncmp(fname, "t              ", sizeof(partsName))) {
1048          if (WITH->PL != 0.0)
1049            x = WITH1->Pirho * WITH->PL * 180.0 / M_PI;
1050          else
1051            x = WITH1->Pirho * 180.0 / M_PI;
1052        } else if (!strncmp(fname, "t1             ", sizeof(partsName)))
1053          x = WITH1->PTx1;
1054        else if (!strncmp(fname, "t2             ", sizeof(partsName)))
1055          x = WITH1->PTx2;
1056        else if (!strncmp(fname, "gap            ", sizeof(partsName)))
1057          x = WITH1->Pgap;
1058        else if (!strncmp(fname, "roll           ", sizeof(partsName))) {
1059          if (WITH->Pkind == Mpole)
1060            x = WITH1->PdTpar;
1061          else if (WITH->Pkind == Wigl)
1062            x = WITH1->PdTpar;
1063        } else if (!strncmp(fname, "n              ", sizeof(partsName))) {
1064          if (((1 << ((long)WITH->Pkind)) &
1065               ((1 << ((long)Mpole)) | (1 << ((long)Wigl)))) != 0)
1066            x = WITH1->PN;
1067        } else if (!strncmp(fname, "b              ", sizeof(partsName)))
1068          x = GetKparm(1, &V);
1069        else if (!strncmp(fname, "a              ", sizeof(partsName)))
1070          x = GetKparm(2L, &V);
1071        else {
1072          /* error detected */
1073          getest(P_expset(SET1, 0), "  illegal extension...",
1074                 LINK->LINK->LINK);
1075        }
1076        PUSH(x, LINK->LINK->LINK);
1077        GetSym(LINK->LINK->LINK);
1078      } else {
1079        V.i = CheckBLOCKStable(LINK->LINK->LINK->id, LINK->LINK->LINK->LINK);
1080        if (V.i != 0) {
1081          getest(P_addset(P_expset(SET, 0), (long)period_), "<.> expected",
1082                 LINK->LINK->LINK);
1083          getest(P_addset(P_expset(SET2, 0), (long)lsym), "illegal component",
1084                 LINK->LINK->LINK);
1085          x = BlockLength(V.i, LINK->LINK->LINK);
1086          PUSH(x, LINK->LINK->LINK);
1087          GetSym(LINK->LINK->LINK);
1088        } else {  /*4: function ?*/
1089          memcpy(fname, LINK->LINK->LINK->id, sizeof(alfa_));
1090          memset(fname + sizeof(alfa_), ' ',
1091                 sizeof(partsName) - sizeof(alfa_));
1092          GetSym(LINK->LINK->LINK);
1093          switch (*LINK->LINK->LINK->sym) {   /*5*/
1094
1095          case semicolon:
1096            GetSym(LINK->LINK->LINK);
1097            break;
1098
1099          case lparent:  /*6: of lparent*/
1100            GetSym(LINK->LINK->LINK);
1101            Expression(LINK->LINK->LINK);
1102            if (!strncmp(fname, "sin            ", sizeof(partsName)))
1103              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1104                sin(LINK->LINK->LINK->S[LINK->LINK->LINK->t]);
1105            else if (!strncmp(fname, "cos            ", sizeof(partsName)))
1106              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1107                cos(LINK->LINK->LINK->S[LINK->LINK->LINK->t]);
1108            else if (!strncmp(fname, "tan            ", sizeof(partsName)))
1109              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1110                tan(LINK->LINK->LINK->S[LINK->LINK->LINK->t]);
1111            else if (!strncmp(fname, "arcsin         ", sizeof(partsName)))
1112              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1113                ArcSin(LINK->LINK->LINK->S[LINK->LINK->LINK->t],
1114                       LINK->LINK->LINK);
1115            else if (!strncmp(fname, "arccos         ", sizeof(partsName)))
1116              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1117                ArcCos(LINK->LINK->LINK->S[LINK->LINK->LINK->t],
1118                       LINK->LINK->LINK);
1119            else if (!strncmp(fname, "arctan         ", sizeof(partsName)))
1120              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1121                atan(LINK->LINK->LINK->S[LINK->LINK->LINK->t]);
1122            else if (!strncmp(fname, "sqrt           ", sizeof(partsName)))
1123              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1124                sqrt(LINK->LINK->LINK->S[LINK->LINK->LINK->t]);
1125            else if (!strncmp(fname, "log            ", sizeof(partsName)))
1126              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1127                log(LINK->LINK->LINK->S[LINK->LINK->LINK->t]);
1128            else if (!strncmp(fname, "exp            ", sizeof(partsName)))
1129              LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1130                exp(LINK->LINK->LINK->S[LINK->LINK->LINK->t]);
1131            writes(LINK->LINK->LINK);
1132            if (*LINK->LINK->LINK->sym == rparent)
1133              GetSym(LINK->LINK->LINK);
1134            else
1135              longjmp(LINK->LINK->LINK->_JL999, 1);
1136            break;
1137            /*6:of lparent*/
1138          default:
1139            break;
1140          }/*5: of case */
1141        }  /*4: of function?*/
1142      }
1143    }
1144  }  /*1: of ident*/
1145  else if (*LINK->LINK->LINK->sym == realcon) {
1146    PUSH(*LINK->LINK->LINK->rnum, LINK->LINK->LINK);
1147    GetSym(LINK->LINK->LINK);
1148  } else if (*LINK->LINK->LINK->sym == intcon) {
1149    x = *LINK->LINK->LINK->inum;
1150    PUSH(x, LINK->LINK->LINK);
1151    GetSym(LINK->LINK->LINK);
1152  } else if (*LINK->LINK->LINK->sym == lparent) {
1153    GetSym(LINK->LINK->LINK);
1154    Expression(LINK->LINK->LINK);
1155    if (*LINK->LINK->LINK->sym == rparent)
1156      GetSym(LINK->LINK->LINK);
1157    else
1158      longjmp(LINK->LINK->LINK->_JL999, 1);
1159  } else
1160    longjmp(LINK->LINK->LINK->_JL999, 1);
1161
1162  if (*LINK->LINK->LINK->sym != pwrsym)
1163    return;
1164
1165  GetSym(LINK->LINK->LINK);
1166  if (*LINK->LINK->LINK->sym != intcon)
1167    longjmp(LINK->LINK->LINK->_JL999, 1);
1168  LINK->LINK->LINK->S[LINK->LINK->LINK->t] =
1169    pow(LINK->LINK->LINK->S[LINK->LINK->LINK->t],
1170        (double)(*LINK->LINK->LINK->inum));
1171  GetSym(LINK->LINK->LINK);
1172}
1173
1174
1175static void Term(struct LOC_Expression *LINK)
1176{
1177  struct LOC_Term V;
1178  Lat_symbol mulop;
1179
1180  V.LINK = LINK;
1181  /* term */
1182  Factor(&V);
1183  while ((unsigned int)(*LINK->LINK->sym) < 32 &&
1184         ((1 << ((long)(*LINK->LINK->sym))) &
1185          ((1 << ((long)times)) | (1 << ((long)rdiv)))) != 0) {
1186    mulop = *LINK->LINK->sym;
1187    GetSym(LINK->LINK);
1188    Factor(&V);
1189    if (mulop == times) {
1190      LINK->LINK->S[LINK->LINK->t - 1] *= LINK->LINK->S[LINK->LINK->t];
1191      LINK->LINK->t--;
1192      writes(LINK->LINK);
1193    } else {
1194      if (mulop == rdiv) {
1195        LINK->LINK->S[LINK->LINK->t - 1] /= LINK->LINK->S[LINK->LINK->t];
1196        LINK->LINK->t--;
1197        writes(LINK->LINK);
1198      }
1199    }
1200  }
1201}
1202
1203
1204static void Expression(struct LOC_Lat_EVAL *LINK)
1205{
1206  struct LOC_Expression V;
1207  Lat_symbol addop;
1208
1209  V.LINK = LINK;
1210  /* Expression */
1211  if ((unsigned int)(*LINK->sym) < 32 &&
1212      ((1 << ((long)(*LINK->sym))) &
1213       ((1 << ((long)plus_)) | (1 << ((long)minus_)))) != 0) {
1214    addop = *LINK->sym;
1215    GetSym(LINK);
1216    Term(&V);
1217    if (addop == minus_)
1218      LINK->S[LINK->t] = -LINK->S[LINK->t];
1219  } else
1220    Term(&V);
1221
1222  while ((unsigned int)(*LINK->sym) < 32 &&
1223         ((1 << ((long)(*LINK->sym))) &
1224          ((1 << ((long)plus_)) | (1 << ((long)minus_)))) != 0) {
1225    addop = *LINK->sym;
1226    GetSym(LINK);
1227    Term(&V);
1228    if (addop == plus_) {
1229      LINK->S[LINK->t - 1] += LINK->S[LINK->t];
1230      LINK->t--;
1231      writes(LINK);
1232    } else {
1233      if (addop == minus_) {
1234        LINK->S[LINK->t - 1] -= LINK->S[LINK->t];
1235        LINK->t--;
1236        writes(LINK);
1237      }
1238    }
1239  }
1240}
1241
1242/******************************
1243 *                           *
1244 *        E V A L            *
1245 *                           *
1246 ******************************/
1247
1248static double Lat_EVAL(FILE **fi_, FILE **fo_, long *cc_, long *ll_,
1249                       long *errpos_,
1250                       long *lc_, long *nkw_, long *inum_, long emax__,
1251                       long emin__,
1252                       long kmax__, long nmax__, char *chin_, char *id_,
1253                       double *rnum_,
1254                       bool *skipflag_, bool *rsvwd_, char *line_,
1255                       Lat_symbol *sym_,
1256                       alfa_ *key_, Lat_symbol *ksy_, Lat_symbol *sps_,
1257                       struct LOC_Lattice_Read *LINK)
1258{  /* eval */
1259  struct LOC_Lat_EVAL V;
1260  double Result = 0.0;
1261  symset SET;
1262
1263  V.LINK = LINK;
1264  V.fi = fi_;
1265  V.fo = fo_;
1266  V.cc = cc_;
1267  V.ll = ll_;
1268  V.errpos = errpos_;
1269  V.lc = lc_;
1270  V.nkw = nkw_;
1271  V.inum = inum_;
1272  V.emax_ = emax__;
1273  V.emin_ = emin__;
1274  V.kmax_ = kmax__;
1275  V.nmax_ = nmax__;
1276  V.chin = chin_;
1277  V.id = id_;
1278  V.rnum = rnum_;
1279  V.skipflag = skipflag_;
1280  V.rsvwd = rsvwd_;
1281  V.line = line_;
1282  V.sym = sym_;
1283  V.key = key_;
1284  V.ksy = ksy_;
1285  V.sps = sps_;
1286  if (setjmp(V._JL999))
1287    goto _L999;
1288  P_addset(P_expset(V.facbegsys, 0), (long)intcon);
1289  P_addset(V.facbegsys, (long)realcon);
1290  P_addset(V.facbegsys, (long)ident);
1291  P_addset(V.facbegsys, (long)lparent);
1292  GetSym(&V);
1293  V.t = 0;
1294  Expression(&V);
1295  if (V.t != 1)
1296    goto _L999;
1297  Result = V.S[1];
1298  goto _L888;
1299
1300 _L999:
1301  ErrFlag = true;
1302  test(P_expset(SET, 0), "** Lat_Eval: error",
1303       &V);
1304 _L888:   /* exit */
1305
1306  return Result;
1307}
1308
1309
1310/* Local variables for Lat_ProcessBlockInput: */
1311struct LOC_Lat_ProcessBlockInput
1312{
1313  struct LOC_Lattice_Read *LINK;
1314  FILE **fi, **fo;
1315  long *cc, *ll, *errpos, *lc, *nkw, *inum, emax_, emin_, kmax_, nmax_;
1316  char *chin;
1317  char *id;
1318  double *rnum;
1319  bool *skipflag, *rsvwd;
1320  char *line;
1321  Lat_symbol *sym;
1322  alfa_ *key;
1323  Lat_symbol *ksy;
1324  Lat_symbol *sps;
1325};
1326
1327static void GetBlock(struct LOC_Lat_ProcessBlockInput *LINK);
1328
1329
1330static void errorm_(const char *cmnt, struct LOC_Lat_ProcessBlockInput *LINK)
1331{
1332  Lat_errorm(cmnt, LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos,
1333             LINK->lc, LINK->chin, LINK->skipflag, LINK->line, LINK->LINK);
1334}
1335
1336
1337static void GetSym_(struct LOC_Lat_ProcessBlockInput *LINK)
1338{
1339  Lat_GetSym(LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos, LINK->lc,
1340             LINK->nkw, LINK->inum, LINK->emax_, LINK->emin_, LINK->kmax_,
1341             LINK->nmax_, LINK->chin, LINK->id, LINK->rnum, LINK->skipflag,
1342             LINK->rsvwd, LINK->line, LINK->sym, LINK->key, LINK->ksy,
1343             LINK->sps, LINK->LINK);
1344}
1345
1346
1347static void test_(long *s1, const char *cmnt, struct LOC_Lat_ProcessBlockInput *LINK)
1348{
1349  /*test*/
1350  if (!P_inset(*LINK->sym, s1))
1351    errorm_(cmnt, LINK);
1352}
1353
1354
1355static void getest_(long *s1, const char *cmnt,
1356                    struct LOC_Lat_ProcessBlockInput *LINK)
1357{
1358  /*test*/
1359  GetSym_(LINK);
1360  if (!P_inset(*LINK->sym, s1))
1361    errorm_(cmnt, LINK);
1362}
1363
1364
1365double EVAL(struct LOC_Lat_ProcessBlockInput *LINK)
1366{
1367  return (Lat_EVAL(LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos,
1368                   LINK->lc, LINK->nkw, LINK->inum, LINK->emax_, LINK->emin_,
1369                   LINK->kmax_, LINK->nmax_, LINK->chin, LINK->id, LINK->rnum,
1370                   LINK->skipflag, LINK->rsvwd, LINK->line, LINK->sym,
1371                   LINK->key, LINK->ksy, LINK->sps, LINK->LINK));
1372}
1373
1374
1375static void DeBlock(long ii, long k4, struct LOC_Lat_ProcessBlockInput *LINK)
1376{
1377  long  k1, k2, k3, k5;
1378
1379  k2 = LINK->LINK->BlockS[ii - 1].BSTART;
1380  k3 = LINK->LINK->BlockS[ii - 1].BOWARI;
1381  for (k5 = 1; k5 <= k4; k5++) {
1382    for (k1 = k2 - 1; k1 < k3; k1++) {  /*11*/
1383      LINK->LINK->Bpointer++;
1384      if (LINK->LINK->Bpointer >= NoBEmax) {
1385        printf("** DeBlock: NoBEmax exceeded %ld (%d)\n",
1386               LINK->LINK->Bpointer, NoBEmax);
1387        exit(1);
1388      }
1389      LINK->LINK->Bstack[LINK->LINK->Bpointer - 1] = LINK->LINK->Bstack[k1];
1390    }  /*11*/
1391  }
1392  GetSym_(LINK);
1393  if (*LINK->sym == comma)
1394    GetSym_(LINK);
1395}
1396
1397/* Local variables for GetBlock: */
1398struct LOC_GetBlock
1399{
1400  struct LOC_Lat_ProcessBlockInput *LINK;
1401};
1402
1403
1404static void InsideParent(long k4, struct LOC_GetBlock *LINK)
1405{
1406  long    b, b1, b2, k1;
1407  symset  SET;
1408
1409  b1 = LINK->LINK->LINK->Bpointer + 1;
1410  GetSym_(LINK->LINK);
1411  GetBlock(LINK->LINK);
1412  b2 = LINK->LINK->LINK->Bpointer;
1413  if (k4 >= 2) {
1414    for (k1 = 2; k1 <= k4; k1++) {
1415      for (b = b1 - 1; b < b2; b++) {
1416        LINK->LINK->LINK->Bpointer++;
1417        if (LINK->LINK->LINK->Bpointer >= NoBEmax) {
1418          printf("** InsideParent: NoBEmax exceeded %ld (%d)\n",
1419                 LINK->LINK->LINK->Bpointer, NoBEmax);
1420          exit(1);
1421        }
1422        LINK->LINK->LINK->Bstack[LINK->LINK->LINK->Bpointer - 1] = LINK->
1423          LINK->LINK->Bstack[b];
1424      }
1425    }
1426  }
1427  test_(P_expset(SET, 1 << ((long)rparent)), "<)> expected",
1428        LINK->LINK);
1429  getest_(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon)) |
1430                   (1 << ((long)rparent))), "<, > or <;> expected",
1431          LINK->LINK);
1432  if (*LINK->LINK->sym == comma)
1433    GetSym_(LINK->LINK);
1434}
1435
1436
1437static void Doinverse(struct LOC_GetBlock *LINK)
1438{
1439  long    b, b1, b2, b3, k1, k4;
1440  symset  SET;
1441  long    FORLIM;
1442
1443  getest_(P_expset(SET, 1 << ((long)lparent)), "<(> expected after INV",
1444          LINK->LINK);
1445  b1 = LINK->LINK->LINK->Bpointer + 1;
1446  GetSym_(LINK->LINK);
1447  GetBlock(LINK->LINK);
1448  b2 = LINK->LINK->LINK->Bpointer;
1449  /* Bug fix: INV(A, B, ...) for 2 elements
1450  k4 = b2 - b1 */
1451  k4 = b2 - b1 + 1;
1452  if (k4 >= 2) {
1453    k4 /= 2; b3 = b1 + k4; k1 = 0;
1454    FORLIM = b3;
1455  /* Bug fix: INV(A, B, ...) for 2 elements
1456    for (b = b1 - 1; b < FORLIM; b++) { */
1457    for (b = b1 - 1; b < FORLIM-1; b++) {
1458      b3 = LINK->LINK->LINK->Bstack[b];
1459      LINK->LINK->LINK->Bstack[b] = LINK->LINK->LINK->Bstack[b2-k1-1];
1460      LINK->LINK->LINK->Bstack[b2-k1-1] = b3;
1461      k1++;
1462    }
1463  }
1464  test_(P_expset(SET, 1 << ((long)rparent)), "<)> expected", LINK->LINK);
1465  getest_(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon)) |
1466                   (1 << ((long)rparent))), "<, > or <;> expected",
1467          LINK->LINK);
1468  if (*LINK->LINK->sym == comma) GetSym_(LINK->LINK);
1469}
1470
1471
1472static void GetBlock(struct LOC_Lat_ProcessBlockInput *LINK)
1473{
1474  struct  LOC_GetBlock V;
1475  long    i, ii, k1, k4 = 0;
1476  long    SET[(long)invsym / 32 + 2];
1477  symset  SET1;
1478
1479  V.LINK = LINK;
1480  do {   /*7*/
1481    P_addset(P_expset(SET, 0), (long)ident);
1482    P_addset(SET, (long)intcon);
1483    P_addset(SET, (long)lparent);
1484    test_(P_addset(SET, (long)invsym),
1485          "<Element/Block name>, <integer>, <INV> or <(> expected",
1486          LINK);
1487    if (*LINK->sym == lparent)   /*7*/
1488      InsideParent(1, &V);
1489    else {
1490      if (*LINK->sym == invsym)
1491        Doinverse(&V);
1492      else {
1493        if (*LINK->sym == ident) {  /*8*/
1494          i = CheckElementtable(LINK->id, LINK->LINK);
1495          if (i != 0) {  /*9*/
1496            LINK->LINK->Bpointer++;
1497            if (LINK->LINK->Bpointer >= NoBEmax) {
1498              printf("** GetBlock: NoBEmax exceeded %ld (%d)\n",
1499                     LINK->LINK->Bpointer, NoBEmax);
1500              exit(1);
1501            }
1502            LINK->LINK->Bstack[LINK->LINK->Bpointer - 1] = i;
1503            GetSym_(LINK);
1504            if (*LINK->sym == comma)
1505              GetSym_(LINK);
1506          } else {  /*9*/
1507            ii = CheckBLOCKStable(LINK->id, LINK->LINK);
1508            if (ii != 0) {  /*10*/
1509              DeBlock(ii, 1, LINK);
1510            } else {  /*10*/
1511              ii = CheckUDItable(LINK->id, LINK->LINK);
1512              if (ii != 0) {  /*11*/
1513                k4 = (long)floor(LINK->LINK->UDItable[ii - 1].Uvalue + 0.5);
1514                GetSym_(LINK);
1515              } else
1516                test_(P_expset(SET1, 0), "invalid identifier",
1517                      LINK);
1518              /*11*/
1519              test_(P_expset(SET1, 1 << ((long)times)), "<*> expected",
1520                    LINK);
1521              if (*LINK->sym == times) {  /*11*/
1522                P_addset(P_expset(SET, 0), (long)ident);
1523                P_addset(SET, (long)lparent);
1524                getest_(P_addset(SET, (long)invsym),
1525                        "<element/Block name>, <INV> or <(> expected",
1526                        LINK);
1527                if (*LINK->sym == lparent)
1528                  InsideParent(k4, &V);
1529                else {
1530                  if (*LINK->sym == invsym)
1531                    Doinverse(&V);
1532                  else {
1533                    if (*LINK->sym == ident) {  /*12*/
1534                      i = CheckElementtable(LINK->id, LINK->LINK);
1535                      if (i != 0) {  /*13*/
1536                        for (k1 = 1; k1 <= k4; k1++) {  /*14*/
1537                          LINK->LINK->Bpointer++;
1538                          if (LINK->LINK->Bpointer >= NoBEmax) {
1539                            printf("** GetBlock: NoBEmax exceeded %ld (%d)\n",
1540                                   LINK->LINK->Bpointer, NoBEmax);
1541                            exit(1);
1542                          }
1543                          LINK->LINK->Bstack[LINK->LINK->Bpointer - 1] = i;
1544                        }  /*14*/
1545                        GetSym_(LINK);
1546                        if (*LINK->sym == comma)
1547                          GetSym_(LINK);
1548                      } else {  /*13*/
1549                        ii = CheckBLOCKStable(LINK->id, LINK->LINK);
1550                        if (ii == 0)
1551                          test_(P_expset(SET1, 0), "invalid name",
1552                                LINK);
1553                        DeBlock(ii, k4, LINK);
1554                      }  /*13*/
1555                    }  /*12*/
1556                  }
1557                }
1558              }  /*11*/
1559            }  /*10*/
1560          }  /*9*/
1561        } else {
1562          if (*LINK->sym == intcon) {  /*8*/
1563            k4 = *LINK->inum;
1564            GetSym_(LINK);
1565            test_(P_expset(SET1, 1 << ((long)times)), "<*> expected",
1566                  LINK);
1567            if (*LINK->sym == times) {  /*9*/
1568              GetSym_(LINK);
1569              P_addset(P_expset(SET, 0), (long)ident);
1570              P_addset(SET, (long)lparent);
1571              test_(P_addset(SET, (long)invsym),
1572                    "<element/Block name>, <INV> or <(> expected",
1573                    LINK);
1574              if (*LINK->sym == lparent)
1575                InsideParent(k4, &V);
1576              else {
1577                if (*LINK->sym == invsym)
1578                  Doinverse(&V);
1579                else {
1580                  if (*LINK->sym == ident) {  /*10*/
1581                    i = CheckElementtable(LINK->id, LINK->LINK);
1582                    if (i != 0) {  /*11*/
1583                      for (k1 = 1; k1 <= k4; k1++) {  /*12*/
1584                        LINK->LINK->Bpointer++;
1585                        if (LINK->LINK->Bpointer >= NoBEmax) {
1586                          printf("** GetBlock: NoBEmax exceeded %ld (%d)\n",
1587                                 LINK->LINK->Bpointer, NoBEmax);
1588                          exit(1);
1589                        }
1590                        LINK->LINK->Bstack[LINK->LINK->Bpointer - 1] = i;
1591                      }  /*12*/
1592                      GetSym_(LINK);
1593                      if (*LINK->sym == comma)
1594                        GetSym_(LINK);
1595                    } else {  /*11*/
1596                      ii = CheckBLOCKStable(LINK->id, LINK->LINK);
1597                      if (ii == 0)
1598                        test_(P_expset(SET1, 0), "invalid name",
1599                              LINK);
1600                      DeBlock(ii, k4, LINK);
1601                    }  /*11*/
1602                  }  /*10*/
1603                }
1604              }
1605            }  /*9*/
1606          } else {
1607            if (*LINK->sym == minus_) {  /*8*/
1608              GetSym_(LINK);
1609              i = CheckElementtable(LINK->id, LINK->LINK);
1610              if (i != 0) {  /*9*/
1611                LINK->LINK->Bpointer++;
1612                if (LINK->LINK->Bpointer >= NoBEmax) {
1613                  printf("** GetBlock: NoBEmax exceeded %ld (%d)\n",
1614                         LINK->LINK->Bpointer, NoBEmax);
1615                  exit(1);
1616                }
1617                LINK->LINK->Bstack[LINK->LINK->Bpointer - 1] = -i;
1618                GetSym_(LINK);
1619                if (*LINK->sym == comma)
1620                  GetSym_(LINK);
1621              } else
1622                test_(P_expset(SET1, 0), "<element name> expected.",
1623                      LINK);
1624            }  /*8*/
1625          }
1626        }
1627      }
1628    }
1629  } while (*LINK->sym == (long)invsym || *LINK->sym == (long)minus_ ||
1630           *LINK->sym == (long)intcon || *LINK->sym == (long)ident);
1631}
1632
1633
1634static void Lat_ProcessBlockInput(FILE **fi_, FILE **fo_, long *cc_, long *ll_,
1635                                  long *errpos_, long *lc_, long *nkw_,
1636                                  long *inum_, long emax__, long emin__,
1637                                  long kmax__, long nmax__, char *chin_,
1638                                  char *id_, char *BlockName,
1639                                  double *rnum_, bool *skipflag_, bool *rsvwd_,
1640                                  char *line_,
1641                                  Lat_symbol *sym_, alfa_ *key_,
1642                                  Lat_symbol *ksy_, Lat_symbol *sps_,
1643                                  struct LOC_Lattice_Read *LINK)
1644{
1645  struct LOC_Lat_ProcessBlockInput V;
1646  long i;
1647  symset SET;
1648  _REC_BlockStype *WITH;
1649
1650  V.LINK = LINK;
1651  V.fi = fi_;
1652  V.fo = fo_;
1653  V.cc = cc_;
1654  V.ll = ll_;
1655  V.errpos = errpos_;
1656  V.lc = lc_;
1657  V.nkw = nkw_;
1658  V.inum = inum_;
1659  V.emax_ = emax__;
1660  V.emin_ = emin__;
1661  V.kmax_ = kmax__;
1662  V.nmax_ = nmax__;
1663  V.chin = chin_;
1664  V.id = id_;
1665  V.rnum = rnum_;
1666  V.skipflag = skipflag_;
1667  V.rsvwd = rsvwd_;
1668  V.line = line_;
1669  V.sym = sym_;
1670  V.key = key_;
1671  V.ksy = ksy_;
1672  V.sps = sps_;
1673  i = CheckElementtable(BlockName, LINK);
1674  if (i != 0) {
1675    test_(P_expset(SET, 0), "<Block name>: conflict with Element name", &V);
1676    return;
1677  }
1678  /* Increment number of defined blocks */
1679  LINK->NoB++;
1680  if (LINK->NoB > NoBmax) {
1681    printf("** NoBmax exhausted: %ld(%d)\n", LINK->NoB, NoBmax);
1682    return;
1683  }
1684  WITH = &LINK->BlockS[LINK->NoB - 1];
1685  memcpy(WITH->Bname, BlockName, sizeof(partsName));
1686  WITH->BSTART = LINK->Bpointer + 1;
1687  GetBlock(&V);
1688  test_(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", &V);
1689  GetSym_(&V);
1690  WITH->BOWARI = LINK->Bpointer;
1691}  /* ProcessBlockInput */
1692
1693
1694static bool Lat_CheckWiggler(FILE **fo, long i, struct LOC_Lattice_Read *LINK)
1695{
1696  bool         Result;
1697  double       a, Lambda, L, diff;
1698  long         NN;
1699  ElemFamType  *WITH;
1700  elemtype     *WITH1;
1701  WigglerType  *WITH2;
1702
1703  Result = false;
1704  WITH = &ElemFam[i-1]; WITH1 = &WITH->ElemF; WITH2 = WITH1->W;
1705  Lambda = WITH2->lambda;
1706  L = WITH1->PL; a = L/Lambda;
1707  NN = (long)floor(a+0.01+0.5);
1708  diff = fabs((L-NN*Lambda)/L);
1709  if (diff < 1e-5) return true;
1710  printf("\n");
1711  printf(">>> Incorrect definition of %.*s\n\n", NameLength, WITH1->PName);
1712  printf("    L      ( total length ) =%20.12f [m]\n", L);
1713  printf("    Lambda ( wave  length ) =%20.12f [m]\n", Lambda);
1714  printf("    # of Period = L/Lambda  =%20.12f ?????\n\n", L / Lambda);
1715  return true;
1716}
1717
1718/* Local variables for Lat_DealElement: */
1719struct LOC_Lat_DealElement
1720{
1721  struct      LOC_Lattice_Read *LINK;
1722  FILE        **fi, **fo;
1723  long        *cc, *ll, *errpos, *lc, *nkw, *inum, emax_, emin_, kmax_, nmax_;
1724  char        *chin;
1725  char        *id;
1726  char        *BlockName;
1727  double      *rnum;
1728  bool        *skipflag, *rsvwd;
1729  char        *line;
1730  Lat_symbol  *sym;
1731  alfa_       *key;
1732  Lat_symbol  *ksy;
1733  Lat_symbol  *sps;
1734  jmp_buf     _JL9999;
1735  double      B[HOMmax+HOMmax+1];
1736  bool        BA[HOMmax+HOMmax+1];
1737  int         n_harm, harm[n_harm_max];
1738  double      kxV[n_harm_max], BoBrhoV[n_harm_max];
1739  double      kxH[n_harm_max], BoBrhoH[n_harm_max], phi[n_harm_max];
1740  long        DBNsavemax;
1741  DBNameType  DBNsave[nKidMax];
1742};
1743
1744
1745static void errorm__(const char *cmnt, struct LOC_Lat_DealElement *LINK)
1746{
1747  Lat_errorm(cmnt, LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos,
1748             LINK->lc, LINK->chin, LINK->skipflag, LINK->line, LINK->LINK);
1749}
1750
1751
1752static void GetSym__(struct LOC_Lat_DealElement *LINK)
1753{
1754  Lat_GetSym(LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos, LINK->lc,
1755             LINK->nkw, LINK->inum, LINK->emax_, LINK->emin_, LINK->kmax_,
1756             LINK->nmax_, LINK->chin, LINK->id, LINK->rnum, LINK->skipflag,
1757             LINK->rsvwd, LINK->line, LINK->sym, LINK->key, LINK->ksy,
1758             LINK->sps, LINK->LINK);
1759}
1760
1761static void test__(long *s1, const char *cmnt, struct LOC_Lat_DealElement *LINK)
1762{
1763  /*test*/
1764  if (!P_inset(*LINK->sym, s1))
1765    errorm__(cmnt, LINK);
1766}
1767
1768
1769static void getest__(long *s1, const char *cmnt, struct LOC_Lat_DealElement *LINK)
1770{
1771  /*test*/
1772  GetSym__(LINK);
1773  if (!P_inset(*LINK->sym, s1))
1774    errorm__(cmnt, LINK);
1775}
1776
1777
1778static double EVAL_(struct LOC_Lat_DealElement *LINK)
1779{
1780  return (Lat_EVAL(LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos,
1781                   LINK->lc, LINK->nkw, LINK->inum, LINK->emax_, LINK->emin_,
1782                   LINK->kmax_, LINK->nmax_, LINK->chin, LINK->id, LINK->rnum,
1783                   LINK->skipflag, LINK->rsvwd, LINK->line, LINK->sym,
1784                   LINK->key, LINK->ksy, LINK->sps, LINK->LINK));
1785}
1786
1787static void ProcessBlockInput(struct LOC_Lat_DealElement *LINK)
1788{
1789  Lat_ProcessBlockInput(LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos,
1790                        LINK->lc, LINK->nkw, LINK->inum, LINK->emax_,
1791                        LINK->emin_, LINK->kmax_, LINK->nmax_, LINK->chin,
1792                        LINK->id, LINK->BlockName, LINK->rnum, LINK->skipflag,
1793                        LINK->rsvwd, LINK->line, LINK->sym, LINK->key,
1794                        LINK->ksy, LINK->sps, LINK->LINK);
1795}  /* ProcessBlockInput */
1796
1797
1798static void CheckWiggler(long i, struct LOC_Lat_DealElement *LINK)
1799{
1800  if (!Lat_CheckWiggler(LINK->fo, i, LINK->LINK))
1801    longjmp(LINK->_JL9999, 1);
1802}
1803
1804static void GetHOM(struct LOC_Lat_DealElement *LINK)
1805{
1806  long    i;
1807  double  x, y;
1808  symset  SET;
1809
1810  getest__(P_expset(SET, 1 << ((long)lparent)), "<(> expected", LINK);
1811  do {
1812    i = (long)floor(EVAL_(LINK) + 0.5);
1813    if (i < 1 || HOMmax < i)
1814      getest__(P_expset(SET, 0), "invalid value detected", LINK);
1815    test__(P_expset(SET, 1 << ((long)comma)), "<, > expected", LINK);
1816    x = EVAL_(LINK);
1817    test__(P_expset(SET, 1 << ((long)comma)), "<, > expected", LINK);
1818    y = EVAL_(LINK);
1819    LINK->B[i+HOMmax] = x; LINK->B[HOMmax-i] = y;
1820    LINK->BA[i+HOMmax] = true; LINK->BA[HOMmax-i] = true;
1821    test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)rparent))),
1822           "<, > or <)> expected", LINK);
1823  } while (*LINK->sym != rparent);
1824  GetSym__(LINK);
1825}
1826
1827
1828static void ClearHOMandDBN(struct LOC_Lat_DealElement *LINK)
1829{
1830  long i;
1831
1832  for (i = -HOMmax; i <= HOMmax; i++) {
1833    LINK->B[i + HOMmax] = 0.0;
1834    LINK->BA[i + HOMmax] = false;
1835  }
1836  LINK->DBNsavemax = 0;
1837}
1838
1839
1840static void AssignHOM(long elem, struct LOC_Lat_DealElement *LINK)
1841{
1842  long       i;
1843  MpoleType  *M;
1844
1845  M = ElemFam[elem-1].ElemF.M;
1846  for (i = -HOMmax; i <= HOMmax; i++) {
1847    if (LINK->BA[i+HOMmax]) {
1848      M->PBpar[i+HOMmax] = LINK->B[i+HOMmax];
1849      M->Porder = max(abs(i), M->Porder);
1850    }
1851  }
1852}
1853
1854/*********************
1855 get the high order field
1856**********************/
1857static void GetHarm(struct LOC_Lat_DealElement *LINK)
1858{
1859  long    i, n;
1860  symset  SET;
1861
1862  getest__(P_expset(SET, 1 << ((long)lparent)), "<(> expected", LINK);
1863  LINK->n_harm = 0;
1864  do {
1865    LINK->n_harm++; n = LINK->n_harm;
1866    i = (long)floor(EVAL_(LINK)+0.5);
1867    if (i < 1 || n_harm_max < LINK->n_harm+1)
1868      getest__(P_expset(SET, 0), "invalid value detected", LINK);
1869    LINK->harm[n-1] = i;
1870    test__(P_expset(SET, 1 << ((long)comma)), "<, > expected", LINK);
1871    LINK->kxV[n-1] = EVAL_(LINK);
1872    test__(P_expset(SET, 1 << ((long)comma)), "<, > expected", LINK);
1873    LINK->BoBrhoV[n-1] = EVAL_(LINK);
1874    test__(P_expset(SET, 1 << ((long)comma)), "<, > expected", LINK);
1875    LINK->kxH[n-1] = EVAL_(LINK);
1876    test__(P_expset(SET, 1 << ((long)comma)), "<, > expected", LINK);
1877    LINK->BoBrhoH[n-1] = EVAL_(LINK);
1878    test__(P_expset(SET, 1 << ((long)comma)), "<, > expected", LINK);
1879    LINK->phi[n-1] = EVAL_(LINK);
1880    test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)rparent))),
1881           "<, > or <)> expected", LINK);
1882  } while (*LINK->sym != rparent);
1883  GetSym__(LINK);
1884}
1885
1886
1887static void AssignHarm(long elem, struct LOC_Lat_DealElement *LINK)
1888{
1889  long         i;
1890  WigglerType  *W;
1891
1892  W = ElemFam[elem-1].ElemF.W; W->n_harm += LINK->n_harm;
1893  // the fundamental is stored in harm[0], etc.
1894  for (i = 1; i < W->n_harm; i++) {
1895    W->harm[i] = LINK->harm[i-1];
1896    W->kxV[i] = LINK->kxV[i-1]; W->BoBrhoV[i] = LINK->BoBrhoV[i-1];
1897    W->kxH[i] = LINK->kxH[i-1]; W->BoBrhoH[i] = LINK->BoBrhoH[i-1];
1898    W->phi[i] = LINK->phi[i-1];
1899  }
1900}
1901
1902
1903static void GetDBN_(struct LOC_Lat_DealElement *LINK)
1904
1905{
1906  long i;
1907  symset SET;
1908  long SET1[(long)squote / 32 + 2];
1909
1910  getest__(P_expset(SET, 1 << ((long)lparent)), "<(> expected:GetDBN", LINK);
1911  do {
1912    getest__(P_addset(P_expset(SET1, 0), (long)squote),
1913             "<'> expected:GetDBN", LINK);
1914    LINK->DBNsavemax++;
1915    for (i = 0; i < DBNameLen; i++)
1916      LINK->DBNsave[LINK->DBNsavemax - 1][i] = ' ';
1917    i = 0;
1918    while (*LINK->chin != '\'' && i < DBNameLen) {
1919      i++;
1920      LINK->DBNsave[LINK->DBNsavemax - 1][i - 1] = *LINK->chin;
1921      Lat_Nextch(LINK->fi, LINK->fo, LINK->cc, LINK->ll, LINK->errpos,
1922                 LINK->lc, LINK->chin, LINK->skipflag, LINK->line,
1923                 LINK->LINK);
1924    }
1925    getest__(P_addset(P_expset(SET1, 0), (long)squote),
1926             "<'> expected:GetDBN", LINK);
1927    getest__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)rparent))),
1928             "<, > or <)> expected:GetDBN", LINK);
1929  }
1930  while (*LINK->sym != rparent);
1931  GetSym__(LINK);
1932}
1933
1934
1935static void adjdbname(char *DBname1, char *DBname2)
1936{
1937  long i, j, k, offset;
1938  bool first, blank;
1939
1940  blank = true;
1941  for (j = 0; j < DBNameLen; j++) {
1942    if (DBname1[j] != ' ')
1943      blank = false;
1944  }
1945  first = true;
1946  j = 0;
1947  offset = 0;
1948  if (blank)
1949    return;
1950  do {
1951    j++;
1952    DBname2[j + offset - 1] = DBname1[j - 1];
1953    if (first && DBname1[j - 1] == ' ') {
1954      first = false;
1955      k = -1;
1956      do {
1957        k++;
1958      } while (DBname1[j+k] == ' ' && j+k+1 != DBNameLen);
1959      for (i = j+k; i <= 7; i++)
1960        DBname2[i] = ' ';
1961      offset = 8 - j;
1962    }
1963  }
1964  while (j < DBNameLen-offset);
1965}
1966
1967
1968static void SetDBN(struct LOC_Lat_DealElement *LINK)
1969{
1970  long         i;
1971  ElemFamType  *WITH;
1972  long         FORLIM;
1973
1974  if (globval.Elem_nFam > Elem_nFamMax) {
1975    printf("Elem_nFamMax exceeded: %ld(%ld)\n",
1976           globval.Elem_nFam, (long)Elem_nFamMax);
1977    exit_(1);
1978  }
1979  WITH = &ElemFam[globval.Elem_nFam-1]; WITH->NoDBN = LINK->DBNsavemax;
1980  if (WITH->NoDBN > 0) {
1981    FORLIM = WITH->NoDBN;
1982    for (i = 0; i < FORLIM; i++)
1983      adjdbname(LINK->DBNsave[i], WITH->DBNlist[i]);
1984  }
1985}
1986
1987/**************************************************************************
1988 static bool Lat_DealElement(FILE **fi_, FILE **fo_, long *cc_, long *ll_,
1989                          long *errpos_, long *lc_, long *nkw_, long *inum_,
1990                          long emax__, long emin__,
1991                          long kmax__, long nmax__, char *chin_, char *id_,
1992                          char *ElementName,
1993                          char *BlockName_, double *rnum_, bool *skipflag_,
1994                          bool *rsvwd_,
1995                          char *line_, Lat_symbol *sym_, alfa_ *key_,
1996                          Lat_symbol *ksy_,
1997                          Lat_symbol *sps_, struct LOC_Lattice_Read *LINK)
1998
1999   Purpose:
2000             Read the lattice from file fi_,
2001              first, read the element parameters,
2002              then, read the lattice
2003              element type: drift, quadrupole,sextupole,bending,cavity,.....
2004
2005    *************************************************************************/
2006static bool Lat_DealElement(FILE **fi_, FILE **fo_, long *cc_, long *ll_,
2007                            long *errpos_, long *lc_, long *nkw_, long *inum_,
2008                            long emax__, long emin__,
2009                            long kmax__, long nmax__, char *chin_, char *id_,
2010                            char *ElementName,
2011                            char *BlockName_, double *rnum_, bool *skipflag_,
2012                            bool *rsvwd_,
2013                            char *line_, Lat_symbol *sym_, alfa_ *key_,
2014                            Lat_symbol *ksy_,
2015                            Lat_symbol *sps_, struct LOC_Lattice_Read *LINK)
2016{
2017  struct LOC_Lat_DealElement V;
2018  bool           Result = false;
2019  double         t, t1, t2, gap, QL = 0.0, QK;
2020  double         QKick;  /* kick angle of the corrector [rad]*/
2021  int            Kplane; /* kick plane of the corrector, 1 for H plane, -1 for V plane*/
2022  double         QKV, QKH, QKxV, QKxH, QPhi, QKS;
2023  double         dt, Frf, Vrf;
2024  long k1, k2, harnum, FF1, FF2;
2025  double FFscaling;
2026  Lat_symbol     sym1;
2027  symset         mysys, SET;
2028  long           SET1[(long)lsym / 32 + 2];
2029  ElemFamType    *WITH;
2030  elemtype       *WITH1;
2031  MpoleType      *WITH2;
2032  symset         SET2;
2033  CavityType     *WITH3;
2034  long           SET3[(long)possym / 32 + 2];
2035  long           SET4[(long)dbnsym / 32 + 2];
2036  WigglerType    *WITH4;
2037  FieldMapType   *WITH6;
2038  InsertionType  *WITH5;   /* ID Laurent */
2039  SolenoidType   *WITH7;
2040  char str1[100] = "";
2041  char str2[100] = "";
2042  bool firstflag  = false; // flag for first kick input
2043  bool secondflag = false; // flag for second kick input
2044  long           i;
2045  double         scaling1, scaling2;
2046
2047  V.LINK = LINK;
2048  V.fi = fi_; V.fo = fo_; V.cc = cc_; V.ll = ll_; V.errpos = errpos_;
2049  V.lc = lc_; V.nkw = nkw_; V.inum = inum_; V.emax_ = emax__; V.emin_ = emin__;
2050  V.kmax_ = kmax__; V.nmax_ = nmax__; V.chin = chin_; V.id = id_;
2051  V.BlockName = BlockName_; V.rnum = rnum_; V.skipflag = skipflag_;
2052  V.rsvwd = rsvwd_; V.line = line_; V.sym = sym_; V.key = key_; V.ksy = ksy_;
2053  V.sps = sps_;
2054  if (setjmp(V._JL9999)) goto _L9999;
2055  Result = false;
2056
2057  switch (*V.sym) {
2058   
2059    /**************************************************************************
2060       Drift
2061    **************************************************************************
2062
2063    <name>: Drift,
2064            L=<length>; [m]
2065
2066    Example
2067
2068      L1: Drift, L=0.30;
2069
2070    *************************************************************************/
2071  case drfsym:
2072    getest__(P_expset(SET, 1 << ((long)comma)), "<comma> expected", &V);
2073    getest__(P_addset(P_expset(SET1, 0), (long)lsym), "<L> expected", &V);
2074    getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2075    *V.rnum = EVAL_(&V);
2076    test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", &V);
2077    GetSym__(&V);
2078    globval.Elem_nFam++;
2079    if (globval.Elem_nFam <= Elem_nFamMax) {
2080      WITH = &ElemFam[globval.Elem_nFam-1];
2081      WITH1 = &WITH->ElemF;
2082      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2083      WITH1->PL = *V.rnum;
2084      WITH1->Pkind = PartsKind(drift);
2085      Drift_Alloc(&WITH->ElemF);
2086    } else {
2087      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2088             globval.Elem_nFam, (long)Elem_nFamMax);
2089      exit_(1);
2090    }
2091    break;
2092
2093    /**************************************************************************
2094      bending
2095    **************************************************************************
2096
2097    <name>: Bending,
2098            L      = <length>, ( [m] )
2099            T      = <bending angle>, ( [degree] )
2100            T1     = <entrance angle>, ( [degree] )
2101            T2     = <exit angle>, ( [degree] )
2102            gap    = <total magnet gap>, ( [m] )
2103            K      = <K-value>, ( [m^-2] )
2104                     ( K > 0 : focusing in horizontal plane )
2105                     ( K < 0 : defocusing in horizontal plane )
2106            N      = <# of kicks>,
2107            method = <method>, ( 2 or 4. The method to divide Q into slices.)
2108                           ( The detail of <method> will be discussed later.)
2109            Default value is 2.
2110            Roll   = <roll angle>, ( [deg], design roll angle )
2111            HOM    = (i, <Bi>, <Ai>, ( multipole compoments (power series) )
2112                      j, <Bj>, <Aj>, ( Systematic error Only )
2113                      ............   ( Random errors are assigned )
2114                     n, <Bn>, <An>); ( in a Program File using procedures )
2115
2116    Example
2117
2118      B: bending, L=0.70, T=10.0, T1:=5.0, T2:=5.0, K=-1.0, N=8, Method=2;
2119
2120    *************************************************************************/
2121  case bndsym:  /*4*/
2122    getest__(P_expset(SET, 1 << ((long)comma)), "<, > expected", &V);
2123    GetSym__(&V);
2124    QL = 0.0;   /* L */
2125    QK = 0.0;   /* K, quadrupole components*/
2126    k1 = 0;   /* N */
2127    t  = 0.0;   /* T */
2128    t1 = 0.0;   /* T1 */
2129    t2 = 0.0;   /* T2 */
2130    gap = 0.0;   /* gap */
2131    k2 = Meth_Linear;   /* method */
2132    dt = 0.0;
2133    ClearHOMandDBN(&V);
2134    P_addset(P_expset(mysys, 0), (long)lsym);
2135    P_addset(mysys, (long)ksym);
2136    P_addset(mysys, (long)nsym);
2137    P_addset(mysys, (long)mthsym);
2138    P_addset(mysys, (long)tsym);
2139    P_addset(mysys, (long)rollsym);
2140    P_addset(mysys, (long)t1sym);
2141    P_addset(mysys, (long)t2sym);
2142    P_addset(mysys, (long)gapsym);
2143    P_addset(mysys, (long)homsym);
2144    P_addset(mysys, (long)dbnsym);
2145    do {
2146      test__(mysys, "illegal parameter", &V);
2147      sym1 = *V.sym;
2148      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2149
2150      switch (sym1) {
2151
2152      case lsym:
2153        QL = EVAL_(&V);
2154        break;
2155
2156      case ksym:
2157        QK = EVAL_(&V);
2158        break;
2159
2160      case nsym:
2161        k1 = (long)floor(EVAL_(&V) + 0.5);
2162        break;
2163
2164      case tsym:
2165        t = EVAL_(&V);
2166        break;
2167
2168      case rollsym:
2169        dt = EVAL_(&V);
2170        break;
2171
2172      case t1sym:
2173        t1 = EVAL_(&V);
2174        break;
2175
2176      case t2sym:
2177        t2 = EVAL_(&V);
2178        break;
2179
2180      case gapsym:
2181        gap = EVAL_(&V);
2182        break;
2183
2184      case mthsym:
2185        k2 = (long)floor(EVAL_(&V) + 0.5);
2186        if (k2 != Meth_Linear) globval.MatMeth = false;
2187        if ((unsigned int)k2 >= 32 ||
2188            ((1 << k2) & ((1 << Meth_Linear) | (1 << Meth_Second) |
2189                          (1 << Meth_Fourth))) == 0)
2190          getest__(P_expset(SET, 0), "Check integrator..", &V);
2191        break;
2192
2193      case homsym:
2194        GetHOM(&V);
2195        break;
2196
2197      case dbnsym:
2198        GetDBN_(&V);
2199        break;
2200      default:
2201        break;
2202      }
2203      test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
2204             "<, > or <;> expected", &V);
2205      if (*V.sym == comma)
2206        GetSym__(&V);
2207    } while (P_inset(*V.sym, mysys));
2208    test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected.", &V);
2209    GetSym__(&V);
2210    globval.Elem_nFam++;
2211    if (globval.Elem_nFam <= Elem_nFamMax) {
2212      WITH = &ElemFam[globval.Elem_nFam-1];
2213      WITH1 = &WITH->ElemF;
2214      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2215      WITH1->PL = QL;
2216      WITH1->Pkind = Mpole;
2217      Mpole_Alloc(&WITH->ElemF);
2218      WITH2 = WITH1->M;
2219      WITH2->Pmethod = k2;
2220      WITH2->PN = k1;
2221      if (WITH1->PL != 0.0)
2222        WITH2->Pirho = t * M_PI / 180.0 / WITH1->PL;
2223      else
2224        WITH2->Pirho = t * M_PI / 180.0;
2225      WITH2->PTx1 = t1; WITH2->PTx2 = t2; WITH2->Pgap = gap;
2226      WITH2->n_design = Dip;
2227      AssignHOM(globval.Elem_nFam, &V);
2228      SetDBN(&V);
2229      WITH2->PBpar[HOMmax+2] = QK; WITH2->PdTpar = dt;
2230    } else {
2231      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2232             globval.Elem_nFam, (long)Elem_nFamMax);
2233      exit_(1);
2234    }
2235    break;
2236   
2237    /**************************************************************************
2238      Quadrupole
2239    **************************************************************************
2240
2241    <name>: Quadrupole,
2242            L=<length>, ( [m] )
2243            K =<K-value>, ( [m-2] )
2244            N =<# of kicks>,
2245            Method=<method>,
2246            FF1, FF1 = 1 means scaling factor entrance, exit
2247            FFscaling = Fringe field scaling factor, default 1.0,
2248            Roll or tilt =<roll angle>, ( [deg], design roll angle )
2249            HOM=(i, <Bi>, <Ai>, ( higher order component in USA notation )
2250                 j, <Bj>, <Aj>, ( Systematic error Only )
2251                 ............    ( Random errors are assigned )
2252                 n, <Bn>, <An>); ( in a Program File using procedures )
2253
2254    Example
2255
2256      Q1: Quadrupole, L=0.5, K=2.213455, N=4, Method=4;
2257      Q1 : Quadrupole, L=0.5, K=2.213455, N=4, FF1 =1, FF2, FFscaling=1, Method=4;
2258
2259    **************************************************************/
2260  case qdsym:  /*4*/
2261    getest__(P_expset(SET, 1 << ((long)comma)), "<, > expected", &V);
2262    GetSym__(&V);
2263    QL = 0.0;   /* L */
2264    QK = 0.0;   /* K */
2265    k1 = 0;   /* N */
2266    k2 = Meth_Fourth;   /* method */
2267    dt = 0.0;
2268    FF1 = 0;     /* Entrance Fringe field */
2269    FF2 = 0;     /* Exit Fringe field */
2270    FFscaling = 1.0; /*Fringe field scaling factor */
2271    ClearHOMandDBN(&V);
2272    P_addset(P_expset(mysys, 0), (long)lsym);
2273    P_addset(mysys, (long)ksym);
2274    P_addset(mysys, (long)nsym);
2275    P_addset(mysys, (long)mthsym);
2276    P_addset(mysys, (long)rollsym);
2277    P_addset(mysys, (long)tiltsym);
2278    P_addset(mysys, (long)homsym);
2279    P_addset(mysys, (long)dbnsym);
2280    P_addset(mysys, (long)ff1sym);
2281    P_addset(mysys, (long)ff2sym);
2282    P_addset(mysys, (long)ffscalingsym);
2283    do {   /*5: read L, K, N, T, T1, T2 */
2284      test__(mysys, "illegal parameter", &V);
2285      sym1 = *V.sym;
2286      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2287      switch (sym1) {   /*6*/
2288
2289      case lsym:
2290        QL = EVAL_(&V);
2291        break;
2292
2293      case ffscalingsym:     // read quad ff scaling factor, default 1.0
2294      FFscaling = EVAL_(&V);
2295      break;
2296
2297      case ff1sym:
2298           FF1 = (long)EVAL_(&V);
2299           break;
2300
2301      case ff2sym:
2302           FF2 = (long)EVAL_(&V);
2303           break;
2304
2305      case ksym:
2306        QK = EVAL_(&V);
2307        break;
2308
2309      case nsym:
2310        k1 = (long)floor(EVAL_(&V) + 0.5);
2311        break;
2312
2313      case mthsym:
2314        k2 = (long)floor(EVAL_(&V) + 0.5);
2315        if (k2 != Meth_Linear) globval.MatMeth = false;
2316        if ((unsigned int)k2 >= 32 ||
2317            ((1 << k2) & ((1 << Meth_Linear) | (1 << Meth_First) |
2318                          (1 << Meth_Second) | (1 << Meth_Fourth))) == 0)
2319          getest__(P_expset(SET, 0), "Check integrator..", &V);
2320        break;
2321
2322      case rollsym:
2323        dt = EVAL_(&V);
2324        break;
2325
2326      case tiltsym:
2327        dt = EVAL_(&V);
2328        break;
2329
2330      case homsym:
2331        GetHOM(&V);
2332        break;
2333
2334      case dbnsym:
2335        GetDBN_(&V);
2336        break;
2337      default:
2338        break;
2339      }
2340      test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
2341             "<, > or <;> expected", &V);
2342      if (*V.sym == comma)
2343        GetSym__(&V);
2344    } while (P_inset(*V.sym, mysys));   /*5*/
2345    test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected.", &V);
2346    GetSym__(&V);
2347    globval.Elem_nFam++;
2348    if (globval.Elem_nFam <= Elem_nFamMax) {
2349      WITH = &ElemFam[globval.Elem_nFam-1];
2350      WITH1 = &WITH->ElemF;
2351      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2352      WITH1->PL = QL; WITH1->Pkind = Mpole;
2353      Mpole_Alloc(&WITH->ElemF);
2354      WITH2 = WITH1->M;
2355      WITH2->Pmethod = k2; WITH2->PN = k1; WITH2->PdTpar = dt;
2356      WITH2->quadFF1 = FF1; /* entrance fringe field flag */
2357      WITH2->quadFF2 = FF2; /* exit fringe field flag */
2358      WITH2->quadFFscaling = FFscaling; /*fringe field scaling factor flag */
2359          WITH2->sextFF1 = 0; /* entrance sextupole fringe field flag */
2360       WITH2->sextFF2 = 0; /* exit sextupole fringe field flag */
2361
2362      AssignHOM(globval.Elem_nFam, &V);
2363      SetDBN(&V);
2364      WITH2->n_design = Quad; WITH2->PBpar[HOMmax+2] = QK;
2365    } else {
2366      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2367             globval.Elem_nFam, (long)Elem_nFamMax);
2368      exit_(1);
2369    }
2370    break;
2371
2372    /**************************************************************************
2373      Sextupole
2374    ***************************************************************************
2375
2376    <name>: Sextupole,
2377            L=<length>, ( [m] )
2378            K =<K-value>, ( [m-3] )
2379            Roll=<roll angle>, ( [degree], design roll angle )
2380            HOM=(i, <Bi>, <Ai>, ( higher order component in USA notation )
2381                 j, <Bj>, <Aj>, ( Systematic error Only )
2382                 ............    ( Random errors are assigned )
2383                 n, <Bn>, <An>); ( in a Program File using procedures )
2384
2385    Example
2386
2387      SF: Sextupole, K=-10.236345;
2388
2389    **************************************************************************/
2390  case sexsym:  /*4*/
2391    QL = 0.0;   /* L */
2392    QK = 0.0;   /* K */
2393    k1 = 0;   /* N */
2394    k2 = Meth_Fourth;   /* method */
2395    FF1 = 0;     /* Entrance Fringe field */
2396    FF2 = 0;     /* Exit Fringe field */
2397    dt = 0.0;
2398    ClearHOMandDBN(&V);
2399    getest__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
2400             "<, > or <;> expected", &V);
2401    if (*V.sym == comma) {
2402      GetSym__(&V);
2403      P_addset(P_expset(mysys, 0), (long)lsym);
2404      P_addset(mysys, (long)ksym);
2405      P_addset(mysys, (long)nsym);
2406      P_addset(mysys, (long)mthsym);
2407      P_addset(mysys, (long)rollsym);
2408      P_addset(mysys, (long)homsym);
2409      P_addset(mysys, (long)dbnsym);
2410      P_addset(mysys, (long)ff1sym);
2411      P_addset(mysys, (long)ff2sym);
2412      do {   /*5: read L, K, N, T, T1, T2 */
2413        test__(mysys, "illegal parameter", &V);
2414        sym1 = *V.sym;
2415        getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2416        switch (sym1)
2417          {   /*6*/
2418          case lsym:
2419            QL = EVAL_(&V);
2420            break;
2421
2422          case ksym:
2423            QK = EVAL_(&V);
2424            break;
2425
2426          case nsym:
2427            k1 = (long)floor(EVAL_(&V) + 0.5);
2428            break;
2429
2430          case ff1sym:
2431        FF1 = (long)EVAL_(&V);
2432        break;
2433
2434      case ff2sym:
2435        FF2 = (long)EVAL_(&V);
2436        break;
2437
2438          case mthsym:
2439            k2 = (long)floor(EVAL_(&V) + 0.5);
2440            if (k2 != Meth_Linear) globval.MatMeth = false;
2441            if ((unsigned int)k2 >= 32 ||
2442                ((1 << k2) & ((1 << Meth_Linear) | (1 << Meth_Second) |
2443                              (1 << Meth_Fourth))) == 0)
2444              getest__(P_expset(SET, 0), "Check integrator..", &V);
2445            break;
2446
2447          case rollsym:
2448            dt = EVAL_(&V);
2449            break;
2450
2451          case homsym:
2452            GetHOM(&V);
2453            break;
2454
2455          case dbnsym:
2456            GetDBN_(&V);
2457            break;
2458          default:
2459            break;
2460          }
2461
2462        test__(P_expset(SET,
2463                        (1 << ((long)comma)) | (1 << ((long)semicolon))),
2464               "<, > or <;> expected", &V);
2465
2466        if (*V.sym == comma)
2467          GetSym__(&V);
2468
2469      } while (P_inset(*V.sym, mysys));   /*5*/
2470      test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected.", &V);
2471    }
2472    GetSym__(&V);
2473    globval.Elem_nFam++;
2474    if (globval.Elem_nFam <= Elem_nFamMax) {
2475      WITH = &ElemFam[globval.Elem_nFam-1];
2476      WITH1 = &WITH->ElemF;
2477      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2478      WITH1->PL = QL;
2479      WITH1->Pkind = Mpole;
2480      Mpole_Alloc(&WITH->ElemF);
2481      WITH2 = WITH1->M;
2482      WITH2->Pmethod = k2;
2483      WITH2->PN = k1;
2484      WITH2->quadFF1 = 0; /* entrance fringe field flag */
2485      WITH2->quadFF2 = 0; /* exit fringe field flag */
2486      WITH2->sextFF1 = FF1; /* entrance fringe field flag */
2487      WITH2->sextFF2 = FF2; /* exit fringe field flag */
2488      if (WITH1->PL != 0.0)
2489        WITH2->Pthick = pthicktype(thick);
2490      else
2491        WITH2->Pthick = pthicktype(thin);
2492      WITH2->PdTpar = dt; WITH2->n_design = Sext;
2493      AssignHOM(globval.Elem_nFam, &V);
2494      SetDBN(&V);
2495      WITH2->PBpar[HOMmax + 3] = QK;
2496    } else {
2497      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2498             globval.Elem_nFam, (long)Elem_nFamMax);
2499      exit_(1);
2500    }
2501    break;
2502
2503    /**************************************************************************
2504      Cavity
2505    ***************************************************************************
2506
2507    <name>: Cavity,
2508            Frequency = <Frf>,   ( [Hz] )
2509            Voltage   = <Vrf>,   ( [V]  )
2510            Phase     = <phi_rf> (degrees)
2511            harnum    = <h>
2512
2513    Example
2514
2515      CAV: Cavity, Frequency = 499.95e6, Voltage=1.22e6, harnum=328;
2516
2517    **************************************************************************/
2518  case cavsym:
2519    ClearHOMandDBN(&V);
2520    getest__(P_expset(SET, 1 << ((long)comma)), "<, > expected", &V);
2521    GetSym__(&V);
2522    Frf = 0.0;   /* Frf */
2523    Vrf = 0.0;   /* Vrf */
2524    QPhi = 0.0;
2525    harnum = 0;   /* Voff */
2526    P_addset(P_expset(mysys, 0), (long)frqsym);
2527    P_addset(mysys, (long)vrfsym);
2528    P_addset(mysys, (long)phisym);
2529    P_addset(mysys, (long)harnumsym);
2530    P_addset(mysys, (long)dbnsym);
2531    do {
2532      test__(mysys, "illegal parameter", &V);
2533      sym1 = *V.sym;
2534      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2535      switch (sym1) {
2536
2537      case frqsym:
2538        Frf = EVAL_(&V);
2539        break;
2540
2541      case vrfsym:
2542        Vrf = EVAL_(&V);
2543        break;
2544
2545      case phisym:
2546        QPhi = EVAL_(&V);
2547        break;
2548
2549      case harnumsym:
2550        harnum = (long)floor(EVAL_(&V) + 0.5);
2551        break;
2552
2553      case dbnsym:
2554        GetDBN_(&V);
2555        break;
2556      default:
2557        break;
2558      }
2559      test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
2560             "<, > or <;> expected", &V);
2561      if (*V.sym == comma)
2562        GetSym__(&V);
2563    } while (P_inset(*V.sym, mysys));
2564    test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected.", &V);
2565    GetSym__(&V);
2566    globval.Elem_nFam++;
2567    if (globval.Elem_nFam <= Elem_nFamMax) {
2568      WITH = &ElemFam[globval.Elem_nFam-1];
2569      WITH1 = &WITH->ElemF;
2570      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2571      WITH1->Pkind = Cavity;
2572      Cav_Alloc(&WITH->ElemF);
2573      WITH3 = WITH1->C;
2574      WITH3->Pvolt = Vrf;   /* Voltage [V] */
2575      WITH3->Pfreq = Frf;   /* Frequency in Hz */
2576      WITH3->phi = QPhi*M_PI/180.0;
2577      WITH3->Ph = harnum;  /* RF harmonic number */
2578      SetDBN(&V);
2579    } else {
2580      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2581             globval.Elem_nFam, (long)Elem_nFamMax);
2582      exit_(1);
2583    }
2584    break;
2585
2586    /**************************************************************************
2587      Corrector
2588    ***************************************************************************
2589   
2590      Kickers specific for orbit correction.
2591   
2592    <name>: Corrector, <direction>, L=<length>, kick = <kick angle [rad]>;
2593
2594    <name> :== Alphanumeric string. Up to NameLength character length.
2595              BEGIN with an alphabet.
2596    <direction> :== 'horizontal'|'vertical'
2597
2598    Example
2599
2600      COH: Corrector, horizontal;
2601
2602    **************************************************************************/
2603
2604  case corsym:  /*4*/
2605    QL = 0.0;   /* L */
2606    QKick = 0.0; /* kick angle of the corrector [rad]*/
2607    Kplane = 0;   /* 1 is horizontal corrector, -1 is vertical corrector */
2608    k1 = 0;     /* N */
2609    k2 = Meth_Linear;   /* method */
2610    dt = 0.0;
2611    ClearHOMandDBN(&V);
2612    getest__(P_expset(SET, 1 << ((long)comma)), "<, > expected", &V);
2613    if (*V.sym == comma) {
2614      GetSym__(&V);
2615      P_addset(P_expset(mysys, 0), (long)lsym);
2616      P_addset(mysys, (long)nsym);
2617      P_addset(mysys, (long)mthsym);
2618      P_addset(mysys, (long)horsym);
2619      P_addset(mysys, (long)versym);
2620      P_addset(mysys, (long)corkicksym);
2621      P_addset(mysys, (long)rollsym);
2622      P_addset(mysys, (long)dbnsym);
2623      do {   /*5: read L, K, N, T, T1, T2 */
2624        test__(mysys, "illegal parameter", &V); sym1 = *V.sym;
2625        if (*V.sym == (long)dbnsym || *V.sym == (long)rollsym ||
2626            *V.sym == (long)mthsym || *V.sym == (long)nsym || 
2627            *V.sym == (long)lsym || *V.sym == (long)corkicksym)
2628          {
2629            getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2630            if (*V.sym == eql) {
2631              switch (sym1) {   /*6*/
2632
2633              case lsym:
2634                QL = EVAL_(&V);
2635                break;
2636
2637                 case corkicksym:
2638                QKick = EVAL_(&V);
2639                break; 
2640               
2641              case nsym:
2642                k1 = (long)floor(EVAL_(&V) + 0.5);
2643                break;
2644
2645              case mthsym:
2646                k2 = (long)floor(EVAL_(&V) + 0.5);
2647                if (k2 != Meth_Linear) globval.MatMeth = false;
2648                if ((unsigned int)k2 >= 32 ||
2649                    ((1 << k2) & ((1 << Meth_Linear) | (1 << Meth_Second) |
2650                                  (1 << Meth_Fourth))) == 0)
2651                  getest__(P_expset(SET2, 0), "Check integrator..", &V);
2652                break;
2653
2654              case rollsym:
2655                dt = EVAL_(&V);
2656                break;
2657              case dbnsym:
2658                GetDBN_(&V);
2659                break;
2660              default:
2661                break;
2662              }
2663            }
2664          } else {
2665            if (sym1 == horsym)
2666              Kplane = 1;
2667            else if (sym1 == versym)
2668              Kplane = -1;
2669            GetSym__(&V);
2670          }
2671        test__(P_expset(SET,
2672                        (1 << ((long)comma)) | (1 << ((long)semicolon))),
2673               "<, > or <;> expected", &V);
2674        if (*V.sym == comma)
2675          GetSym__(&V);
2676      } while (P_inset(*V.sym, mysys));   /*5*/
2677
2678      test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected.", &V);
2679    }
2680    GetSym__(&V);
2681    globval.Elem_nFam++;
2682    if (globval.Elem_nFam <= Elem_nFamMax) {
2683      WITH = &ElemFam[globval.Elem_nFam-1];
2684      WITH1 = &WITH->ElemF;
2685      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2686      WITH1->PL = QL;
2687      WITH1->Pkind = Mpole;
2688      Mpole_Alloc(&WITH->ElemF);
2689      WITH2 = WITH1->M;
2690      SetDBN(&V);
2691      if (WITH1->PL != 0.0)
2692        WITH2->Pthick = pthicktype(thick);
2693      else
2694        WITH2->Pthick = pthicktype(thin);
2695      WITH2->Pmethod = k2;
2696      WITH2->PN = k1;
2697      WITH2->PdTpar = dt;
2698
2699        if(Kplane == 0){
2700                cout << "t2lat: Error! Must specify the type of the corrector, Horizontal or vertical!" << endl;
2701                exit_(1);
2702            }
2703      WITH2->PBpar[Kplane*Dip + HOMmax] = -1*QKick;  //assign the kick angle [rad]
2704 
2705   
2706   
2707    } else {
2708      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2709             globval.Elem_nFam, (long)Elem_nFamMax);
2710      exit_(1);
2711    }
2712    break;
2713
2714    /**************************************************************************
2715      Beam Position Monitor
2716    ***************************************************************************
2717
2718    <name>: Beam Position Monitor;
2719
2720    <name>:== Alphanumeric string. Up to NameLength character length.
2721              BEGIN with an alphabet.
2722
2723    Example
2724
2725      BPM: Beam Position Monitor;
2726
2727    **************************************************************************/
2728
2729  case bemsym:
2730    ClearHOMandDBN(&V);
2731    getest__(P_addset(P_expset(SET3, 0), (long)possym),
2732             "<position> expected", &V);
2733    getest__(P_expset(SET, 1 << ((long)monsym)), "<monitor> expected", &V);
2734    getest__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
2735             "<, > or <;> expected", &V);
2736    if (*V.sym == comma) {
2737      getest__(P_addset(P_expset(SET4, 0), (long)dbnsym),
2738               "illegal parameter", &V);
2739      sym1 = *V.sym;
2740      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2741      if (sym1 == dbnsym)
2742        GetDBN_(&V);
2743      test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", &V);
2744    }
2745    GetSym__(&V);
2746    globval.Elem_nFam++;
2747    if (globval.Elem_nFam <= Elem_nFamMax) {
2748      WITH = &ElemFam[globval.Elem_nFam-1];
2749      WITH1 = &WITH->ElemF;
2750      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2751      WITH1->Pkind = Mpole;
2752      Mpole_Alloc(&WITH->ElemF);
2753      WITH2 = WITH1->M;
2754      WITH2->Pthick = pthicktype(thin);
2755      SetDBN(&V);
2756    } else {
2757      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2758             globval.Elem_nFam, (long)Elem_nFamMax);
2759      exit_(1);
2760    }
2761    break;
2762
2763
2764    /**************************************************************************
2765      Marker
2766    ***************************************************************************
2767
2768    <name>: Marker;
2769
2770    <name>:== Alphanumeric string. Up to NameLength character length.
2771              BEGIN with an alphabet.
2772
2773    Example
2774
2775      SYM: Marker;
2776
2777    **************************************************************************/
2778
2779  case mrksym:
2780    ClearHOMandDBN(&V);
2781    getest__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
2782             "<, > or <;> expected", &V);
2783    if (*V.sym == comma) {
2784      getest__(P_addset(P_expset(SET4, 0), (long)dbnsym),
2785               "illegal parameter", &V);
2786      sym1 = *V.sym;
2787      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2788      if (sym1 == dbnsym)
2789        GetDBN_(&V);
2790      test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", &V);
2791    }
2792    GetSym__(&V);
2793    globval.Elem_nFam++;
2794    if (globval.Elem_nFam <= Elem_nFamMax) {
2795      WITH = &ElemFam[globval.Elem_nFam-1];
2796      WITH1 = &WITH->ElemF;
2797      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2798      WITH1->PL = 0.0;
2799      WITH1->Pkind = PartsKind(marker);
2800      SetDBN(&V);
2801    } else {
2802      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2803             globval.Elem_nFam, (long)Elem_nFamMax);
2804      exit_(1);
2805    }
2806    break;
2807
2808   
2809    /**************************************************************************
2810      Ghost
2811    ***************************************************************************
2812
2813     <name>: Ghost;
2814
2815     <name>:== Alphanumeric string. Up to NameLength character length.
2816               BEGIN with an alphabet.
2817
2818     Example
2819
2820       OBAKE : Ghost;
2821
2822    **************************************************************************/
2823    /*----------->>>
2824      GstSym:BEGIN
2825      getest([comma], '<, > expexted');
2826      getest([typsym], '<type> expected');
2827      getest([eql], '<=> expected');
2828      QL:=Eval;
2829      test([semicolon], '<;> expected');
2830      getsym;
2831      if sym=DBNsym then GetDBN;
2832      globval.Elem_nFam := globval.Elem_nFam + 1;
2833      if globval.Elem_nFam <= Elem_nFamMax then
2834      begin
2835      with ElemFam[globval.Elem_nFam].ElemF do
2836      with ElementT[globval.Elem_nFam] do
2837      BEGIN
2838      Pname:=ElementName; Pkind:=Ghost; PN:=round(QL);
2839      SetDBN;
2840      END;
2841      end
2842      else
2843      writeln('Elem_nFamMax exceeded: ', globval.Elem_nFam:1,
2844      '(', Elem_nFamMax:1, ')');
2845      END;
2846      <<-----------------------------*/
2847
2848
2849    /**************************************************************************
2850      Multipole
2851    ***************************************************************************
2852
2853    <name>: Multipole,
2854            L=<length>, ( [m] )
2855            T =<bending angle>, ( [degree] )
2856            T1=<entrance angle>, ( [degree] )
2857            T2=<exit angle>, ( [degree] )
2858            Roll=<roll angle>, ( [deg], design roll angle )
2859            N =<# of kicks>,
2860            method=<method>, ( 2 or 4. The method to divide Q into slices.)
2861                    ( The detail of <method> will be discussed later.)
2862            Default value is 2.
2863            HOM=(i, <Bi>, <Ai>, ( higher order component in USA notation )
2864                 j, <Bj>, <Aj>, ( Systematic error Only )
2865                 ............    ( Random errors are assigned )
2866                 n, <Bn>, <An>); ( in a Program File using procedures )
2867
2868    Example
2869
2870      B: multipole, L=0.70, T=10.0, T1=5.0, T2=5.0,
2871         HOM=(2, -1.0, 0), N=8, Method=2;
2872
2873
2874      QF: multipole, L=0.70,
2875          HOM=(2, 2.50, 0.0,
2876               4, 1.01e7, 0.0),
2877          N=8, Method=2;
2878
2879    **************************************************************************/
2880
2881  case mpsym:  /*4*/
2882    getest__(P_expset(SET, 1 << ((long)comma)), "<, > expected", &V);
2883    GetSym__(&V);
2884    QL = 0.0;   /* L */
2885    QK = 0.0;   /* K */
2886    k1 = 0;   /* N */
2887    t = 0.0;   /* T */
2888    t1 = 0.0;   /* T1 */
2889    t2 = 0.0;   /* T2 */
2890    gap = 0.0;   /* gap */
2891    k2 = Meth_Linear;   /* method */
2892    dt = 0.0;
2893    ClearHOMandDBN(&V);
2894    P_addset(P_expset(mysys, 0), (long)lsym);
2895    P_addset(mysys, (long)nsym);
2896    P_addset(mysys, (long)mthsym);
2897    P_addset(mysys, (long)tsym);
2898    P_addset(mysys, (long)t1sym);
2899    P_addset(mysys, (long)t2sym);
2900    P_addset(mysys, (long)gapsym);
2901    P_addset(mysys, (long)rollsym);
2902    P_addset(mysys, (long)homsym);
2903    P_addset(mysys, (long)dbnsym);
2904    do {   /* read L, K, N */
2905      test__(mysys, "illegal parameter", &V);
2906      sym1 = *V.sym;
2907      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
2908      switch (sym1) {
2909
2910      case lsym:
2911        QL = EVAL_(&V);
2912        break;
2913
2914      case nsym:
2915        k1 = (long)floor(EVAL_(&V) + 0.5);
2916        break;
2917
2918      case tsym:
2919        t = EVAL_(&V);
2920        break;
2921
2922      case rollsym:
2923        dt = EVAL_(&V);
2924        break;
2925
2926      case t1sym:
2927        t1 = EVAL_(&V);
2928        break;
2929
2930      case t2sym:
2931        t2 = EVAL_(&V);
2932        break;
2933
2934      case gapsym:
2935        gap = EVAL_(&V);
2936        break;
2937
2938      case mthsym:
2939        k2 = (long)floor(EVAL_(&V) + 0.5);
2940        if (k2 != Meth_Linear) globval.MatMeth = false;
2941        if ((unsigned int)k2 >= 32 ||
2942            ((1 << k2) & ((1 << Meth_Linear) | (1 << Meth_Second) |
2943                          (1 << Meth_Fourth))) == 0)
2944          getest__(P_expset(SET, 0), "Check integrator..", &V);
2945        break;
2946
2947      case homsym:
2948        GetHOM(&V);
2949        break;
2950
2951      case dbnsym:
2952        GetDBN_(&V);
2953        break;
2954      default:
2955        break;
2956      }
2957      test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
2958             "<, > or <;> expected", &V);
2959      if (*V.sym == comma)
2960        GetSym__(&V);
2961    } while (P_inset(*V.sym, mysys));
2962    test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected.", &V);
2963    GetSym__(&V);
2964    globval.Elem_nFam++;
2965    if (globval.Elem_nFam <= Elem_nFamMax) {
2966      WITH = &ElemFam[globval.Elem_nFam-1];
2967      WITH1 = &WITH->ElemF;
2968      Mpole_Alloc(&WITH->ElemF);
2969      WITH2 = WITH1->M;
2970      memcpy(WITH1->PName, ElementName, sizeof(partsName));
2971      WITH1->Pkind = Mpole;
2972      WITH1->PL = QL;
2973      if (WITH1->PL != 0e0) {
2974        WITH2->Pthick = pthicktype(thick);
2975        WITH2->Pirho = t * M_PI / 180.0 / WITH1->PL;
2976      } else {
2977        WITH2->Pthick = pthicktype(thin);
2978        WITH2->Pirho = t * M_PI / 180.0;
2979      }
2980      WITH2->PN = k1; WITH2->Pmethod = k2;
2981      WITH2->PTx1 = t1; WITH2->PTx2 = t2; WITH2->Pgap = gap;
2982      WITH2->PdTpar = dt;
2983      AssignHOM(globval.Elem_nFam, &V);
2984      WITH2->n_design = WITH2->Porder;
2985      SetDBN(&V);
2986    } else {
2987      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
2988             globval.Elem_nFam, (long)Elem_nFamMax);
2989      exit_(1);
2990    }
2991    break;
2992
2993    /************************************************************************
2994      Wiggler
2995    *************************************************************************
2996
2997    <name>: Wiggler,
2998            L       = <length [m]>,
2999            BoBrhoV = <B/Brho [1/m]>,
3000            BoBrhoH = <B/Brho [1/m]>,
3001            Lambda  = <period [m]>,
3002            kxV     = <[m]>,
3003            kxH     = <[m]>,
3004            phi     = <phase [deg]>,
3005            harm(n, kxV, BoBrhoV, kxH, BoBrhoH, phi)
3006            ...
3007            N       = <no of integration steps>,
3008            Method  = <method>,
3009
3010    Example
3011
3012      U143: wiggler, L=4.80, K=0.5, Lambda=0.15, N=20, Method=0;
3013
3014      EPU:  wiggler, L=4.80, Lambda=0.15, N=20, Method=0,
3015            harm=(3, kxV_3, BoBrhoV_3, kxH_3, BoBrhoH_3, phi_3,
3016                  ...
3017                  5, kxV_5, BoBrhoV_5, kxH_5, BoBrhoH_5, phi_5);
3018 
3019    **************************************************************************/
3020
3021  case wglsym:
3022    getest__(P_expset(SET, 1 << ((long)comma)), "<, > expected", &V);
3023    GetSym__(&V);
3024    QL = 0e0; QK = 0e0; QKV = 0e0; QKH = 0e0; QKxV = 0e0; QKxH = 0e0;
3025    QPhi = 0e0; QKS = 0e0; k1 = 0; k2 = Meth_Linear; dt = 0e0;
3026    ClearHOMandDBN(&V);
3027    P_addset(P_expset(mysys, 0), (long)lsym);
3028    P_addset(mysys, (long)lmdsym);
3029    P_addset(mysys, (long)bobrhovsym);
3030    P_addset(mysys, (long)bobrhohsym);
3031    P_addset(mysys, (long)kxvsym);
3032    P_addset(mysys, (long)kxhsym);
3033    P_addset(mysys, (long)phisym);
3034    P_addset(mysys, (long)harmsym);
3035    P_addset(mysys, (long)nsym);
3036    P_addset(mysys, (long)mthsym);
3037    P_addset(mysys, (long)rollsym);
3038    P_addset(mysys, (long)dbnsym);
3039    do {
3040      test__(mysys, "illegal parameter", &V);
3041      sym1 = *V.sym;
3042      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
3043      switch (sym1) {   /*6*/
3044           
3045      case lsym:
3046        QL = EVAL_(&V);
3047        break;
3048
3049      case bobrhovsym:
3050        QKV = EVAL_(&V);
3051        break;
3052
3053      case bobrhohsym:
3054        QKH = EVAL_(&V);
3055        break;
3056
3057      case kxvsym:
3058        QKxV = EVAL_(&V);
3059        break;
3060
3061      case kxhsym:
3062        QKxH = EVAL_(&V);
3063        break;
3064
3065      case phisym:
3066        QPhi = EVAL_(&V);
3067        break;
3068
3069      case nsym:
3070        k1 = (long)floor(EVAL_(&V) + 0.5);
3071        break;
3072
3073      case mthsym:
3074        k2 = (long)floor(EVAL_(&V) + 0.5);
3075        if (k2 != Meth_Linear) globval.MatMeth = false;
3076        if ((unsigned int)k2 >= 32 ||
3077            ((1 << k2) &
3078             ((1 << Meth_Linear) | (1 << Meth_First) | (1 << Meth_Second) |
3079              (1 << Meth_Fourth) | (1 << Meth_genfun))) == 0)
3080          getest__(P_expset(SET, 0), "Check integrator..", &V);
3081        break;
3082
3083      case lmdsym:
3084        QKS = EVAL_(&V);
3085        break;
3086
3087      case rollsym:
3088        dt = EVAL_(&V);
3089        break;
3090
3091      case harmsym:
3092        GetHarm(&V);
3093        break;
3094
3095      case dbnsym:
3096        GetDBN_(&V);
3097        break;
3098
3099      default:
3100        break;
3101      }
3102      test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
3103             "<, > or <;> expected", &V);
3104      if (*V.sym == comma)
3105        GetSym__(&V);
3106    } while (P_inset(*V.sym, mysys));   /*5*/
3107    test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", &V);
3108    GetSym__(&V);
3109    globval.Elem_nFam++;
3110    if (globval.Elem_nFam <= Elem_nFamMax) {
3111      WITH = &ElemFam[globval.Elem_nFam-1]; WITH1 = &WITH->ElemF;
3112      memcpy(WITH1->PName, ElementName, sizeof(partsName));
3113      WITH1->PL = QL; WITH1->Pkind = Wigl;
3114      Wiggler_Alloc(&WITH->ElemF); WITH4 = WITH1->W;
3115      WITH4->Pmethod = k2; WITH4->PN = k1;
3116      WITH4->PdTpar = dt;
3117      SetDBN(&V);
3118      WITH4->lambda = QKS; WITH4->n_harm = 1; WITH4->harm[0] = 1;
3119      WITH4->kxV[0] = QKxV; WITH4->BoBrhoV[0] = QKV;
3120      WITH4->kxH[0] = QKxH; WITH4->BoBrhoH[0] = QKH;
3121      WITH4->phi[0] = QPhi;
3122      AssignHarm(globval.Elem_nFam, &V);
3123      /* Equivalent vertically focusing gradient */
3124      WITH4->PBW[HOMmax+2] = -QK*QK/2e0;
3125      CheckWiggler(globval.Elem_nFam, &V);
3126    } else {
3127      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
3128             globval.Elem_nFam, (long)Elem_nFamMax);
3129      exit_(1);
3130    }
3131    break;
3132
3133    /************************************************************************
3134      Field Map
3135    *************************************************************************
3136
3137    <name> : Fieldmap,
3138             L     = <length [m]>,
3139             N     = <no of integration steps>,
3140             file1 = <file name (lower case)>
3141 
3142    Example
3143
3144      FM: Fieldmap, L = 1.0, N = 20, file1 = "U19_Bxyz.dat";
3145
3146    **************************************************************************/
3147
3148  case fmsym:
3149    getest__(P_expset(SET, 1 << ((long)comma)), "<, > expected", &V);
3150    GetSym__(&V);
3151    QL = 0.0; k1 = 0; strcpy(str1, ""); strcpy(str2, "");
3152    P_addset(P_expset(mysys, 0), (long)lsym);
3153    P_addset(mysys, (long)nsym);
3154    P_addset(mysys, (long)fnamesym1);
3155    P_addset(mysys, (long)fnamesym2);
3156    do {
3157      test__(mysys, "illegal parameter", &V);
3158      sym1 = *V.sym;
3159      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
3160      switch (sym1) {   /*6*/
3161           
3162      case lsym:
3163        QL = EVAL_(&V);
3164        break;
3165
3166      case nsym:
3167        k1 = (long)floor(EVAL_(&V) + 0.5);
3168        break;
3169
3170      case fnamesym1:
3171        GetSym__(&V);
3172        for (i = 1; i < (signed)strlen(id_); i++) {
3173          if (id_[i] == '"') break;
3174          strncat(str1, &id_[i], 1);
3175        }
3176        GetSym__(&V);
3177        break;
3178
3179      default:
3180        break;
3181      }
3182      test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
3183             "<, > or <;> expected", &V);
3184      if (*V.sym == comma)
3185        GetSym__(&V);
3186    } while (P_inset(*V.sym, mysys));   /*5*/
3187    test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", &V);
3188    GetSym__(&V);
3189    globval.Elem_nFam++;
3190    if (globval.Elem_nFam <= Elem_nFamMax) {
3191      WITH = &ElemFam[globval.Elem_nFam-1]; WITH1 = &WITH->ElemF;
3192      memcpy(WITH1->PName, ElementName, sizeof(partsName));
3193      WITH1->PL = QL; WITH1->Pkind = FieldMap;
3194      FieldMap_Alloc(WITH1, true);
3195      WITH6 = WITH1->FM; WITH6->n_step = k1;
3196      if (CheckUDItable("energy         ", LINK) != 0) { 
3197        RefUDItable("energy         ", &globval.Energy, LINK);
3198        if (strcmp(str1, "") != 0) get_B(str1, WITH6);
3199      } else {
3200        cout << "Fieldmap: energy not defined" << endl;
3201        exit_(1);
3202      }
3203    } else {
3204      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
3205             globval.Elem_nFam, (long)Elem_nFamMax);
3206      exit_(1);
3207    }
3208    break;
3209
3210    /**********************************************************************
3211      Insertion introduced for SOLEIL using Radia Maps
3212    ***********************************************************************
3213         
3214    <name> : insertion,             
3215             N = <number of thin lenses>,
3216             scaling1 or 2 = scaling factor: should be 1. Default value
3217             file1 = <filename>, in lowercase (first order defaults)
3218             file2 = <filename>, in lowercase (second order defaults)
3219             method = <method>, ( 1 or 3. The method to divide Q into slices.)
3220             ( The detail of <method> will be discussed later.)
3221   
3222    Example
3223         
3224      ID1 : insertion, scaling2 = 1, N=10, file2="hu80_lh";
3225      ID2 : insertion, scaling1 = 1, N=10, file1="hu80_lh_bdl";
3226      ID3 : insertion, N=10, file1="hu80_lh_dbl"; file2="hu80_lh";
3227         
3228    Notes
3229      file1 and file2 must have the same structures and meshing
3230      optional parameter scaling must be at first (weird bug otherwise)
3231    **************************************************************************/
3232
3233  case idsym:
3234    getest__(P_expset(SET, 1L << ((long) comma)), "<, > expected", &V);
3235        GetSym__(&V);
3236        QK   = 0e0;
3237        QKxV = 0e0;
3238        QKS  = 0e0;
3239        k1 = 1; // number of slices of the lattice element
3240        k2 = 3; // 1 linear interpolation, 3 spline interpolation
3241        dt = 0e0;
3242        scaling1 = 1.0; // scaling factor
3243        scaling2 = 1.0; // scaling factor
3244        P_addset(P_expset(mysys, 0), (long) nsym);
3245        P_addset(mysys, (long) fnamesym1);
3246        P_addset(mysys, (long) fnamesym2);
3247        P_addset(mysys, (long) scalingsym1);
3248        P_addset(mysys, (long) scalingsym2);
3249        P_addset(mysys, (long) mthsym);
3250        do {
3251            test__(mysys, "illegal parameter", &V);
3252            sym1 = *V.sym;
3253            getest__(P_expset(SET, 1L << ((long) eql)), "<=> expected", &V);
3254           
3255            //read the parameters setting from the lattice
3256            switch (sym1) {
3257             
3258            case nsym: /* Read number of slices */
3259                k1 = abs((long) floor(EVAL_(&V)));
3260                GetSym__(&V);
3261                break;
3262
3263            case scalingsym1: /* read scaling factor for debugging purpose*/
3264                scaling1 = EVAL_(&V);
3265                break;
3266
3267            case scalingsym2: /* read scaling factor for debugging purpose*/
3268                scaling2 = EVAL_(&V);
3269                break;
3270
3271            case fnamesym1: /* Read filename for insertion device first order kicks*/
3272                firstflag = true;
3273                GetSym__(&V);
3274                for (i = 1; i < (signed) strlen(id_); i++) {
3275                    if (id_[i] == '"')
3276                        break;
3277                    strncat(str1, &id_[i], 1);
3278                }
3279                GetSym__(&V);
3280                break;
3281
3282            case fnamesym2: /* Read filename for insertion
3283             device second order kicks */
3284                secondflag = true;
3285                GetSym__(&V);
3286                for (i = 1; i < (signed) strlen(id_); i++) {
3287                    if (id_[i] == '"')
3288                        break;
3289                    strncat(str2, &id_[i], 1);
3290                }
3291                GetSym__(&V);
3292                break;
3293
3294            case mthsym: // method for interpolation: 1 means linear 3 spline
3295                k2 = (long) floor(EVAL_(&V));
3296                if (k2 != Meth_Linear)
3297                    globval.MatMeth = false;
3298                break;
3299            default:
3300                break;
3301            }
3302            if (*V.sym == comma)
3303                GetSym__(&V);
3304        } while (P_inset(*V.sym, mysys)); /*5*/
3305        GetSym__(&V);
3306        globval.Elem_nFam++;
3307
3308        /* Fills up the ID */
3309        if (globval.Elem_nFam <= Elem_nFamMax) {
3310            WITH = &ElemFam[globval.Elem_nFam - 1];
3311            WITH1 = &WITH->ElemF;
3312            memcpy(WITH1->PName, ElementName, sizeof(partsName));
3313            WITH1->Pkind = Insertion;
3314            Insertion_Alloc(&WITH->ElemF);
3315            WITH5 = WITH1->ID;
3316            WITH5->Pmethod = k2;
3317            WITH5->PN = k1;
3318            WITH5->scaling1 = scaling1;
3319            WITH5->scaling2 = scaling2;
3320
3321            // Check if filename given for first order kicks
3322            if (firstflag) {
3323                if (strcmp(str1, "") == 0)
3324                    strcpy(WITH5->fname1, "/*No_Filename1_Given*/");
3325                strcpy(WITH5->fname1, str1);
3326                // Read Id file for first order kicks
3327                WITH5->firstorder = true;
3328                Read_IDfile(WITH5->fname1, &WITH1->PL, &WITH5->nx1,
3329                        &WITH5->nz1, WITH5->tabx1, WITH5->tabz1, WITH5->thetax1,
3330                        WITH5->thetaz1);
3331            } else {
3332                strcpy(WITH5->fname1, "/*No_Filename1_Given*/");
3333            }
3334
3335            // Check if filename given for Second order kicks
3336            if (secondflag) {
3337                if (strcmp(str2, "") != 0)
3338                    strcpy(WITH5->fname2, "/*No_Filename2_Given*/");
3339                strcpy(WITH5->fname2, str2);
3340                WITH5->secondorder = secondflag;
3341                // Read Id file for second order kicks
3342                Read_IDfile(WITH5->fname2, &WITH1->PL, &WITH5->nx2,
3343                        &WITH5->nz2, WITH5->tabx2, WITH5->tabz2,
3344                        WITH5->thetax2, WITH5->thetaz2);
3345            } else {
3346                strcpy(WITH5->fname2, "/*No_Filename2_Given*/");
3347            }
3348
3349            // check whether no Radia filename read: something is wrong
3350            if (!firstflag && !secondflag) {
3351                printf("Error no Insertion filename found as"
3352                    " an input in lattice file\n");
3353                exit_(-1);
3354            }
3355
3356            if (k2 == 3 | k2 == 2) { // cubic interpolation
3357                WITH5->linear = false;
3358            } else { // linear interpolation
3359                WITH5->linear = true;
3360            }
3361
3362            // stuff for spline interpolation
3363            if (!WITH5->linear) {
3364              if (firstflag){
3365                WITH5->tx1 = dmatrix(1, WITH5->nz1, 1, WITH5->nx1);
3366                WITH5->tz1 = dmatrix(1, WITH5->nz1, 1, WITH5->nx1);
3367                WITH5->TabxOrd1 = (double *) malloc((WITH5->nx1) * sizeof(double));
3368                WITH5->TabzOrd1 = (double *) malloc((WITH5->nz1) * sizeof(double));
3369                WITH5->f2x1 = dmatrix(1, WITH5->nz1, 1, WITH5->nx1);
3370                WITH5->f2z1 = dmatrix(1, WITH5->nz1, 1, WITH5->nx1);
3371                Matrices4Spline(WITH5,1);}
3372
3373              if (secondflag){
3374                WITH5->tx2 = dmatrix(1, WITH5->nz2, 1, WITH5->nx2);
3375                WITH5->tz2 = dmatrix(1, WITH5->nz2, 1, WITH5->nx2);
3376                WITH5->TabxOrd2 = (double *) malloc((WITH5->nx2) * sizeof(double));
3377                WITH5->TabzOrd2 = (double *) malloc((WITH5->nz2) * sizeof(double));
3378                WITH5->f2x2 = dmatrix(1, WITH5->nz2, 1, WITH5->nx2);
3379                WITH5->f2z2 = dmatrix(1, WITH5->nz2, 1, WITH5->nx2);
3380                Matrices4Spline(WITH5,2);}
3381
3382            }
3383            // to put somewhere
3384            //      /** Free memory **/
3385            //      free(tab1);
3386            //      free(tab2);
3387            //
3388            //      free_matrix(tx,1,nz,1,nx);
3389            //      free_matrix(tz,1,nz,1,nx);
3390            //      free_matrix(f2x,1,nz,1,nx);
3391            //      free_matrix(f2z,1,nz,1,nx);
3392
3393        } else {
3394            printf("Elem_nFamMax exceeded: %ld(%ld)\n", globval.Elem_nFam,
3395                    (long) Elem_nFamMax);
3396            exit_(1);
3397        }
3398        break;
3399
3400    /**************************************************************************
3401       Spreader
3402    **************************************************************************
3403
3404    <name>: Spreader,
3405
3406    Example
3407
3408      S1: Spreader;
3409
3410    *************************************************************************/
3411
3412  case sprsym:
3413    getest__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
3414             "<, > or <;> expected", &V);
3415    if (*V.sym == comma) {
3416      getest__(P_addset(P_expset(SET4, 0), (long)dbnsym),
3417               "illegal parameter", &V);
3418      sym1 = *V.sym;
3419      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
3420      if (sym1 == dbnsym)
3421        GetDBN_(&V);
3422      test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", &V);
3423    }
3424    GetSym__(&V);
3425    globval.Elem_nFam++;
3426    if (globval.Elem_nFam <= Elem_nFamMax) {
3427      WITH = &ElemFam[globval.Elem_nFam-1];
3428      WITH1 = &WITH->ElemF;
3429      memcpy(WITH1->PName, ElementName, sizeof(partsName));
3430      WITH1->PL = *V.rnum;
3431      WITH1->Pkind = PartsKind(Spreader);
3432      Spreader_Alloc(&WITH->ElemF);
3433    } else {
3434      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
3435             globval.Elem_nFam, (long)Elem_nFamMax);
3436      exit_(1);
3437    }
3438    break;
3439
3440    /**************************************************************************
3441       Recombiner
3442    **************************************************************************
3443
3444    <name>: Recombiner,
3445
3446    Example
3447
3448      S1: Recombiner;
3449
3450    *************************************************************************/
3451
3452  case recsym:
3453    getest__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
3454             "<, > or <;> expected", &V);
3455    if (*V.sym == comma) {
3456      getest__(P_addset(P_expset(SET4, 0), (long)dbnsym),
3457               "illegal parameter", &V);
3458      sym1 = *V.sym;
3459      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
3460      if (sym1 == dbnsym)
3461        GetDBN_(&V);
3462      test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", &V);
3463    }
3464    GetSym__(&V);
3465    globval.Elem_nFam++;
3466    if (globval.Elem_nFam <= Elem_nFamMax) {
3467      WITH = &ElemFam[globval.Elem_nFam-1];
3468      WITH1 = &WITH->ElemF;
3469      memcpy(WITH1->PName, ElementName, sizeof(partsName));
3470      WITH1->PL = *V.rnum;
3471      WITH1->Pkind = PartsKind(Recombiner);
3472      Recombiner_Alloc(&WITH->ElemF);
3473    } else {
3474      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
3475             globval.Elem_nFam, (long)Elem_nFamMax);
3476      exit_(1);
3477    }
3478    break;
3479
3480    /**************************************************************************
3481      Solenoid
3482    ***************************************************************************
3483
3484    <name>: Solenoid,
3485            L=<length>, ( [m] )
3486            BoBrho = <B/Brho [1/m]>,
3487            N =<# of kicks>,
3488            method=<method>
3489
3490    Example
3491
3492      B: solenoid, L=0.70, BoBrho=10.0;
3493
3494    **************************************************************************/
3495
3496  case solsym:
3497    getest__(P_expset(SET, 1 << ((long)comma)), "<, > expected", &V);
3498    GetSym__(&V);
3499    QL = 0.0; /* L */
3500    QK = 0.0; /* K */
3501    k1 = 0;   /* N */
3502    P_addset(P_expset(mysys, 0), (long)lsym);
3503    P_addset(mysys, (long)bobrhosym);
3504    P_addset(mysys, (long)nsym);
3505    do { /* read L, K, N */
3506      test__(mysys, "illegal parameter", &V);
3507      sym1 = *V.sym;
3508      getest__(P_expset(SET, 1 << ((long)eql)), "<=> expected", &V);
3509      switch (sym1) {
3510
3511      case lsym:
3512        QL = EVAL_(&V);
3513        break;
3514
3515      case bobrhosym:
3516        QK = EVAL_(&V);
3517        break;
3518
3519      case nsym:
3520        k1 = (long)floor(EVAL_(&V) + 0.5);
3521        break;
3522
3523      default:
3524        cout << "Solenoid: undef. case" << endl;
3525        exit_(1);
3526        break;
3527      }
3528      test__(P_expset(SET, (1 << ((long)comma)) | (1 << ((long)semicolon))),
3529             "<, > or <;> expected", &V);
3530      if (*V.sym == comma)
3531        GetSym__(&V);
3532    } while (P_inset(*V.sym, mysys));
3533    test__(P_expset(SET, 1 << ((long)semicolon)), "<;> expected.", &V);
3534    GetSym__(&V);
3535    globval.Elem_nFam++;
3536    if (globval.Elem_nFam <= Elem_nFamMax) {
3537      WITH = &ElemFam[globval.Elem_nFam-1];
3538      WITH1 = &WITH->ElemF;
3539      Solenoid_Alloc(&WITH->ElemF);
3540      WITH7 = WITH1->Sol;
3541      memcpy(WITH1->PName, ElementName, sizeof(partsName));
3542      WITH1->Pkind = Solenoid;
3543      WITH1->PL = QL; WITH7->N = k1; WITH7->BoBrho = QK;
3544    } else {
3545      printf("Elem_nFamMax exceeded: %ld(%ld)\n",
3546             globval.Elem_nFam, (long)Elem_nFamMax);
3547      exit_(1);
3548    }
3549    break;
3550
3551    /**************************************************************************
3552      BLOCK DEFINITION
3553    **************************************************************************/
3554
3555  case ident:
3556  case intcon:
3557  case invsym:   /* Block Definition */
3558    ProcessBlockInput(&V);
3559    break;
3560  default:
3561    break;
3562  }/*3.5:of CASE*/
3563
3564  Result = true;
3565
3566 _L9999:
3567
3568  return Result;
3569}  /*of procedure Lat_DealElement*/
3570
3571
3572static void errorm___(const char *cmnt, struct LOC_Lattice_Read *LINK)
3573{
3574  /*write(fo, ' ****')*/
3575  /*error*/
3576  if (LINK->cc > LINK->errpos) {
3577    fprintf(*LINK->fo, "%*c^%.80s", (int)(LINK->cc - LINK->errpos),
3578            ' ', cmnt);
3579    LINK->errpos = LINK->cc + 3;
3580  }
3581  while (!P_eof(*LINK->fi))
3582    Lat_Nextch(LINK->fi, LINK->fo, &LINK->cc, &LINK->ll, &LINK->errpos,
3583               &LINK->lc, &LINK->chin, &LINK->skipflag, LINK->line, LINK);
3584  ErrFlag = true;
3585  longjmp(LINK->_JL9999, 1);
3586}
3587
3588
3589static void GetSym___(struct LOC_Lattice_Read *LINK)
3590{
3591  /* reads next symbol  */
3592  /*GetSym*/
3593  Lat_GetSym(LINK->fi, LINK->fo, &LINK->cc, &LINK->ll, &LINK->errpos,
3594             &LINK->lc, &LINK->nkw, &LINK->inum, (long)emax, (long)emin,
3595             (long)kmax, nmax, &LINK->chin, LINK->id, &LINK->rnum,
3596             &LINK->skipflag, &LINK->rsvwd, LINK->line, &LINK->sym, LINK->key,
3597             LINK->ksy, LINK->sps, LINK);
3598}
3599
3600
3601static void test___(long *s1, const char *cmnt, struct LOC_Lattice_Read *LINK)
3602{
3603  /*test*/
3604  if (!P_inset(LINK->sym, s1))
3605    errorm___(cmnt, LINK);
3606}
3607
3608
3609static void getest___(long *s1, const char *cmnt, struct LOC_Lattice_Read *LINK)
3610{
3611  /*test*/
3612  GetSym___(LINK);
3613  if (!P_inset(LINK->sym, s1))
3614    errorm___(cmnt, LINK);
3615}
3616
3617
3618/* Local variables for init_reserved_words: */
3619struct LOC_init_reserved_words
3620{
3621  struct LOC_Lattice_Read *LINK;
3622};
3623
3624
3625static void Reg(const char *name, Lat_symbol ks,
3626                struct LOC_init_reserved_words *LINK)
3627{
3628  LINK->LINK->nkw++;  /* incrementation of the number of keywords */
3629  if (LINK->LINK->nkw > Lat_nkw_max) {
3630    cout << "Reg: Lat_nkw_max exceeded " << LINK->LINK->nkw
3631         << "(" << Lat_nkw_max << ")" << endl;
3632  }
3633  memcpy(LINK->LINK->key[LINK->LINK->nkw - 1], name, sizeof(alfa_));
3634  LINK->LINK->ksy[LINK->LINK->nkw - 1] = ks;
3635}
3636
3637
3638static void init_reserved_words(struct LOC_Lattice_Read *LINK)
3639{
3640  struct LOC_init_reserved_words V;
3641
3642  V.LINK = LINK;
3643  LINK->nkw = 0; /* Number of keywords equals zero */
3644  /*-------------------------------------------------------------
3645  To define reserved symbol in the lattice reading,
3646  MUST follow alphabetical list !!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3647  --------------------------------------------------------------*/
3648  Reg("and            ", andsym, &V);
3649  Reg("beam           ", bemsym, &V);
3650  Reg("bending        ", bndsym, &V);
3651  Reg("cavity         ", cavsym, &V);
3652  Reg("cell           ", celsym, &V);
3653  Reg("chromaticity   ", chmsym, &V);
3654  Reg("corrector      ", corsym, &V);  /* corrector */
3655  Reg("dbname         ", dbnsym, &V);
3656  Reg("define         ", defsym, &V);
3657  Reg("dispersion     ", dspsym, &V);
3658  Reg("drift          ", drfsym, &V);
3659  Reg("dt             ", dtsym, &V);
3660  Reg("end            ", endsym, &V);
3661  Reg("ff1            ", ff1sym, &V);     /* Laurent */
3662  Reg("ff2            ", ff2sym, &V);     /* Laurent */
3663  Reg("ffscaling      ", ffscalingsym, &V);/* quadff scaling  J.Zhang */
3664  Reg("fieldmap       ", fmsym, &V);
3665  Reg("file1          ", fnamesym1, &V);   /* ID Laurent */
3666  Reg("file2          ", fnamesym2, &V);   /* ID Laurent */
3667  Reg("focusing       ", fcssym, &V);
3668  Reg("frequency      ", frqsym, &V);
3669  Reg("fringe         ", frgsym, &V);
3670  Reg("galilean       ", xytsym, &V);
3671  Reg("gap            ", gapsym, &V);
3672  Reg("ghost          ", gstsym, &V);
3673  Reg("harm           ", harmsym, &V);
3674  Reg("harnum         ", harnumsym, &V);
3675  Reg("hom            ", homsym, &V);
3676  Reg("horizontal     ", horsym, &V); /* with "corrector", define the horizontal corrector in the lattice*/
3677  Reg("insertion      ", idsym, &V);   /* ID Laurent */
3678  Reg("inv            ", invsym, &V);
3679  Reg("kick           ", corkicksym, &V); /* with "corrector", define the kick angle of the corrector ,  Jianfeng Zhang*/
3680  Reg("kicker         ", kicksym, &V);
3681  Reg("ks             ", kssym, &V);
3682  Reg("lambda         ", lmdsym, &V);
3683  Reg("lattice        ", latsym, &V);
3684  Reg("marker         ", mrksym, &V);
3685  Reg("matrix         ", matsym, &V);
3686  Reg("method         ", mthsym, &V);
3687  Reg("monitor        ", monsym, &V);
3688  Reg("multipole      ", mpsym, &V);
3689  Reg("nonlinear      ", nbdsym, &V);
3690  Reg("parameter      ", prmsym, &V);
3691  Reg("position       ", possym, &V);
3692  Reg("print          ", prnsym, &V);
3693  Reg("quadrupole     ", qdsym, &V);
3694  Reg("recombiner     ", recsym, &V);
3695  Reg("roll           ", rollsym, &V);
3696  Reg("scaling1       ", scalingsym1, &V); /* ID Laurent */
3697  Reg("scaling2       ", scalingsym2, &V); /* ID Laurent */
3698  Reg("sextupole      ", sexsym, &V);
3699  Reg("solenoid       ", solsym, &V);
3700  Reg("spreader       ", sprsym, &V);
3701  Reg("symmetry       ", symsym, &V);
3702  Reg("t1             ", t1sym, &V);
3703  Reg("t2             ", t2sym, &V);
3704  Reg("table          ", tblsym, &V);
3705  Reg("task           ", tsksym, &V);
3706  Reg("tilt           ", tiltsym, &V); // added for compatibility with Tracy II
3707  Reg("type           ", typsym, &V);
3708  Reg("use            ", usesym, &V);
3709  Reg("vertical       ", versym, &V); /* with "corrector", define the vertical corrector in the lattice*/
3710  Reg("voltage        ", vrfsym, &V);
3711  Reg("wiggler        ", wglsym, &V);
3712
3713  if (trace) fprintf(stdout,"Nb of keywords = %ld (%d)\n",
3714                     LINK->nkw, Lat_nkw_max);
3715
3716  LINK->sps[(int)'+'] = plus_;
3717  LINK->sps[(int)'-'] = minus_;
3718  LINK->sps[(int)'('] = lparent;
3719  LINK->sps[(int)')'] = rparent;
3720  LINK->sps[(int)'='] = eql;
3721  LINK->sps[(int)','] = comma;
3722  LINK->sps[(int)'['] = lbrack;
3723  LINK->sps[(int)']'] = rbrack;
3724  LINK->sps[(int)'\''] = squote;
3725  LINK->sps[(int)'&'] = andsy;
3726  LINK->sps[(int)';'] = semicolon;
3727  LINK->sps[(int)'/'] = rdiv;
3728  LINK->sps[(int)':'] = colon;
3729
3730  if (trace)
3731    printf("%d %d %d %d %d %d %d %d %d %d %d %d %d %d \n",
3732           (int)'+', (int)'-', (int)'(', (int)')', (int)'=',
3733           (int)',', (int)'[', (int)'[', (int)']', (int)'\'',
3734           (int)'&', (int)';', (int)'/', (int)':');
3735
3736  LINK->lc = 0;   /* reset line counter */
3737  LINK->ll = 0;   /* reset line length  */
3738  LINK->cc = 0;   /* reset char counter */
3739  LINK->errpos = 0;   /* reset error position */
3740  LINK->chin = ' ';   /* reset current char   */
3741  LINK->skipflag = false;   /* reset skip flag  */
3742  P_addset(P_expset(LINK->defbegsys, 0), (long)ident);
3743  P_addset(P_expset(LINK->elmbegsys, 0), (long)qdsym);
3744  P_addset(LINK->elmbegsys, (long)sexsym); /*link the lattice element name*/
3745  P_addset(LINK->elmbegsys, (long)corsym);
3746  P_addset(LINK->elmbegsys, (long)bemsym);
3747  P_addset(LINK->elmbegsys, (long)gstsym);
3748  P_addset(LINK->elmbegsys, (long)mrksym);
3749  P_addset(LINK->elmbegsys, (long)nbdsym);
3750  P_addset(LINK->elmbegsys, (long)frgsym);
3751  P_addset(LINK->elmbegsys, (long)xytsym);
3752  P_addset(LINK->elmbegsys, (long)drfsym);
3753  P_addset(LINK->elmbegsys, (long)bndsym);
3754  P_addset(LINK->elmbegsys, (long)wglsym);
3755  P_addset(LINK->elmbegsys, (long)mpsym);
3756  P_addset(LINK->elmbegsys, (long)cavsym);
3757  P_addset(LINK->elmbegsys, (long)idsym);      /* ID Laurent */
3758  P_addset(LINK->elmbegsys, (long)fmsym);
3759  P_addset(LINK->elmbegsys, (long)sprsym);
3760  P_addset(LINK->elmbegsys, (long)recsym);
3761  P_addset(LINK->elmbegsys, (long)solsym);
3762  P_addset(LINK->elmbegsys, (long)fnamesym1);  /* ID file name Laurent */
3763  P_addset(LINK->elmbegsys, (long)fnamesym2);  /* ID file name Laurent */
3764//  P_addset(LINK->elmbegsys, (long)scalingsym); /* ID scale factor Laurent */
3765}
3766
3767/* Local variables for DealWithDefns: */
3768struct LOC_DealWithDefns
3769{
3770  struct LOC_Lattice_Read *LINK;
3771};
3772
3773static double EVAL__(struct LOC_DealWithDefns *LINK)
3774{
3775  return (Lat_EVAL(LINK->LINK->fi, LINK->LINK->fo, &LINK->LINK->cc,
3776                   &LINK->LINK->ll, &LINK->LINK->errpos, &LINK->LINK->lc,
3777                   &LINK->LINK->nkw, &LINK->LINK->inum, (long)emax,
3778                   (long)emin, (long)kmax, nmax, &LINK->LINK->chin,
3779                   LINK->LINK->id, &LINK->LINK->rnum, &LINK->LINK->skipflag,
3780                   &LINK->LINK->rsvwd, LINK->LINK->line, &LINK->LINK->sym,
3781                   LINK->LINK->key, LINK->LINK->ksy, LINK->LINK->sps,
3782                   LINK->LINK));
3783}
3784
3785
3786/******************************************************
3787 *                                                   *
3788 *                  P A R S E R                      *
3789 *                                                   *
3790 ******************************************************/
3791
3792static void DealWithDefns(struct LOC_Lattice_Read *LINK)
3793{  /*0*/
3794  struct LOC_DealWithDefns V;
3795  partsName idsave, ElementName, BlockName, IdentName;
3796  long i, j, k, k1;
3797  symset SET;
3798  long SET1[(long)ident / 32 + 2];
3799  long SET2[(long)period_ / 32 + 2];
3800  long SET3[(long)invsym / 32 + 2];
3801  _REC_BlockStype *WITH;
3802  long FORLIM;
3803  long SET4[(long)symsym / 32 + 2];
3804  long SET5[(long)endsym / 32 + 2];
3805
3806  /************** DEAL WITH DEFINITIONS *********************************/
3807
3808  V.LINK = LINK;
3809  GetSym___(LINK);
3810  if (LINK->sym != latsym) {  /*1*/
3811    test___(P_expset(SET, 0), "<illegal operand> detected", LINK);
3812    return;
3813  }
3814  /****** The first word must be 'lattice' **********/
3815  getest___(P_expset(SET, 1 << ((long)semicolon)), "<;> expected", LINK);
3816  // Test whether expression exists
3817  getest___(P_addset(P_expset(SET1, 0), (long)ident),
3818            "<identifier> expected", LINK);
3819
3820  /***************************************************************************/
3821  if (LINK->sym == ident) {
3822    do {   /*2*/
3823      if (LINK->sym == ident) {  /*2.5:-----------------------*/
3824        memcpy(idsave, LINK->id, sizeof(alfa_));
3825        memset(idsave + sizeof(alfa_), ' ', sizeof(partsName) - sizeof(alfa_));
3826        memcpy(ElementName, LINK->id, sizeof(alfa_));
3827        memset(ElementName + sizeof(alfa_), ' ',
3828               sizeof(partsName) - sizeof(alfa_));
3829        memcpy(BlockName, LINK->id, sizeof(alfa_));
3830        memset(BlockName + sizeof(alfa_), ' ',
3831               sizeof(partsName) - sizeof(alfa_));
3832        P_addset(P_expset(SET2, 0), (long)colon);
3833        P_addset(SET2, (long)eql);
3834        getest___(P_addset(SET2, (long)period_),
3835                  "<colon>, <=> or <.> expected", LINK);
3836        if (LINK->sym == colon)   /*3:of IF sym=colon*/ {  /*3*/
3837          P_addset(P_expset(SET3, 0), (long)ident);
3838          P_addset(SET3, (long)intcon);
3839          getest___(P_setunion(SET, LINK->elmbegsys,
3840                               P_addset(SET3, (long)invsym)),
3841                    "<identifier>, <integer> or <INV> expected", LINK);
3842          P_addset(P_expset(SET3, 0), (long)ident);
3843          P_addset(SET3, (long)intcon);
3844          if (P_inset(LINK->sym,
3845                      P_setunion(SET, LINK->elmbegsys,
3846                                 P_addset(SET3, (long)invsym)))) {
3847            if (!Lat_DealElement(LINK->fi, LINK->fo, &LINK->cc,
3848                                 &LINK->ll,
3849                                 &LINK->errpos, &LINK->lc, &LINK->nkw,
3850                                 &LINK->inum, (long)emax, (long)emin,
3851                                 (long)kmax, nmax, &LINK->chin,
3852                                 LINK->id,
3853                                 ElementName, BlockName, &LINK->rnum,
3854                                 &LINK->skipflag, &LINK->rsvwd,
3855                                 LINK->line,
3856                                 &LINK->sym, LINK->key, LINK->ksy,
3857                                 LINK->sps, LINK))
3858              longjmp(LINK->_JL9999, 1);
3859          }
3860                 
3861        } else {
3862          /**************************************************************
3863           *                                                           *
3864           *         P A R A M E T E R  A S S I G N M E N T            *
3865           *                                                           *
3866           **************************************************************/
3867
3868          if (LINK->sym == eql) {  /*3:of parameter*/
3869            memcpy(IdentName, idsave, sizeof(partsName));
3870            i = CheckUDItable(IdentName, LINK);
3871            if (i == 0)
3872              EnterUDItable(IdentName, EVAL__(&V), LINK);
3873            else
3874              ModUDItable(i, EVAL__(&V), LINK);
3875            test___(P_expset(SET, 1 << ((long)semicolon)),
3876                    "<;> expected", LINK);
3877            GetSym___(LINK);
3878            /*3:of parameter*/
3879            /*-----
3880              IdentName:=idsave;
3881              i:=CheckElementtable(IdentName);
3882              IF i=0 THEN Test([], '<element name> expected');
3883              getest([lsym, tsym, t1sym, t2sym, gapsym, ksym],
3884                     'illegal component');
3885              sym1:=sym;
3886              getest([eql], '<=> expected');
3887              case sym1 of
3888              lsym:  ElemFam[i].ElemF.PL :=Eval;
3889              ksym:  ElemFam[i].ElemF.Pk :=Eval;
3890              tsym:  ElemFam[i].ElemF.Pt :=Eval;
3891              t1sym: ElemFam[i].ElemF.Pt1:=Eval;
3892              t2sym: ElemFam[i].ElemF.Pt2:=Eval;
3893              gapsym: ElemFam[i].ElemF.Pgap:=Eval;
3894              END;
3895              test([semicolon], '<;> expected');
3896              GetSym;
3897              -----*/
3898            /*3:of parameter */
3899          }  /*3:of parameter */
3900        }
3901
3902        /*****************************
3903         *                           *
3904         *     DEAL WITH ELEMENT     *
3905         *                           *
3906         *****************************/
3907
3908      }  /*2.5*/
3909
3910      /**************************************************************
3911       *                                                           *
3912       *         C E L L    D E F I N I T I O N                    *
3913       *                                                           *
3914       *************************************************************
3915
3916       CELL : <block name>, SYMMETRY=<symmetry>;
3917
3918              <block name>:== name of a block.
3919              <symmetry>:== number of supersymmetry:== number of the block/ring
3920
3921       Example
3922
3923         CELL : BL1, Symmetry=12;
3924
3925      ************************************************************************/
3926
3927      if (LINK->sym == celsym) {  /*3*/
3928        getest___(P_expset(SET, 1 << ((long)colon)),
3929                  "<colon> expected", LINK);
3930        getest___(P_addset(P_expset(SET1, 0), (long)ident),
3931                  "<Block name> expected", LINK);
3932        i = CheckBLOCKStable(LINK->id, LINK);
3933        if (i == 0)
3934          test___(P_expset(SET, 0),
3935                  "<Block name> expected", LINK);
3936        k = 0;
3937        if (i != 0) {  /*4*/
3938          WITH = &LINK->BlockS[i - 1];
3939          FORLIM = WITH->BOWARI;
3940          for (j = WITH->BSTART - 1; j < FORLIM; j++) {  /*6*/
3941            k++;
3942            if (j < NoBEmax)
3943              k1 = LINK->Bstack[j];
3944            else {
3945              printf("** DealWithDefns: NoBEmax exceeded %ld (%d)\n",
3946                     j+1, NoBEmax);
3947              exit(1);
3948            }
3949            if (k <= Cell_nLocMax)
3950              Cell[k].Fnum = k1;
3951            else {
3952              printf("** Cell_nLocMax exhausted: %ld(%ld)\n",
3953                     k, (long)Cell_nLocMax);
3954              exit_(1);
3955            }
3956          }
3957          /*5*/
3958        }
3959        /*4*/
3960        if (k <= Cell_nLocMax)
3961          globval.Cell_nLoc = k;   /*number of Elements in a cell*/
3962        else {
3963          printf("** Cell_nLocMax exhausted: %ld(%ld)\n",
3964                 k, (long)Cell_nLocMax);
3965          exit_(1);
3966        }
3967        getest___(P_expset(SET, 1 << ((long)comma)),
3968                  "<, > expected", LINK);
3969        getest___(P_addset(P_expset(SET4, 0), (long)symsym),
3970                  "<symmetry> expected", LINK);
3971        getest___(P_expset(SET, 1 << ((long)eql)),
3972                  "<=> expected", LINK);
3973        LINK->Symmetry = (long)floor(EVAL__(&V) + 0.5);
3974        if (LINK->Symmetry >= 1)
3975          LINK->Ring = true;
3976        else {
3977          LINK->Symmetry = 1;
3978          LINK->Ring = false;
3979        }
3980        test___(P_expset(SET, 1 << ((long)semicolon)),
3981                "<;> expected", LINK);
3982        GetSym___(LINK);
3983      }  /*3: of celsym*/
3984
3985
3986      switch (LINK->sym) {   /*2*/
3987
3988        /******************************************
3989
3990                     PRINT element-name
3991                     PRINT block_name
3992                     PRINT parameter
3993
3994        ******************************************/
3995
3996      case prnsym:  /*4*/
3997        getest___(P_addset(P_expset(SET1, 0), (long)ident),
3998                  "<identifiler> expected", LINK);
3999        memcpy(IdentName, LINK->id, sizeof(alfa_));
4000        memset(IdentName + sizeof(alfa_), ' ',
4001               sizeof(partsName) - sizeof(alfa_));
4002        i = CheckElementtable(IdentName, LINK);
4003        if (i == 0) {   /*PrintElementParam(i)*/
4004          i = CheckBLOCKStable(IdentName, LINK);
4005          if (i == 0) {   /*PrintBlockParam(i)*/
4006            i = CheckUDItable(IdentName, LINK);
4007            if (i == 0)
4008              getest___(P_expset(SET, 0),
4009                        " invalid expression", LINK);
4010            /*PrintUDIParam(i)*/
4011          }
4012        }
4013        if (i != 0) {
4014          getest___(P_expset(SET, 1 << ((long)semicolon)),
4015                    "<;> expected", LINK);
4016          GetSym___(LINK);
4017        }
4018        break;
4019        /*4*/
4020      default:
4021        break;
4022      }/*3:of CASE*/
4023
4024    } while (LINK->sym == (long)prnsym || LINK->sym == (long)celsym ||
4025             LINK->sym == (long)dspsym ||
4026             LINK->sym == (long)chmsym || LINK->sym == (long)ident);
4027  }
4028
4029  test___(P_addset(P_expset(SET5, 0), (long)endsym), "<END> expected", LINK);
4030  getest___(P_expset(SET, 1 << ((long)semicolon)), "<;> expexted", LINK);
4031
4032
4033
4034  /*8888888888*/
4035  /*5*/
4036  /*6*/
4037  /*6*/
4038  /*5*/
4039  /*1*/
4040  /*1*/
4041}  /*0*/
4042
4043
4044void GetEnergy(struct LOC_Lattice_Read *LINK)
4045{
4046  long k;
4047
4048  k = CheckUDItable("energy         ", LINK);
4049  if (k == 0) {
4050    printf("> Beam energy is not defined.\n");
4051    printf("  Input beam energy in [GeV] := ");
4052    scanf("%lg%*[^\n]", &globval.Energy);
4053    getchar();
4054    EnterUDItable("energy         ", globval.Energy, LINK);
4055  } else
4056    RefUDItable("energy         ", &globval.Energy, LINK);
4057}
4058
4059
4060void GetRingType(struct LOC_Lattice_Read *LINK)
4061{
4062  long k;
4063
4064  k = CheckUDItable("ringtype       ", LINK);
4065  if (k == 0L) {
4066    fprintf(stdout,"> Ring type is not defined, default is ring.\n");
4067    globval.RingType = 1;
4068  } else {
4069    globval.RingType = (int) LINK->UDItable[k - 1L].Uvalue;
4070    if (globval.RingType != 1 && globval.RingType != 0) {
4071      printf("> ringtype variable is not defined"
4072              " properly in the lattice file\n");
4073      printf("> ringtype set to 1 means ring\n");
4074      printf("> ringtype set to 0 means transfer line\n");
4075      exit_(1);
4076    }
4077  }
4078}
4079
4080
4081/****************************************************************************/
4082/* void GetDP(struct LOC_Lattice_Read *LINK)
4083
4084   Purpose:
4085       Define particle energy offset read from lattice file
4086
4087   Input:
4088       none
4089
4090   Output:
4091       none
4092
4093   Return:
4094       none
4095
4096   Global variables:
4097       none
4098
4099   Specific functions:
4100       none
4101
4102   Comments:
4103       none
4104
4105****************************************************************************/
4106static void GetDP(struct LOC_Lattice_Read *LINK)
4107{
4108  long k;
4109
4110  k = CheckUDItable("dp             ", LINK);
4111  if (k != 0) {
4112      RefUDItable("dp             ", &globval.dPcommon, LINK);
4113      return;
4114    }
4115  printf("> dP/P is not defined.\n");
4116  printf("  Input dP/P := ");
4117  scanf("%lg%*[^\n]", &globval.dPcommon);
4118  getchar();
4119  EnterUDItable("dp             ", globval.dPcommon, LINK);
4120}
4121
4122/****************************************************************************/
4123/* void GetCODEPS(LOC_Lattice_Read *LINK)
4124
4125   Purpose:
4126       Read and assign cod precision read from lattice file
4127
4128   Input:
4129       none
4130
4131   Output:
4132       none
4133
4134   Return:
4135       none
4136
4137   Global variables:
4138       none
4139
4140   Specific functions:
4141       none
4142
4143   Comments:
4144       none
4145
4146****************************************************************************/
4147static void GetCODEPS(struct LOC_Lattice_Read *LINK)
4148{
4149  long k;
4150
4151  k = CheckUDItable("codeps         ", LINK);
4152  if (k != 0) {
4153    RefUDItable("codeps         ", &globval.CODeps, LINK);
4154    return;
4155  }
4156  printf("> CODEPS is not defined.\n");
4157  printf("  Input CODEPS := ");
4158  scanf("%lg%*[^\n]", &globval.CODeps);
4159  getchar();
4160  EnterUDItable("codeps         ", globval.Energy, LINK);
4161}
4162
4163/****************************************************************************/
4164/* Local double Circumference(struct LOC_Lattice_Read *LINK)
4165
4166   Purpose:
4167
4168
4169   Input:
4170       none
4171
4172   Output:
4173       none
4174
4175   Return:
4176       none
4177
4178   Global variables:
4179       none
4180
4181   Specific functions:
4182       none
4183
4184   Comments:
4185       none
4186
4187****************************************************************************/
4188static double Circumference(struct LOC_Lattice_Read *LINK)
4189{
4190  long i;
4191  double S;
4192  long FORLIM;
4193
4194  S = 0.0;
4195  FORLIM = globval.Cell_nLoc;
4196  for (i = 1; i <= FORLIM; i++)
4197    S += ElemFam[Cell[i].Fnum - 1].ElemF.PL;
4198  return S;
4199}
4200
4201
4202/****************************************************************************/
4203/* Local void RegisterKids(struct LOC_Lattice_Read *LINK)
4204
4205   Purpose:
4206
4207
4208   Input:
4209       none
4210
4211   Output:
4212       none
4213
4214   Return:
4215       none
4216
4217   Global variables:
4218       none
4219
4220   Specific functions:
4221       none
4222
4223   Comments:
4224       none
4225 ****************************************************************************/
4226
4227static void RegisterKids(struct LOC_Lattice_Read *LINK)
4228{
4229  long i, FORLIM;
4230  ElemFamType *WITH;
4231
4232  if (globval.Elem_nFam <= Elem_nFamMax) {
4233    FORLIM = globval.Elem_nFam;
4234    for (i = 0; i < FORLIM; i++)
4235      ElemFam[i].nKid = 0;
4236  } else {
4237    printf("Elem_nFamMax exceeded: %ld(%d)\n",
4238           globval.Elem_nFam, Elem_nFamMax);
4239    exit_(1);
4240  }
4241
4242  FORLIM = globval.Cell_nLoc;
4243  for (i = 1; i <= FORLIM; i++) {
4244    WITH = &ElemFam[Cell[i].Fnum - 1];
4245    WITH->nKid++;
4246    if (WITH->nKid <= nKidMax) {
4247      WITH->KidList[WITH->nKid - 1] = i;
4248      Cell[i].Knum = WITH->nKid;
4249    } else
4250      printf("nKidMax exceeded: %d(%d)\n", WITH->nKid, nKidMax);
4251  }
4252}
4253
4254
4255/****************************************************************************/
4256/* void PrintResult(struct LOC_Lattice_Read *LINK)
4257
4258   Purpose:
4259       Print Lattice statistics
4260
4261   input:
4262       LINK
4263
4264   output:
4265       none
4266
4267   return:
4268       none
4269
4270   global variables:
4271       none
4272
4273   specific functions:
4274       none
4275
4276   comments
4277       none
4278
4279****************************************************************************/
4280void PrintResult(struct LOC_Lattice_Read *LINK)
4281{
4282  long j, nKid, FORLIM;
4283  struct tm *newtime;
4284
4285  /* Get time and date */
4286  newtime = GetTime();
4287
4288  printf("\n");
4289  printf("  TRACY III v. 3.5 compiled on %s\n",__DATE__);
4290  printf("\n");
4291  printf("  LATTICE Statistics for today %s \n\n", asctime2(newtime));
4292  printf("  Number of constants: UDIC                 =%5ld"
4293         ", UDImax          =%5d\n",
4294         LINK->UDIC, UDImax);
4295  printf("  Number of keywords : nkw                  =%5ld"
4296         ", Lat_nkw_max     =%5d\n",
4297         LINK->nkw, Lat_nkw_max);
4298  printf("  Number of Families : globval.Elem_nFam    =%5ld"
4299         ", Elem_nFamMax    =%5d\n",
4300         globval.Elem_nFam, Elem_nFamMax);
4301  nKid = 0L;
4302  FORLIM = globval.Elem_nFam;
4303  for (j = 0L; j < FORLIM; j++) {
4304    if (ElemFam[j].nKid > nKid)
4305      nKid = ElemFam[j].nKid;
4306  }
4307  printf("  Max number of Kids : nKidMax              =%5ld"
4308         ", nKidMax         =%5d\n",
4309         nKid, nKidMax);
4310  printf("  Number of Blocks   : NoB                  =%5ld"
4311         ", NoBmax          =%5d\n",
4312         LINK->NoB, NoBmax);
4313  printf("  Max Block size     : NoBE                 =%5ld"
4314         ", NoBEmax         =%5d\n",
4315         LINK->Bpointer, NoBEmax);
4316  printf("  Number of Elements : globval.Cell_nLoc    =%5ld"
4317         ", Cell_nLocmax    =%5d\n",
4318         globval.Cell_nLoc, Cell_nLocMax);
4319  printf("  Circumference      : %12.7f [m]\n", Circumference(LINK));
4320  printf("\n");
4321  printf("\n");
4322}
4323
4324/****************************************************************************/
4325/* long ElemIndex(const char *name)
4326
4327   Purpose:
4328       return family index of the element .
4329       Note: in the PASCAL version the element family index could be
4330       comfortably accessed using the element name.
4331       This is no longer possible because we gave up on the interpretive PASCAL-S.
4332   Input:
4333       name    Family name
4334
4335   Output:
4336       none
4337
4338   Return:
4339       none
4340
4341   Global variables:
4342       none
4343
4344   Specific functions:
4345       none
4346
4347   Comments:
4348       30-06-2011   Fix the bug to get the correct Cell[i].Elem.PName.
4349****************************************************************************/
4350
4351long ElemIndex(const char *name)
4352{
4353  long       i = 0;
4354  int        n = 0;
4355  partsName  name1, name2;
4356
4357 
4358  const bool  prt = false;
4359
4360  if (prt) printf("\n");
4361
4362while(name[i]!= ' ' && name[i]!= '\0'){
4363  name1[i] = tolower(name[i]);
4364  i++;
4365  }
4366name1[i] = '\0'; 
4367
4368
4369  if (globval.Elem_nFam > Elem_nFamMax) {
4370    printf("ElemIndex: Elem_nFamMax exceeded: %ld(%d)\n",
4371           globval.Elem_nFam, Elem_nFamMax);
4372    exit_(1);
4373  }
4374
4375  for (i = 1; i <= globval.Elem_nFam; i++) 
4376  {
4377    n = 0;
4378    while ((ElemFam[i-1].ElemF.PName[n] != ' ') && (ElemFam[i-1].ElemF.PName[n] != '\0')) 
4379    {
4380      name2[n] = ElemFam[i-1].ElemF.PName[n]; 
4381      n++;
4382    }
4383    name2[n] = '\0';
4384
4385    if (prt)
4386      printf("%d %d %s |%s|%s|\n",
4387             n, strcmp(name2, name1) == 0, name1, ElemFam[i-1].ElemF.PName, name2);
4388    if (strcmp(name2, name1) == 0) 
4389      return i;
4390  }
4391
4392  printf("ElemIndex: undefined element %s\n", name);
4393  exit_(1);
4394  return 0;
4395}
4396
4397/****************************************************************************/
4398/* boolean Lattice_Read(FILE **fi_, FILE **fo_)
4399
4400   Purpose:
4401       Low level routine for reading lattice file
4402
4403   Input:
4404       none
4405
4406   Output:
4407       none
4408
4409   Return:
4410       none
4411
4412   Global variables:
4413       globval
4414
4415   Specific functions:
4416       setjmp f2c
4417       init_reserved_words
4418
4419   Comments:
4420       none
4421
4422****************************************************************************/
4423
4424bool Lattice_Read(FILE **fi_, FILE **fo_)
4425{
4426  struct LOC_Lattice_Read V;
4427
4428  if (trace)
4429    printf("t2lat: dndsym = %d, solsym = %d, max_set = %d, SETBITS = %u\n",
4430           bndsym, solsym, max_set, (unsigned)(4*sizeof(long int)));
4431
4432  V.fi = fi_; /* input lattice file */
4433  V.fo = fo_; /* output lattice file */
4434  if (setjmp(V._JL9999))
4435    goto _L9999;
4436  V.UDIC = 0; 
4437  globval.Cell_nLoc = 0; 
4438  globval.Elem_nFam = 0; 
4439  V.NoB = 0;
4440  V.Symmetry = 0; 
4441  V.Bpointer = 0;
4442
4443  globval.CODeps   = 0.0; 
4444  globval.dPcommon = 0.0; 
4445  globval.Energy   = 0.0;
4446
4447  ErrFlag = false;
4448
4449  init_reserved_words(&V);
4450
4451  GetSym___(&V);
4452
4453  if (V.sym == defsym) 
4454  DealWithDefns(&V);
4455
4456  if (V.Symmetry != 0) 
4457  {
4458    GetRingType(&V);/* define whether a ring or a transfer line */
4459    GetEnergy(&V);  /* define particle energy */
4460    GetCODEPS(&V);  /* define COD precision */
4461    GetDP(&V);      /* define energy offset */
4462  }
4463
4464  if (*V.fi != NULL)
4465    fclose(*V.fi);  /* Close lat file */
4466  *V.fi = NULL;
4467  if (*V.fo != NULL)
4468    fclose(*V.fo);  /* Close lax file */
4469  *V.fo = NULL;
4470  RegisterKids(&V); /* Check whether too many elements */
4471  PrintResult(&V);  /* Print lattice statistics */
4472 _L9999:
4473  return (!ErrFlag);
4474}
4475
4476#undef NoBmax
4477#undef NoBEmax
4478#undef UDImax
4479#undef LatLLng
4480#undef Lat_nkw_max
4481#undef emax
4482#undef emin
4483#undef kmax
4484#undef nmax
4485
4486#undef nn              // added by nsrl-ii
4487#undef tmax            // added by nsrl-ii
4488
4489#undef smax
4490#undef xmax
4491
Note: See TracBrowser for help on using the repository browser.