ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/trunk/smi23d/src/mengine/readprm.c
Revision: 25
Committed: Tue Jul 8 15:38:24 2008 UTC (13 years, 3 months ago) by tjod
File size: 118460 byte(s)
Log Message:
Incoporate the parameter files mmff94.prm and mmxconst.prm
into the mengine code by turning them into header files.
Each line becomes a string array element that is parsed in
place of lines in the .prm files.  The .prm files are no
longer needed, but are not deleted from the repository.

Line User Rev File contents
1 tjod 3 #define EXTERN extern
2    
3     #include "pcwin.h"
4     #include "pcmod.h"
5     #include "pot.h"
6     #include "field.h"
7     #include "atom_k.h"
8     #include "pdb.h"
9    
10     #define MAXOPLS 600
11    
12     FILE * fopen_path ( char * , char * , char * ) ;
13     void zero_data(void);
14 tjod 25 void readpid(char *);
15 tjod 3 void read_datafiles(char *);
16     void readmmxdata(char *);
17     void torphase(int, float *, int *);
18 tjod 25 void read_parameterfile(char **);
19 tjod 3 void numeral(int, char *, int);
20     void message_alert(char *, char *);
21     int FetchRecord(FILE *, char *);
22     void read_added_const(void);
23     int check_dipoleconst(int,int,float);
24     int check_torsionconst(int,int,int,int,float,float,float,int,int,int);
25     int check_torsion4const(int,int,int,int,float,float,float,int,int,int);
26     int check_torsion5const(int,int,int,int,float,float,float,int,int,int);
27     int check_angleconst(int,int,int,float,float,float,float);
28     int check_angle5const(int,int,int,float,float,float,float);
29     int check_angle4const(int,int,int,float,float,float,float);
30     int check_angle3const(int,int,int,float,float,float,float);
31     int check_bondconst(int,int,float,float);
32     int check_bond5const(int,int,float,float);
33     int check_bond4const(int,int,float,float);
34     int check_bond3const(int,int,float,float);
35     int check_opbend(int,int,int,int,float);
36     int check_strbnd(int,int,int,float,float,float);
37     int check_vdwpr(int,int,float,float);
38    
39     EXTERN struct t_minim_control {
40     int type, method, field, added_const;
41     char added_path[256],added_name[256];
42     } minim_control;
43     EXTERN struct t_minim_values {
44     int iprint, ndc, nconst;
45     float dielc;
46     } minim_values;
47    
48     struct t_units {
49     double bndunit, cbnd, qbnd;
50     double angunit, cang, qang, pang, sang, aaunit;
51     double stbnunit, ureyunit, torsunit, storunit, v14scale;
52     double aterm, bterm, cterm, dielec, chgscale;
53     } units;
54    
55     struct t_bondk1 {
56     int use_ring3, use_ring4, use_ring5;
57     int nbnd, nbnd3, nbnd4, nbnd5, ndeloc;
58     char kb[MAXBONDCONST][7], kb3[MAXBOND3CONST][7],kb4[MAXBOND4CONST][7],kb5[MAXBOND5CONST][7];
59     char kbdel[MAXBONDDELOC][7];
60     float s[MAXBONDCONST], t[MAXBONDCONST];
61     float s3[MAXBOND3CONST], t3[MAXBOND3CONST];
62     float s4[MAXBOND4CONST], t4[MAXBOND4CONST];
63     float s5[MAXBOND5CONST], t5[MAXBOND5CONST];
64     float sdel[MAXBONDDELOC], tdel[MAXBONDDELOC];
65     } bondk1;
66    
67     struct t_electroneg {
68     int nelecti, nelectj;
69     int itype[200],ibond[200],iattach[200];
70     int jtype[200],jbond[200],jattach[200];
71     float icorr[200],jcorr[200];
72     } electroneg;
73    
74     struct t_piatomk {
75     int npiatom;
76     int kat[MAXPIATOM], qp[MAXPIATOM];
77     float q[MAXPIATOM], ion[MAXPIATOM], emz[MAXPIATOM], zor[MAXPIATOM], zz[MAXPIATOM];
78     float w1[MAXPIATOM], w2[MAXPIATOM], w3[MAXPIATOM];
79     } piatomk;
80    
81     struct t_bondpk {
82     int npibond,npibond44i,npibond45i,npibond46i,npibond55i,npibond56i,npibond66i;
83     int npibond44o,npibond45o,npibond46o,npibond55o,npibond56o,npibond66o;
84     int npibond40,npibond50,npibond60;
85     int use_pibond,use_pibond44i,use_pibond45i,use_pibond46i,use_pibond55i,use_pibond56i,use_pibond66i;
86     int use_pibond44o,use_pibond45o,use_pibond46o,use_pibond55o,use_pibond56o,use_pibond66o;
87     int use_pibond40,use_pibond50,use_pibond60;
88     char kb[MAXPIBONDCONST][7],kb44i[20][7],kb45i[20][7],kb46i[20][7],kb55i[20][7],kb56i[20][7],kb66i[20][7];
89     char kb44o[20][7],kb45o[20][7],kb46o[20][7],kb55o[20][7],kb56o[20][7],kb66o[20][7];
90     char kb40[20][7],kb50[20][7],kb60[20][7];
91     float bk[MAXPIBONDCONST], bl[MAXPIBONDCONST], bmom[MAXPIBONDCONST],sslop[MAXPIBONDCONST], tslop[MAXPIBONDCONST], tslop2[MAXPIBONDCONST];
92     float bk44i[20], bl44i[20], bmom44i[20],sslop44i[20], tslop44i[20], tslop244i[20];
93     float bk45i[20], bl45i[20], bmom45i[20],sslop45i[20], tslop45i[20], tslop245i[20];
94     float bk46i[20], bl46i[20], bmom46i[20],sslop46i[20], tslop46i[20], tslop246i[20];
95     float bk55i[20], bl55i[20], bmom55i[20],sslop55i[20], tslop55i[20], tslop255i[20];
96     float bk56i[20], bl56i[20], bmom56i[20],sslop56i[20], tslop56i[20], tslop256i[20];
97     float bk66i[20], bl66i[20], bmom66i[20],sslop66i[20], tslop66i[20], tslop266i[20];
98     float bk44o[20], bl44o[20], bmom44o[20],sslop44o[20], tslop44o[20], tslop244o[20];
99     float bk45o[20], bl45o[20], bmom45o[20],sslop45o[20], tslop45o[20], tslop245o[20];
100     float bk46o[20], bl46o[20], bmom46o[20],sslop46o[20], tslop46o[20], tslop246o[20];
101     float bk55o[20], bl55o[20], bmom55o[20],sslop55o[20], tslop55o[20], tslop255o[20];
102     float bk56o[20], bl56o[20], bmom56o[20],sslop56o[20], tslop56o[20], tslop256o[20];
103     float bk66o[20], bl66o[20], bmom66o[20],sslop66o[20], tslop66o[20], tslop266o[20];
104     float bk40[20], bl40[20], bmom40[20],sslop40[20], tslop40[20], tslop240[20];
105     float bk50[20], bl50[20], bmom50[20],sslop50[20], tslop50[20], tslop250[20];
106     float bk60[20], bl60[20], bmom60[20],sslop60[20], tslop60[20], tslop260[20];
107     } bondpk;
108    
109     struct t_picorrect {
110     float fudge1, fudge2, fudge3, fudge4, beta_fudgef,beta_fudgeh;
111     } picorrect;
112    
113     struct t_angk1 {
114     int use_ang3, use_ang4, use_ang5;
115     int nang, nang3, nang4, nang5;
116     int ndel, ndel3, ndel4;
117     char ktype[MAXANGCONST][10],ktype3[MAXANG3CONST][10],ktype4[MAXANG4CONST][10],ktype5[MAXANG5CONST][10];
118     char kdel[MAXANGDEL][10],kdel3[MAXANG3DEL][10],kdel4[MAXANG4DEL][10];
119     float con[MAXANGCONST], ang[MAXANGCONST][3];
120     float con3[MAXANG3CONST], ang3[MAXANG3CONST][3];
121     float con4[MAXANG4CONST], ang4[MAXANG4CONST][3];
122     float con5[MAXANG5CONST], ang5[MAXANG5CONST][3];
123     float condel[MAXANGDEL], angdel[MAXANGDEL][3];
124     float condel3[MAXANG3DEL], angdel3[MAXANG3DEL][3];
125     float condel4[MAXANG4DEL], angdel4[MAXANG4DEL][3];
126     } angk1;
127    
128     struct t_angkp1 {
129     int npiang;
130     char kpa[MAXPIANGCONST][10];
131     float pacon[MAXPIANGCONST], panat[MAXPIANGCONST][3], piobp[MAXPIANGCONST];
132     } angkp1;
133    
134     struct t_angf {
135     int use_angf, nfang;
136     char kftype[MAXANGCONST][10];
137     float fcon[MAXANGCONST],fc0[MAXANGCONST],fc1[MAXANGCONST],fc2[MAXANGCONST];
138     float fc3[MAXANGCONST],fc4[MAXANGCONST],fc5[MAXANGCONST],fc6[MAXANGCONST];
139     } angf;
140    
141     struct t_ureybrad_k {
142     int nurey_brad;
143     char kubang[MAXUREY][10];
144     float ubconst[MAXUREY],ubdist[MAXUREY];
145     } ureybrad_k;
146    
147     struct t_crossterm_k {
148     int nangang, nstrbnd, nstrtor;
149     int ang_ang[MAXAA], stbnindex[MAXSTBN] ;
150     char str_tor[MAXSTRTOR][7], stbn[MAXSTBN][10];
151     float aacon[MAXAA][3], stbncon[MAXSTBN][3], str_torcon[MAXSTRTOR];
152     } crossterm_k;
153    
154     struct t_ooplane_k {
155     int nopbend;
156     char iopb[MAXOOP][13];
157     float copb[MAXOOP];
158     } ooplane_k;
159    
160     struct t_improptor_k {
161     int nimptor, nimprop;
162     char kv[MAXIMP][13];
163     float cimptor[MAXIMP], tdi[MAXIMP];
164     float v1[MAXIMP], v2[MAXIMP], v3[MAXIMP];
165     int ph1[MAXIMP], ph2[MAXIMP], ph3[MAXIMP];
166     } improptor_k;
167    
168     struct t_torkn1 {
169     int use_tor4, use_tor5;
170     int ntor, ntor4, ntor5, ntordel, torindex[MAXTORDEL];
171     char kv[MAXTORCONST][13], kv4[MAXTOR4CONST][13], kv5[MAXTOR5CONST][13];
172     char kvdel[MAXTORDEL][13];
173     float tv1[MAXTORCONST], tv2[MAXTORCONST], tv3[MAXTORCONST];
174     float tv4[MAXTORCONST], tv5[MAXTORCONST], tv6[MAXTORCONST];
175     int phase1[MAXTORCONST], phase2[MAXTORCONST], phase3[MAXTORCONST];
176     int phase4[MAXTORCONST], phase5[MAXTORCONST], phase6[MAXTORCONST];
177     float tv41[MAXTOR4CONST], tv42[MAXTOR4CONST], tv43[MAXTOR4CONST];
178     int phase41[MAXTOR4CONST], phase42[MAXTOR4CONST], phase43[MAXTOR4CONST];
179     float tv51[MAXTOR5CONST], tv52[MAXTOR5CONST], tv53[MAXTOR5CONST];
180     int phase51[MAXTOR5CONST], phase52[MAXTOR5CONST], phase53[MAXTOR5CONST];
181     float tvdel1[MAXTORDEL], tvdel2[MAXTORDEL], tvdel3[MAXTORDEL];
182     int phasedel1[MAXTORDEL], phasedel2[MAXTORDEL], phasedel3[MAXTORDEL];
183     } torkn1;
184    
185     struct t_torknp {
186     int npitor;
187     char kv[MAXPITORCONST][13];
188     int ph1[MAXPITORCONST], ph2[MAXPITORCONST], ph3[MAXPITORCONST];
189     float tv1[MAXPITORCONST], tv2[MAXPITORCONST], tv3[MAXPITORCONST];
190     } torknp;
191    
192     struct t_vdw1 {
193     int nvdw;
194     float rad[MAXVDWCONST], eps[MAXVDWCONST];
195     int lpd[MAXVDWCONST], ihtyp[MAXVDWCONST], ihdon[MAXVDWCONST];
196     float alpha[MAXVDWCONST],n[MAXVDWCONST],a[MAXVDWCONST],g[MAXVDWCONST];
197     char da[MAXVDWCONST][2];
198     } vdw1;
199    
200     struct t_vdwpr_k {
201     int nvdwpr;
202     int ia1[MAXBONDCONST],ia2[MAXBONDCONST];
203     char kv[MAXBONDCONST][7];
204     float radius[MAXBONDCONST], eps[MAXBONDCONST];
205     } vdwpr_k;
206    
207     struct t_charge_k {
208     int ncharge, nbndchrg, nbndchrgdel;
209     int type[MAXATOMTYPE], btype[MAXBONDCONST], btypedel[MAXBONDCONST];
210     float charge[MAXATOMTYPE], bcharge[MAXBONDCONST], formchrg[MAXATOMTYPE], bchargedel[MAXBONDCONST];
211     float typechrg[MAXATOMTYPE];
212     } charge_k;
213    
214     struct t_amberchrg_k {
215     int ncharge, type[MAXATOMTYPE], res_type[MAXATOMTYPE];
216     char symbol[MAXATOMTYPE][4];
217     float chrg[MAXATOMTYPE];
218     } amberchrg_k;
219    
220     struct t_oplschrg_k {
221     int ncharge, type[MAXOPLS];
222     char chrgstring[MAXOPLS][25];
223     float chrg[MAXOPLS];
224     } oplschrg_k;
225    
226     struct t_dipole_k {
227     int ndipole,ndipole3,ndipole4,ndipole5;
228     char kb[MAXBONDCONST][7], kb3[MAXBOND3CONST][7], kb4[MAXBOND4CONST][7], kb5[MAXBOND5CONST][7];
229     float bmom[MAXBONDCONST],bmom3[MAXBOND3CONST],bmom4[MAXBOND4CONST],bmom5[MAXBOND5CONST];
230     } dipole_k;
231    
232     struct t_ehpara {
233     int neheat, nevhf;
234     int nn[220];
235     float ee[220], estr[220], et[220]; // bond data
236     char cc[220][21], mm[155][21];
237     int cent[155];
238     char kheat[155][13];
239     float ss[155], ssless[155]; // structural data
240     } ehpara;
241    
242     struct t_epiheat {
243     int npihf;
244     int nnp[125];
245     float eep[125], aa[125], bb[125], ccc[125];
246     char ccp[125][21];} epiheat;
247    
248     struct t_metaldata {
249     int nmetal, type[200];
250     char name[200][3];
251     float radius[200], eps[200];
252     } metaldata;
253    
254     EXTERN struct t_user {
255     int dielec;
256     } user;
257    
258     void zero_data()
259     {
260     int ii,j;
261    
262     bondk1.use_ring3 = FALSE;
263     bondk1.use_ring4 = FALSE;
264     bondk1.use_ring5 = FALSE;
265     bondpk.use_pibond = FALSE;
266     bondpk.use_pibond44i = FALSE;
267     bondpk.use_pibond45i = FALSE;
268     bondpk.use_pibond46i = FALSE;
269     bondpk.use_pibond55i = FALSE;
270     bondpk.use_pibond56i = FALSE;
271     bondpk.use_pibond66i = FALSE;
272     bondpk.use_pibond44o = FALSE;
273     bondpk.use_pibond45o = FALSE;
274     bondpk.use_pibond46o = FALSE;
275     bondpk.use_pibond55o = FALSE;
276     bondpk.use_pibond56o = FALSE;
277     bondpk.use_pibond66o = FALSE;
278     bondpk.use_pibond40 = FALSE;
279     bondpk.use_pibond50 = FALSE;
280     bondpk.use_pibond60 = FALSE;
281    
282     angk1.use_ang3 = FALSE;
283     angk1.use_ang4 = FALSE;
284     angk1.use_ang5 = FALSE;
285     angf.use_angf = FALSE;
286    
287     torkn1.use_tor4 = FALSE;
288     torkn1.use_tor5 = FALSE;
289    
290     units.bndunit = 1.0;
291     units.cbnd = 0.0;
292     units.qbnd = 0.0;
293     units.angunit = 1.0 / (radian*radian);
294     units.cang = 0.0;
295     units.qang = 0.0;
296     units.pang = 0.0;
297     units.sang = 0.0;
298     units.aaunit = 1.0 / (radian*radian);
299     units.stbnunit = 1.0;
300     units.ureyunit = 1.0;
301     units.torsunit = 1.0;
302     units.storunit = 1.0;
303     units.v14scale = 1.0;
304     units.aterm = 0.0;
305     units.bterm = 0.0;
306     units.cterm = 0.0;
307     // units.dielec = 1.0;
308     units.chgscale = 1.0;
309     // atom type constants
310     atom_k.natomtype = 0;
311     charge_k.ncharge = 0;
312     charge_k.nbndchrg = 0;
313     charge_k.nbndchrgdel = 0;
314     amberchrg_k.ncharge = 0;
315     oplschrg_k.ncharge = 0;
316     for (ii = 0; ii < MAXATOMTYPE; ii ++)
317     {
318     atom_k.type[ii]= 0;
319     atom_k.valency[ii]= 0;
320     atom_k.tclass[ii] = 0;
321     atom_k.tclass1[ii] = 0;
322     atom_k.tclass2[ii] = 0;
323     atom_k.number[ii] = 0;
324     atom_k.ligands[ii] = 0;
325     atom_k.weight[ii] = 0.00F;
326     strcpy(atom_k.symbol[ii]," ");
327     strcpy(atom_k.description[ii]," ");
328     charge_k.type[ii] = 0;
329     charge_k.charge[ii] = 0.00F;
330     charge_k.formchrg[ii] = 0.00F;
331     charge_k.typechrg[ii] = 0.00F;
332     amberchrg_k.type[ii] = 0;
333     amberchrg_k.res_type[ii] = 0;
334     amberchrg_k.chrg[ii] = 0.0;
335     strcpy(amberchrg_k.symbol[ii]," ");
336     }
337     for (ii = 0; ii < MAXOPLS; ii++)
338     {
339     oplschrg_k.type[ii] = 0;
340     oplschrg_k.chrg[ii] = 0.0;
341     strcpy(oplschrg_k.chrgstring[ii]," ");
342     }
343     // metals
344     metaldata.nmetal = 0;
345     for (ii = 0; ii < 200; ii++)
346     {
347     metaldata.type[ii] = 0;
348     metaldata.radius[ii] = 0.0;
349     metaldata.eps[ii] = 0.0;
350     strcpy(metaldata.name[ii]," ");
351     }
352     // bond constant data
353     bondk1.nbnd = 0;
354     vdwpr_k.nvdwpr = 0;
355     dipole_k.ndipole = 0;
356     electroneg.nelecti = 0;
357     electroneg.nelectj = 0;
358    
359     for( ii = 0; ii < MAXBONDCONST; ii++ )
360     {
361     strcpy(bondk1.kb[ii]," ");
362     bondk1.s[ii] = 0.;
363     bondk1.t[ii] = 0.;
364     strcpy(vdwpr_k.kv[ii]," ");
365     vdwpr_k.ia1[ii] = 0;
366     vdwpr_k.ia2[ii] = 0;
367     vdwpr_k.radius[ii] = 0.0;
368     vdwpr_k.eps[ii] = 0.0;
369     strcpy(dipole_k.kb[ii]," ");
370     dipole_k.bmom[ii] = 0.;
371     charge_k.btype[ii] = 0;
372     charge_k.bcharge[ii] = 0.0;
373     charge_k.btypedel[ii] = 0;
374     charge_k.bchargedel[ii] = 0.0;
375     }
376     // bond 3 constant
377     bondk1.nbnd3 = 0;
378     dipole_k.ndipole3 = 0;
379     for( ii = 0; ii < MAXBOND3CONST; ii++ )
380     {
381     strcpy(bondk1.kb3[ii]," ");
382     bondk1.s3[ii] = 0.;
383     bondk1.t3[ii] = 0.;
384     strcpy(dipole_k.kb3[ii]," ");
385     dipole_k.bmom3[ii] = 0.;
386     }
387     // bond 4 constant
388     bondk1.nbnd4 = 0;
389     dipole_k.ndipole4 = 0;
390     for( ii = 0; ii < MAXBOND4CONST; ii++ )
391     {
392     strcpy(bondk1.kb4[ii]," ");
393     bondk1.s4[ii] = 0.;
394     bondk1.t4[ii] = 0.;
395     strcpy(dipole_k.kb4[ii]," ");
396     dipole_k.bmom4[ii] = 0.;
397     }
398     // bond 5 constant
399     bondk1.nbnd5 = 0;
400     dipole_k.ndipole5 = 0;
401     for( ii = 0; ii < MAXBOND5CONST; ii++ )
402     {
403     strcpy(bondk1.kb5[ii]," ");
404     bondk1.s5[ii] = 0.;
405     bondk1.t5[ii] = 0.;
406     strcpy(dipole_k.kb5[ii]," ");
407     dipole_k.bmom5[ii] = 0.;
408     }
409     // deloc bond constant
410     bondk1.ndeloc = 0;
411     for( ii = 0; ii < MAXBONDDELOC; ii++ )
412     {
413     strcpy(bondk1.kbdel[ii]," ");
414     bondk1.sdel[ii] = 0.;
415     bondk1.tdel[ii] = 0.;
416     }
417     // VDW constants
418     vdw1.nvdw = 0;
419     for( ii = 0; ii < MAXVDWCONST; ii++ )
420     {
421     vdw1.rad[ii] = 0.;
422     vdw1.eps[ii] = 0.;
423     vdw1.lpd[ii] = 0;
424     vdw1.ihtyp[ii] = 0;
425     vdw1.ihdon[ii] = 0;
426     vdw1.alpha[ii] = 0.0;
427     vdw1.n[ii] = 0.0;
428     vdw1.a[ii] = 0.0;
429     vdw1.g[ii] = 0.0;
430     strcpy(vdw1.da[ii]," ");
431     }
432     // torsion
433     torkn1.ntor = 0;
434     for( ii = 0; ii < MAXTORCONST; ii++ )
435     {
436     strcpy(torkn1.kv[ii]," ");
437     torkn1.tv1[ii] = 0.;
438     torkn1.tv2[ii] = 0.;
439     torkn1.tv3[ii] = 0.;
440     torkn1.tv4[ii] = 0.;
441     torkn1.tv5[ii] = 0.;
442     torkn1.tv6[ii] = 0.;
443     torkn1.phase1[ii] = torkn1.phase2[ii] = torkn1.phase3[ii] = 0;
444     torkn1.phase4[ii] = torkn1.phase5[ii] = torkn1.phase6[ii] = 0;
445     }
446     // torsion 4
447     torkn1.ntor4 = 0;
448     for( ii = 0; ii < MAXTOR4CONST; ii++ )
449     {
450     strcpy(torkn1.kv4[ii]," ");
451     torkn1.tv41[ii] = 0.;
452     torkn1.tv42[ii] = 0.;
453     torkn1.tv43[ii] = 0.;
454     torkn1.phase41[ii] = torkn1.phase42[ii] = torkn1.phase43[ii] = 0;
455     }
456     // torsion 5
457     torkn1.ntor5 = 0;
458     for( ii = 0; ii < MAXTOR5CONST; ii++ )
459     {
460     strcpy(torkn1.kv5[ii]," ");
461     torkn1.tv51[ii] = 0.;
462     torkn1.tv52[ii] = 0.;
463     torkn1.tv53[ii] = 0.;
464     torkn1.phase51[ii] = torkn1.phase52[ii] = torkn1.phase53[ii] = 0;
465     }
466     // torsion delocalized
467     torkn1.ntordel = 0;
468     for( ii = 0; ii < MAXTORDEL; ii++ )
469     {
470     strcpy(torkn1.kvdel[ii]," ");
471     torkn1.torindex[ii] = 0;
472     torkn1.tvdel1[ii] = 0.;
473     torkn1.tvdel2[ii] = 0.;
474     torkn1.tvdel3[ii] = 0.;
475     torkn1.phasedel1[ii] = torkn1.phasedel2[ii] = torkn1.phasedel3[ii] = 0;
476     }
477     // fourier angles
478     angf.nfang = 0;
479     for (ii=0; ii < MAXANGCONST; ii++ )
480     {
481     strcpy(angf.kftype[ii]," ");
482     angf.fcon[ii] = 0.0;
483     angf.fc0[ii] = 0.0;
484     angf.fc1[ii] = 0.0;
485     angf.fc2[ii] = 0.0;
486     angf.fc3[ii] = 0.0;
487     angf.fc4[ii] = 0.0;
488     angf.fc5[ii] = 0.0;
489     angf.fc6[ii] = 0.0;
490     }
491     // angle
492     angk1.nang = 0;
493     for( ii = 0; ii < MAXANGCONST; ii++ )
494     {
495     strcpy(angk1.ktype[ii]," ");
496     angk1.con[ii] = 0.;
497     for( j = 0; j < 3; j++ )
498     angk1.ang[ii][j] = 0.;
499     }
500     // angle 3
501     angk1.nang3 = 0;
502     for( ii = 0; ii < MAXANG3CONST; ii++ )
503     {
504     strcpy(angk1.ktype3[ii]," ");
505     angk1.con3[ii] = 0.;
506     for( j = 0; j < 3; j++ )
507     angk1.ang3[ii][j] = 0.;
508     }
509     // angle 4
510     angk1.nang4 = 0;
511     for( ii = 0; ii < MAXANG4CONST; ii++ )
512     {
513     strcpy(angk1.ktype4[ii]," ");
514     angk1.con4[ii] = 0.;
515     for( j = 0; j < 3; j++ )
516     angk1.ang4[ii][j] = 0.;
517     }
518     // angle 5
519     angk1.nang5 = 0;
520     for( ii = 0; ii < MAXANG5CONST; ii++ )
521     {
522     strcpy(angk1.ktype5[ii]," ");
523     angk1.con5[ii] = 0.;
524     for( j = 0; j < 3; j++ )
525     angk1.ang5[ii][j] = 0.;
526     }
527     // angle delocalized
528     angk1.ndel = 0;
529     for( ii = 0; ii < MAXANGDEL; ii++ )
530     {
531     strcpy(angk1.kdel[ii]," ");
532     angk1.condel[ii] = 0.;
533     for( j = 0; j < 3; j++ )
534     angk1.angdel[ii][j] = 0.;
535     }
536     // angle3 delocalized
537     angk1.ndel3 = 0;
538     for( ii = 0; ii < MAXANG3DEL; ii++ )
539     {
540     strcpy(angk1.kdel3[ii]," ");
541     angk1.condel3[ii] = 0.;
542     for( j = 0; j < 3; j++ )
543     angk1.angdel3[ii][j] = 0.;
544     }
545     // angle4 delocalized
546     angk1.ndel4 = 0;
547     for( ii = 0; ii < MAXANG4DEL; ii++ )
548     {
549     strcpy(angk1.kdel4[ii]," ");
550     angk1.condel4[ii] = 0.;
551     for( j = 0; j < 3; j++ )
552     angk1.angdel4[ii][j] = 0.;
553     }
554     // urey bradley
555     ureybrad_k.nurey_brad = 0;
556     for (ii = 0; ii < MAXUREY; ii++)
557     {
558     strcpy(ureybrad_k.kubang[ii]," ");
559     ureybrad_k.ubconst[ii] = 0.0;
560     ureybrad_k.ubdist[ii] = 0.0;
561     }
562     // out of plane
563     ooplane_k.nopbend = 0;
564     for( ii = 0; ii < MAXOOP; ii++ )
565     {
566     strcpy(ooplane_k.iopb[ii]," ");
567     ooplane_k.copb[ii] = 0.;
568     }
569     // improper torsions
570     improptor_k.nimptor = 0;
571     improptor_k.nimprop = 0;
572     for( ii = 0; ii < MAXIMP; ii++ )
573     {
574     strcpy(improptor_k.kv[ii]," ");
575     improptor_k.cimptor[ii] = 0.;
576     improptor_k.tdi[ii] = 0.;
577     improptor_k.v1[ii] = 0.0;
578     improptor_k.v2[ii] = 0.0;
579     improptor_k.v3[ii] = 0.0;
580     improptor_k.ph1[ii] = 0;
581     improptor_k.ph2[ii] = 0;
582     improptor_k.ph3[ii] = 0;
583     }
584     // cross terms
585     crossterm_k.nangang = 0;
586     crossterm_k.nstrbnd = 0;
587     crossterm_k.nstrtor = 0;
588     for (ii = 0; ii < MAXAA; ii++)
589     {
590     crossterm_k.ang_ang[ii] = 0;
591     crossterm_k.aacon[ii][0] = 0.0;
592     crossterm_k.aacon[ii][1] = 0.0;
593     crossterm_k.aacon[ii][2] = 0.0;
594     }
595     for (ii = 0; ii < MAXSTBN; ii++)
596     {
597     strcpy(crossterm_k.stbn[ii]," ");
598     crossterm_k.stbnindex[ii] = 0;
599     crossterm_k.stbncon[ii][0] = 0.0;
600     crossterm_k.stbncon[ii][1] = 0.0;
601     crossterm_k.stbncon[ii][2] = 0.0;
602     }
603     for (ii = 0; ii < MAXSTRTOR; ii++)
604     {
605     strcpy(crossterm_k.str_tor[ii]," ");
606     crossterm_k.str_torcon[ii] = 0.0;
607     }
608    
609    
610     // heat parameters
611     ehpara.neheat = 0;
612     ehpara.nevhf = 0;
613     for( ii = 0; ii < 220; ii++ )
614     {
615     ehpara.nn[ii] = 0;
616     ehpara.ee[ii] = 0.;
617     ehpara.estr[ii] = 0.;
618     ehpara.et[ii] = 0.;
619     strcpy( ehpara.cc[ii], " " );
620     }
621     for( ii = 0; ii < 155; ii++ )
622     {
623     strcpy( ehpara.mm[ii], " " );
624     strcpy( ehpara.kheat[ii], " " );
625     ehpara.ss[ii] = 0.;
626     ehpara.ssless[ii] = 0.;
627     ehpara.cent[ii] = 0;
628     }
629     // pi heat parameters
630     epiheat.npihf = 0;
631     for( ii = 0; ii < 125; ii++ )
632     {
633     epiheat.nnp[ii] = 0;
634     epiheat.eep[ii] = 0.;
635     epiheat.aa[ii] = 0.;
636     epiheat.bb[ii] = 0.;
637     epiheat.ccc[ii] = 0.;
638     strcpy( epiheat.ccp[ii], " " );
639     }
640     // pi atom data
641     piatomk.npiatom = 0;
642     for (ii=0; ii < MAXPIATOM; ii++)
643     {
644     piatomk.kat[ii] = 0;
645     piatomk.qp[ii] = 0;
646     piatomk.q[ii] = 0.0;
647     piatomk.ion[ii] = 0.0;
648     piatomk.emz[ii] = 0.0;
649     piatomk.zor[ii] = 0.0;
650     piatomk.zz[ii] = 0.0;
651     piatomk.w1[ii] = 0.0;
652     piatomk.w2[ii] = 0.0;
653     piatomk.w3[ii] = 0.0;
654     }
655     // pi bond constants
656     bondpk.npibond = 0;
657     bondpk.npibond44i = 0;
658     bondpk.npibond45i = 0;
659     bondpk.npibond46i = 0;
660     bondpk.npibond55i = 0;
661     bondpk.npibond56i = 0;
662     bondpk.npibond66i = 0;
663     bondpk.npibond44o = 0;
664     bondpk.npibond45o = 0;
665     bondpk.npibond46o = 0;
666     bondpk.npibond55o = 0;
667     bondpk.npibond56o = 0;
668     bondpk.npibond66o = 0;
669     bondpk.npibond40 = 0;
670     bondpk.npibond50 = 0;
671     bondpk.npibond60 = 0;
672     for (ii = 0; ii < 20; ii++)
673     {
674     strcpy(bondpk.kb44i[ii]," ");
675     strcpy(bondpk.kb45i[ii]," ");
676     strcpy(bondpk.kb46i[ii]," ");
677     strcpy(bondpk.kb55i[ii]," ");
678     strcpy(bondpk.kb56i[ii]," ");
679     strcpy(bondpk.kb66i[ii]," ");
680     strcpy(bondpk.kb44o[ii]," ");
681     strcpy(bondpk.kb45o[ii]," ");
682     strcpy(bondpk.kb46o[ii]," ");
683     strcpy(bondpk.kb55o[ii]," ");
684     strcpy(bondpk.kb56o[ii]," ");
685     strcpy(bondpk.kb66o[ii]," ");
686     strcpy(bondpk.kb40[ii]," ");
687     strcpy(bondpk.kb50[ii]," ");
688     strcpy(bondpk.kb60[ii]," ");
689    
690     bondpk.bk44i[ii] = 0.0;
691     bondpk.bl44i[ii] = 0.0;
692     bondpk.bmom44i[ii] = 0.0;
693     bondpk.sslop44i[ii] = 0.0;
694     bondpk.tslop44i[ii] = 0.0;
695     bondpk.tslop244i[ii] = 0.0;
696     bondpk.bk45i[ii] = 0.0;
697     bondpk.bl45i[ii] = 0.0;
698     bondpk.bmom45i[ii] = 0.0;
699     bondpk.sslop45i[ii] = 0.0;
700     bondpk.tslop45i[ii] = 0.0;
701     bondpk.tslop245i[ii] = 0.0;
702     bondpk.bk46i[ii] = 0.0;
703     bondpk.bl46i[ii] = 0.0;
704     bondpk.bmom46i[ii] = 0.0;
705     bondpk.sslop46i[ii] = 0.0;
706     bondpk.tslop46i[ii] = 0.0;
707     bondpk.tslop246i[ii] = 0.0;
708     bondpk.bk55i[ii] = 0.0;
709     bondpk.bl55i[ii] = 0.0;
710     bondpk.bmom55i[ii] = 0.0;
711     bondpk.sslop55i[ii] = 0.0;
712     bondpk.tslop55i[ii] = 0.0;
713     bondpk.tslop255i[ii] = 0.0;
714     bondpk.bk56i[ii] = 0.0;
715     bondpk.bl56i[ii] = 0.0;
716     bondpk.bmom56i[ii] = 0.0;
717     bondpk.sslop56i[ii] = 0.0;
718     bondpk.tslop56i[ii] = 0.0;
719     bondpk.tslop256i[ii] = 0.0;
720     bondpk.bk66i[ii] = 0.0;
721     bondpk.bl66i[ii] = 0.0;
722     bondpk.bmom66i[ii] = 0.0;
723     bondpk.sslop66i[ii] = 0.0;
724     bondpk.tslop66i[ii] = 0.0;
725     bondpk.tslop266i[ii] = 0.0;
726     bondpk.bk44o[ii] = 0.0;
727     bondpk.bl44o[ii] = 0.0;
728     bondpk.bmom44o[ii] = 0.0;
729     bondpk.sslop44o[ii] = 0.0;
730     bondpk.tslop44o[ii] = 0.0;
731     bondpk.tslop244o[ii] = 0.0;
732     bondpk.bk45o[ii] = 0.0;
733     bondpk.bl45o[ii] = 0.0;
734     bondpk.bmom45o[ii] = 0.0;
735     bondpk.sslop45o[ii] = 0.0;
736     bondpk.tslop45o[ii] = 0.0;
737     bondpk.tslop245o[ii] = 0.0;
738     bondpk.bk46o[ii] = 0.0;
739     bondpk.bl46o[ii] = 0.0;
740     bondpk.bmom46o[ii] = 0.0;
741     bondpk.sslop46o[ii] = 0.0;
742     bondpk.tslop46o[ii] = 0.0;
743     bondpk.tslop246o[ii] = 0.0;
744     bondpk.bk55o[ii] = 0.0;
745     bondpk.bl55o[ii] = 0.0;
746     bondpk.bmom55o[ii] = 0.0;
747     bondpk.sslop55o[ii] = 0.0;
748     bondpk.tslop55o[ii] = 0.0;
749     bondpk.tslop255o[ii] = 0.0;
750     bondpk.bk56o[ii] = 0.0;
751     bondpk.bl56o[ii] = 0.0;
752     bondpk.bmom56o[ii] = 0.0;
753     bondpk.sslop56o[ii] = 0.0;
754     bondpk.tslop56o[ii] = 0.0;
755     bondpk.tslop256o[ii] = 0.0;
756     bondpk.bk66o[ii] = 0.0;
757     bondpk.bl66o[ii] = 0.0;
758     bondpk.bmom66o[ii] = 0.0;
759     bondpk.sslop66o[ii] = 0.0;
760     bondpk.tslop66o[ii] = 0.0;
761     bondpk.tslop266o[ii] = 0.0;
762     bondpk.bk40[ii] = 0.0;
763     bondpk.bl40[ii] = 0.0;
764     bondpk.bmom40[ii] = 0.0;
765     bondpk.sslop40[ii] = 0.0;
766     bondpk.tslop40[ii] = 0.0;
767     bondpk.tslop240[ii] = 0.0;
768     bondpk.bk50[ii] = 0.0;
769     bondpk.bl50[ii] = 0.0;
770     bondpk.bmom50[ii] = 0.0;
771     bondpk.sslop50[ii] = 0.0;
772     bondpk.tslop50[ii] = 0.0;
773     bondpk.tslop250[ii] = 0.0;
774     bondpk.bk60[ii] = 0.0;
775     bondpk.bl60[ii] = 0.0;
776     bondpk.bmom60[ii] = 0.0;
777     bondpk.sslop60[ii] = 0.0;
778     bondpk.tslop60[ii] = 0.0;
779     bondpk.tslop260[ii] = 0.0;
780     }
781     for (ii = 0; ii < MAXPIBONDCONST; ii++)
782     {
783     strcpy(bondpk.kb[ii]," ");
784     bondpk.bk[ii] = 0.0;
785     bondpk.bl[ii] = 0.0;
786     bondpk.bmom[ii] = 0.0;
787     bondpk.sslop[ii] = 0.0;
788     bondpk.tslop[ii] = 0.0;
789     bondpk.tslop2[ii] = 0.0;
790     }
791     picorrect.fudge1 = 0.0;
792     picorrect.fudge2 = 0.0;
793     picorrect.fudge3 = 0.0;
794     picorrect.fudge4 = 0.0;
795     picorrect.beta_fudgef = 0.0;
796     picorrect.beta_fudgeh = 0.0;
797    
798     // pi torsions
799     torknp.npitor = 0;
800     for (ii = 0; ii < MAXPITORCONST; ii++)
801     {
802     strcpy(torknp.kv[ii]," ");
803     torknp.ph1[ii] = 0;
804     torknp.ph2[ii] = 0;
805     torknp.ph3[ii] = 0;
806     torknp.tv1[ii] = 0.0;
807     torknp.tv2[ii] = 0.0;
808     torknp.tv3[ii] = 0.0;
809     }
810     // pi angles
811     angkp1.npiang = 0;
812     for (ii = 0; ii < MAXPIANGCONST; ii++)
813     {
814     strcpy(angkp1.kpa[ii]," ");
815     angkp1.pacon[ii] = 0.0;
816     angkp1.panat[ii][0] = 0.0;
817     angkp1.panat[ii][1] = 0.0;
818     angkp1.panat[ii][2] = 0.0;
819     angkp1.piobp[ii] = 0.0;
820     }
821    
822     }
823     /* =========================================== */
824     void get_added_const()
825     {
826     int ia,ib,ic,id, ip1,ip2,ip3,iz,if3,if4,if5;
827     float f1,f2,f3,f4;
828     float v1[6];
829     int se1[6];
830     char pa[3],pb[3],pc[3],pd[3],pt[13];
831     char line[151], iptemp[21], dummy[10];
832     FILE *wfile;
833    
834     wfile = fopen_path(minim_control.added_path,minim_control.added_name,"r");
835     if (wfile == NULL)
836     {
837     message_alert("Error opening added constants file","Error");
838     fprintf(pcmoutfile,"Error reading added constants file: %s\n",minim_control.added_name);
839     return;
840     }
841    
842     fprintf(pcmoutfile,"\nThe Following Added Parameters were read:\n\n");
843     while ( FetchRecord(wfile,line) )
844     {
845     sscanf(line,"%s",iptemp);
846    
847     if (strcmp(iptemp,"bond") == 0)
848     {
849     sscanf(line,"%s %d %d %f %f",dummy,&ia,&ib,&f1,&f2);
850     fprintf(pcmoutfile,"Bond: AtomTypes %d %d Blen %f Bk %f \n",ia,ib,f1,f2);
851     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE)
852     {
853     iz = check_bondconst(ia,ib,f1,f2);
854     if (iz == FALSE)
855     {
856     bondk1.s[bondk1.nbnd] = f1;
857     bondk1.t[bondk1.nbnd] = f2;
858     numeral(ia,pa,3);
859     numeral(ib,pb,3);
860     strcpy(pt,pa);
861     strcat(pt,pb);
862     strcpy(bondk1.kb[bondk1.nbnd],pt);
863     bondk1.nbnd++;
864     }
865     } else
866     {
867     message_alert("Error reading bond atom types in added constants file","Error");
868     }
869     } else if (strcmp(iptemp,"bond5") == 0)
870     {
871     iz = sscanf(line,"%s %d %d %f %f",dummy,&ia,&ib,&f1,&f2);
872     fprintf(pcmoutfile,"Bond5: AtomTypes %d %d Blen %f Bk %f \n",ia,ib,f1,f2);
873     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && iz > 3)
874     {
875     iz = check_bond5const(ia,ib,f1,f2);
876     if (iz == FALSE)
877     {
878     bondk1.s5[bondk1.nbnd5] = f1;
879     bondk1.t5[bondk1.nbnd5] = f2;
880     numeral(ia,pa,3);
881     numeral(ib,pb,3);
882     strcpy(pt,pa);
883     strcat(pt,pb);
884     strcpy(bondk1.kb5[bondk1.nbnd5],pt);
885     bondk1.nbnd5++;
886     }
887     } else
888     {
889     message_alert("Error reading bond atom types in added constants file","Error");
890     }
891     } else if (strcmp(iptemp,"bond4") == 0)
892     {
893     iz = sscanf(line,"%s %d %d %f %f",dummy,&ia,&ib,&f1,&f2);
894     fprintf(pcmoutfile,"Bond4: AtomTypes %d %d Blen %f Bk %f \n",ia,ib,f1,f2);
895     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && iz > 3)
896     {
897     iz = check_bond4const(ia,ib,f1,f2);
898     if (iz == FALSE)
899     {
900     bondk1.s4[bondk1.nbnd4] = f1;
901     bondk1.t4[bondk1.nbnd4] = f2;
902     numeral(ia,pa,3);
903     numeral(ib,pb,3);
904     strcpy(pt,pa);
905     strcat(pt,pb);
906     strcpy(bondk1.kb4[bondk1.nbnd4],pt);
907     bondk1.nbnd4++;
908     }
909     } else
910     {
911     message_alert("Error reading bond atom types in added constants file","Error");
912     }
913     } else if (strcmp(iptemp,"bond3") == 0)
914     {
915     sscanf(line,"%s %d %d %f %f",dummy,&ia,&ib,&f1,&f2);
916     fprintf(pcmoutfile,"Bond3: AtomTypes %d %d Blen %f Bk %f \n",ia,ib,f1,f2);
917     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE)
918     {
919     iz = check_bond3const(ia,ib,f1,f2);
920     if (iz == FALSE)
921     {
922     bondk1.s3[bondk1.nbnd3] = f1;
923     bondk1.t3[bondk1.nbnd3] = f2;
924     numeral(ia,pa,3);
925     numeral(ib,pb,3);
926     strcpy(pt,pa);
927     strcat(pt,pb);
928     strcpy(bondk1.kb3[bondk1.nbnd3],pt);
929     bondk1.nbnd3++;
930     }
931     } else
932     {
933     message_alert("Error reading bond atom types in added constants file","Error");
934     }
935     } else if (strcmp(iptemp,"angle") == 0)
936     {
937     f1 =f2=f3=f4=0.0;
938     sscanf(line,"%s %d %d %d %f %f %f %f",dummy,&ia,&ib,&ic,&f1,&f2,&f3,&f4);
939     fprintf(pcmoutfile,"Angle: AtomTypes %d %d %d Anat %f Acon %f \n",ia,ib,ic,f2,f1);
940     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE)
941     {
942     iz = check_angleconst(ia,ib,ic,f1,f2,f3,f4);
943     if (iz == FALSE)
944     {
945     angk1.con[angk1.nang] = f1;
946     angk1.ang[angk1.nang][0] = f2;
947     angk1.ang[angk1.nang][1] = f3;
948     angk1.ang[angk1.nang][2] = f4;
949     numeral(ia,pa,3);
950     numeral(ib,pb,3);
951     numeral(ic,pc,3);
952     strcpy(pt,pa); strcat(pt,pb); strcat(pt,pc);
953     strcpy(angk1.ktype[angk1.nang],pt);
954     angk1.nang++;
955     }
956     } else
957     {
958     message_alert("Error reading Angle atom types in added constants file","Error");
959     }
960     } else if (strcmp(iptemp,"angle5") == 0)
961     {
962     f1 =f2=f3=f4=0.0;
963     sscanf(line,"%s %d %d %d %f %f %f %f",dummy,&ia,&ib,&ic,&f1,&f2,&f3,&f4);
964     fprintf(pcmoutfile,"Angle5: AtomTypes %d %d %d Anat %f Acon %f \n",ia,ib,ic,f2,f1);
965     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE)
966     {
967     iz = check_angle5const(ia,ib,ic,f1,f2,f3,f4);
968     if (iz == FALSE)
969     {
970     angk1.con5[angk1.nang5] = f1;
971     angk1.ang5[angk1.nang5][0] = f2;
972     angk1.ang5[angk1.nang5][1] = f3;
973     angk1.ang5[angk1.nang5][2] = f4;
974     numeral(ia,pa,3);
975     numeral(ib,pb,3);
976     numeral(ic,pc,3);
977     strcpy(pt,pa); strcat(pt,pb); strcat(pt,pc);
978     strcpy(angk1.ktype5[angk1.nang5],pt);
979     angk1.nang5++;
980     }
981     } else
982     {
983     message_alert("Error reading Angle atom types in added constants file","Error");
984     }
985     } else if (strcmp(iptemp,"angle4") == 0)
986     {
987     f1 =f2=f3=f4=0.0;
988     sscanf(line,"%s %d %d %d %f %f %f %f",dummy,&ia,&ib,&ic,&f1,&f2,&f3,&f4);
989     fprintf(pcmoutfile,"Angle4: AtomTypes %d %d %d Anat %f Acon %f \n",ia,ib,ic,f2,f1);
990     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE)
991     {
992     iz = check_angle4const(ia,ib,ic,f1,f2,f3,f4);
993     if (iz == FALSE)
994     {
995     angk1.con4[angk1.nang4] = f1;
996     angk1.ang4[angk1.nang4][0] = f2;
997     angk1.ang4[angk1.nang4][1] = f3;
998     angk1.ang4[angk1.nang4][2] = f4;
999     numeral(ia,pa,3);
1000     numeral(ib,pb,3);
1001     numeral(ic,pc,3);
1002     strcpy(pt,pa); strcat(pt,pb); strcat(pt,pc);
1003     strcpy(angk1.ktype4[angk1.nang4],pt);
1004     angk1.nang4++;
1005     }
1006     } else
1007     {
1008     message_alert("Error reading Angle atom types in added constants file","Error");
1009     }
1010     } else if (strcmp(iptemp,"angle3") == 0)
1011     {
1012     f1 =f2=f3=f4=0.0;
1013     sscanf(line,"%s %d %d %d %f %f %f %f",dummy,&ia,&ib,&ic,&f1,&f2,&f3,&f4);
1014     fprintf(pcmoutfile,"Angle3: AtomTypes %d %d %d Anat %f Acon %f \n",ia,ib,ic,f2,f1);
1015     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE)
1016     {
1017     iz = check_angle3const(ia,ib,ic,f1,f2,f3,f4);
1018     if (iz == FALSE)
1019     {
1020     angk1.con3[angk1.nang3] = f1;
1021     angk1.ang3[angk1.nang3][0] = f2;
1022     angk1.ang3[angk1.nang3][1] = f3;
1023     angk1.ang3[angk1.nang3][2] = f4;
1024     numeral(ia,pa,3);
1025     numeral(ib,pb,3);
1026     numeral(ic,pc,3);
1027     strcpy(pt,pa); strcat(pt,pb); strcat(pt,pc);
1028     strcpy(angk1.ktype3[angk1.nang3],pt);
1029     angk1.nang3++;
1030     }
1031     } else
1032     {
1033     message_alert("Error reading Angle atom types in added constants file","Error");
1034     }
1035     } else if (strcmp(iptemp,"torsion") == 0)
1036     {
1037     sscanf(line,"%s %d %d %d %d %f %d %f %d %f %d",dummy,&ia,&ib,&ic,&id,&f1,&ip1,&f2,
1038     &ip2,&f3,&ip3);
1039     fprintf(pcmoutfile,"Torsion: AtomTypes %d %d %d %d V1 %f V2 %f V3 %f\n",ia,ib,ic,id,f1,f2,f3);
1040     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE && id < MAXATOMTYPE)
1041     {
1042     iz = check_torsionconst(ia,ib,ic,id,f1,f2,f3,ip1,ip2,ip3);
1043     if (iz == FALSE)
1044     {
1045     v1[0] = f1; v1[1] = f2; v1[2] = f3;
1046     v1[3] = v1[4] = v1[5] = 0.0;
1047     se1[0] = ip1; se1[1] = ip2; se1[2] = ip3;
1048     se1[3] = se1[4] = se1[5] = 0;
1049     torphase(6, v1, se1);
1050     torkn1.tv1[torkn1.ntor] = v1[0];
1051     torkn1.phase1[torkn1.ntor] = se1[0];
1052     torkn1.tv2[torkn1.ntor] = v1[1];
1053     torkn1.phase2[torkn1.ntor] = se1[1];
1054     torkn1.tv3[torkn1.ntor] = v1[2];
1055     torkn1.phase3[torkn1.ntor] = se1[2];
1056     torkn1.tv4[torkn1.ntor] = v1[3];
1057     torkn1.phase4[torkn1.ntor] = se1[3];
1058     torkn1.tv5[torkn1.ntor] = v1[4];
1059     torkn1.phase5[torkn1.ntor] = se1[4];
1060     torkn1.tv6[torkn1.ntor] = v1[5];
1061     torkn1.phase6[torkn1.ntor] = se1[5];
1062     numeral(ia,pa,3);
1063     numeral(ib,pb,3);
1064     numeral(ic,pc,3);
1065     numeral(id,pd,3);
1066     strcpy(pt,pa);
1067     strcat(pt,pb); strcat(pt,pc); strcat(pt,pd);
1068     strcpy(torkn1.kv[torkn1.ntor],pt);
1069     torkn1.ntor++;
1070     }
1071     } else
1072     {
1073     message_alert("Error reading Torsion atom types in added constants file","Error");
1074     }
1075     } else if (strcmp(iptemp,"torsion4") == 0)
1076     {
1077     sscanf(line,"%s %d %d %d %d %f %d %f %d %f %d",dummy,&ia,&ib,&ic,&id,&f1,&ip1,&f2,
1078     &ip2,&f3,&ip3);
1079     fprintf(pcmoutfile,"Torsion4: AtomTypes %d %d %d %d V1 %f V2 %f V3 %f\n",ia,ib,ic,id,f1,f2,f3);
1080     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE && id < MAXATOMTYPE)
1081     {
1082     iz = check_torsion4const(ia,ib,ic,id,f1,f2,f3,ip1,ip2,ip3);
1083     if (iz == FALSE)
1084     {
1085     v1[0] = f1; v1[1] = f2; v1[2] = f3;
1086     v1[3] = v1[4] = v1[5] = 0.0;
1087     se1[0] = ip1; se1[1] = ip2; se1[2] = ip3;
1088     se1[3] = se1[4] = se1[5] = 0;
1089     torphase(6, v1, se1);
1090     torkn1.tv41[torkn1.ntor4] = v1[0];
1091     torkn1.phase41[torkn1.ntor4] = se1[0];
1092     torkn1.tv42[torkn1.ntor4] = v1[1];
1093     torkn1.phase42[torkn1.ntor4] = se1[1];
1094     torkn1.tv43[torkn1.ntor4] = v1[2];
1095     torkn1.phase43[torkn1.ntor4] = se1[2];
1096     numeral(ia,pa,3);
1097     numeral(ib,pb,3);
1098     numeral(ic,pc,3);
1099     numeral(id,pd,3);
1100     strcpy(pt,pa);
1101     strcat(pt,pb); strcat(pt,pc); strcat(pt,pd);
1102     strcpy(torkn1.kv4[torkn1.ntor4],pt);
1103     torkn1.ntor4++;
1104     }
1105     } else
1106     {
1107     message_alert("Error reading Torsion atom types in added constants file","Error");
1108     }
1109     } else if (strcmp(iptemp,"torsion5") == 0)
1110     {
1111     sscanf(line,"%s %d %d %d %d %f %d %f %d %f %d",dummy,&ia,&ib,&ic,&id,&f1,&ip1,&f2,
1112     &ip2,&f3,&ip3);
1113     fprintf(pcmoutfile,"Torsion5: AtomTypes %d %d %d %d V1 %f V2 %f V3 %f\n",ia,ib,ic,id,f1,f2,f3);
1114     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE && id < MAXATOMTYPE)
1115     {
1116     iz = check_torsion5const(ia,ib,ic,id,f1,f2,f3,ip1,ip2,ip3);
1117     if (iz == FALSE)
1118     {
1119     v1[0] = f1; v1[1] = f2; v1[2] = f3;
1120     v1[3] = v1[4] = v1[5] = 0.0;
1121     se1[0] = ip1; se1[1] = ip2; se1[2] = ip3;
1122     se1[3] = se1[4] = se1[5] = 0;
1123     torphase(6, v1, se1);
1124     torkn1.tv51[torkn1.ntor5] = v1[0];
1125     torkn1.phase51[torkn1.ntor5] = se1[0];
1126     torkn1.tv52[torkn1.ntor5] = v1[1];
1127     torkn1.phase52[torkn1.ntor5] = se1[1];
1128     torkn1.tv53[torkn1.ntor5] = v1[2];
1129     torkn1.phase53[torkn1.ntor5] = se1[2];
1130     numeral(ia,pa,3);
1131     numeral(ib,pb,3);
1132     numeral(ic,pc,3);
1133     numeral(id,pd,3);
1134     strcpy(pt,pa);
1135     strcat(pt,pb); strcat(pt,pc); strcat(pt,pd);
1136     strcpy(torkn1.kv5[torkn1.ntor5],pt);
1137     torkn1.ntor5++;
1138     }
1139     } else
1140     {
1141     message_alert("Error reading Torsion atom types in added constants file","Error");
1142     }
1143     } else if (strcmp(iptemp,"vdw") == 0)
1144     {
1145     sscanf(line,"%s %d %f %f %d %d %d",dummy,&ia,&f1,&f2,&if3,&if4,&if5);
1146     fprintf(pcmoutfile,"VDW: AtomTypes %d Rad %f Eps %f \n",ia,f1,f2);
1147     vdw1.rad[ia] = f1;
1148     vdw1.eps[ia] = f2;
1149     vdw1.lpd[ia] = if3;
1150     vdw1.ihdon[ia] = if4;
1151     vdw1.ihtyp[ia] = if5;
1152     } else if (strcmp(iptemp,"metal") == 0)
1153     {
1154     sscanf(line,"%s %d %f %f",dummy,&ia,&f1,&f2);
1155     fprintf(pcmoutfile,"Metal: AtomTypes %d Rad %f Eps %f \n",ia,f1,f2);
1156     metaldata.radius[ia] = f1;
1157     metaldata.eps[ia] = f2;
1158     } else if (strcmp(iptemp,"vdwpr") == 0)
1159     {
1160     sscanf(line,"%s %d %d %f %f",dummy, &ia, &ib, &f1,&f2);
1161     fprintf(pcmoutfile,"VDWPR: AtomTypes %d %d Rad %f Eps %f \n",ia,ib,f1,f2);
1162     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE)
1163     {
1164     iz = check_vdwpr(ia,ib,f1,f2);
1165     if (iz == FALSE)
1166     {
1167     numeral(ia,pa,3);
1168     numeral(ib,pb,3);
1169     strcpy(pt,pa); strcat(pt,pb);
1170     strcpy(vdwpr_k.kv[vdwpr_k .nvdwpr],pt);
1171     vdwpr_k.ia1[vdwpr_k .nvdwpr] = ia;
1172     vdwpr_k.ia2[vdwpr_k .nvdwpr] = ib;
1173     vdwpr_k.radius[vdwpr_k .nvdwpr] = f1;
1174     vdwpr_k.eps[vdwpr_k .nvdwpr] = f2;
1175     vdwpr_k .nvdwpr++;
1176     }
1177     } else
1178     {
1179     message_alert("Error reading VDWPR atom types in added constants file","Error");
1180     }
1181     } else if (strcmp(iptemp,"opbend") == 0)
1182     {
1183     sscanf( line, "%s %d %d %d %d %f", dummy, &ia, &ib, &ic,&id,&f1 );
1184     fprintf(pcmoutfile,"OPBEND: AtomTypes %d %d %d %d Const %f \n",ia,ib,ic,id,f1);
1185     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE)
1186     {
1187     iz = check_opbend(ia,ib,ic,id,f1);
1188     if (iz == FALSE)
1189     {
1190     numeral(ia,pa,3);
1191     numeral(ib,pb,3);
1192     numeral(ic,pc,3);
1193     numeral(id,pd,3);
1194     strcpy(pt,pa); strcat(pt,pb); strcat(pt,pc); strcat(pt,pd);
1195     strcpy(ooplane_k.iopb[ooplane_k.nopbend],pt);
1196     ooplane_k.nopbend++;
1197     }
1198     } else
1199     {
1200     message_alert("Error reading Opbend atom types in added constants file","Error");
1201     }
1202     } else if (strcmp(iptemp,"strbnd") == 0)
1203     {
1204     id = 0;
1205     sscanf(line,"%s %d %d %d %f %f %f",dummy, &ia, &ib, &ic,
1206     &f1,&f2,&f3);
1207     fprintf(pcmoutfile,"STRBEND: AtomTypes %d %d Const %f %f %f\n",ia,ib,f1,f2,f3);
1208     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE && ic < MAXATOMTYPE)
1209     {
1210     iz = check_strbnd(ia,ib,ic,f1,f2,f3);
1211     if (iz == FALSE)
1212     {
1213     numeral(ia,pa,3);
1214     numeral(ib,pb,3);
1215     numeral(ic,pc,3);
1216     strcpy(pt,pa); strcat(pt,pb); strcat(pt,pc);
1217     strcpy(crossterm_k.stbn[crossterm_k.nstrbnd],pt);
1218     crossterm_k.stbncon[crossterm_k.nstrbnd][0] = f1;
1219     crossterm_k.stbncon[crossterm_k.nstrbnd][1] = f2;
1220     crossterm_k.stbncon[crossterm_k.nstrbnd][2] = f3;
1221     crossterm_k.stbnindex[crossterm_k.nstrbnd] = 0;
1222     crossterm_k.nstrbnd++;
1223     }
1224     } else
1225     {
1226     message_alert("Error reading Strbend atom types in added constants file","Error");
1227     }
1228     } else if (strcmp(iptemp,"dipole") == 0)
1229     {
1230     sscanf(line,"%s %d %d %f",dummy,&ia,&ib,&f1);
1231     fprintf(pcmoutfile,"DIPOLE: AtomTypes %d %d Bond Moment %f \n",ia,ib,f1);
1232     if (ia < MAXATOMTYPE && ib < MAXATOMTYPE)
1233     {
1234     iz = check_dipoleconst(ia,ib,f1);
1235     if (iz == FALSE)
1236     {
1237     dipole_k.bmom[dipole_k.ndipole] = f1;
1238     numeral(ia,pa,3);
1239     numeral(ib,pb,3);
1240     strcpy(pt,pa); strcat(pt,pb);
1241     strcpy(dipole_k.kb[dipole_k.ndipole],pt);
1242     dipole_k.ndipole++;
1243     }
1244     } else
1245     {
1246     message_alert("Error reading Dipole atom types in added constants file","Error");
1247     }
1248     }
1249     }
1250     fclose(wfile);
1251     }
1252     /* ========================================== */
1253     int check_vdwpr(int ia, int ib,float f1, float f2)
1254     {
1255     int i;
1256     char pa[3],pb[3],pt[13];
1257    
1258     numeral(ia,pa,3);
1259     numeral(ib,pb,3);
1260     strcpy(pt,pa);
1261     strcat(pt,pb);
1262     for (i=0; i < vdwpr_k .nvdwpr; i++)
1263     {
1264     if (strcmp(pt,vdwpr_k.kv[i]) == 0)
1265     {
1266     vdwpr_k.radius[i] = f1;
1267     vdwpr_k.eps[i] = f2;
1268     return TRUE;
1269     }
1270     }
1271     return FALSE;
1272     }
1273     /* ========================================== */
1274     int check_strbnd(int ia, int ib, int ic,float f1, float f2, float f3)
1275     {
1276     int i;
1277     char pa[3],pb[3],pc[3],pt[13];
1278    
1279     numeral(ia,pa,3);
1280     numeral(ib,pb,3);
1281     numeral(ic,pc,3);
1282     strcpy(pt,pa);
1283     strcat(pt,pb);
1284     strcat(pt,pc);
1285     for (i=0; i < crossterm_k.nstrbnd; i++)
1286     {
1287     if (strcmp(pt,crossterm_k.stbn[i]) == 0)
1288     {
1289     crossterm_k.stbncon[i][0] = f1;
1290     crossterm_k.stbncon[i][1] = f2;
1291     crossterm_k.stbncon[i][2] = f3;
1292     return TRUE;
1293     }
1294     }
1295     return FALSE;
1296     }
1297     /* ========================================== */
1298     int check_opbend(int ia,int ib,int ic,int id, float ftemp)
1299     {
1300     int i;
1301     char pa[3],pb[3],pc[3],pd[3],pt[13];
1302     numeral(ia,pa,3);
1303     numeral(ib,pb,3);
1304     numeral(ic,pc,3);
1305     numeral(id,pd,3);
1306     strcpy(pt,pa);
1307     strcat(pt,pb);
1308     strcat(pt,pc);
1309     strcat(pt,pd);
1310     for (i=0; i < ooplane_k.nopbend; i++)
1311     {
1312     if (strcmp(pt,ooplane_k.iopb[i]) == 0)
1313     {
1314     ooplane_k.copb[i] = ftemp;
1315     return TRUE;
1316     }
1317     }
1318     return FALSE;
1319     }
1320     /* ========================================== */
1321     int check_dipoleconst(int ia,int ib,float bmom)
1322     {
1323     int i;
1324     char pa[3],pb[3],pt[13];
1325    
1326     numeral(ia,pa,3);
1327     numeral(ib,pb,3);
1328     strcpy(pt,pa); strcat(pt,pb);
1329     for (i=0; i < dipole_k.ndipole; i++)
1330     {
1331     if (strcmp(pt,dipole_k.kb[i]) == 0)
1332     {
1333     dipole_k.bmom[i] = bmom;
1334     return TRUE;
1335     }
1336     }
1337     return FALSE;
1338     }
1339     /* ========================================== */
1340     int check_torsionconst(int ia,int ib,int ic,int id,float v1,float v2, float v3,
1341     int ip1, int ip2, int ip3)
1342     {
1343     int i;
1344     char pa[3],pb[3],pc[3],pd[3],pt[13];
1345    
1346     numeral(ia,pa,3);
1347     numeral(ib,pb,3);
1348     numeral(ic,pc,3);
1349     numeral(id,pd,3);
1350     strcpy(pt,pa);
1351     strcat(pt,pb);
1352     strcat(pt,pc);
1353     strcat(pt,pd);
1354     for (i= 0; i < torkn1.ntor; i++)
1355     {
1356     if (strcmp(pt,torkn1.kv[i]) == 0)
1357     {
1358     torkn1.tv1[i] = v1;
1359     torkn1.tv2[i] = v2;
1360     torkn1.tv3[i] = v3;
1361     return TRUE;
1362     }
1363     }
1364     return FALSE;
1365     }
1366     /* ========================================== */
1367     int check_torsion4const(int ia,int ib,int ic,int id,float v1,float v2, float v3,
1368     int ip1, int ip2, int ip3)
1369     {
1370     int i;
1371     char pa[3],pb[3],pc[3],pd[3],pt[13];
1372    
1373     numeral(ia,pa,3);
1374     numeral(ib,pb,3);
1375     numeral(ic,pc,3);
1376     numeral(id,pd,3);
1377     strcpy(pt,pa);
1378     strcat(pt,pb);
1379     strcat(pt,pc);
1380     strcat(pt,pd);
1381     for (i= 0; i < torkn1.ntor4; i++)
1382     {
1383     if (strcmp(pt,torkn1.kv4[i]) == 0)
1384     {
1385     torkn1.tv41[i] = v1;
1386     torkn1.tv42[i] = v2;
1387     torkn1.tv43[i] = v3;
1388     return TRUE;
1389     }
1390     }
1391     return FALSE;
1392     }
1393     /* ========================================== */
1394     int check_torsion5const(int ia,int ib,int ic,int id,float v1,float v2, float v3,
1395     int ip1, int ip2, int ip3)
1396     {
1397     int i;
1398     char pa[3],pb[3],pc[3],pd[3],pt[13];
1399    
1400     numeral(ia,pa,3);
1401     numeral(ib,pb,3);
1402     numeral(ic,pc,3);
1403     numeral(id,pd,3);
1404     strcpy(pt,pa);
1405     strcat(pt,pb);
1406     strcat(pt,pc);
1407     strcat(pt,pd);
1408     for (i= 0; i < torkn1.ntor5; i++)
1409     {
1410     if (strcmp(pt,torkn1.kv5[i]) == 0)
1411     {
1412     torkn1.tv51[i] = v1;
1413     torkn1.tv52[i] = v2;
1414     torkn1.tv53[i] = v3;
1415     return TRUE;
1416     }
1417     }
1418     return FALSE;
1419     }
1420     /* ========================================== */
1421     int check_angleconst(int ia,int ib,int ic,float f1, float f2,float f3,float f4)
1422     {
1423     int i;
1424     char pa[3],pb[3],pc[3],pt[13];
1425    
1426     numeral(ia,pa,3);
1427     numeral(ib,pb,3);
1428     numeral(ic,pc,3);
1429     strcpy(pt,pa);
1430     strcat(pt,pb);
1431     strcat(pt,pc);
1432     for (i=0; i < angk1.nang; i++)
1433     {
1434     if (strcmp(pt,angk1.ktype[i]) == 0)
1435     {
1436     angk1.con[i] = f1;
1437     angk1.ang[i][0] = f2;
1438     angk1.ang[i][1] = f3;
1439     angk1.ang[i][2] = f4;
1440     return TRUE;
1441     }
1442     }
1443     return FALSE;
1444     }
1445     /* ========================================== */
1446     int check_angle5const(int ia,int ib,int ic,float f1, float f2,float f3,float f4)
1447     {
1448     int i;
1449     char pa[3],pb[3],pc[3],pt[13];
1450    
1451     numeral(ia,pa,3);
1452     numeral(ib,pb,3);
1453     numeral(ic,pc,3);
1454     strcpy(pt,pa);
1455     strcat(pt,pb);
1456     strcat(pt,pc);
1457     for (i=0; i < angk1.nang5; i++)
1458     {
1459     if (strcmp(pt,angk1.ktype5[i]) == 0)
1460     {
1461     angk1.con5[i] = f1;
1462     angk1.ang5[i][0] = f2;
1463     angk1.ang5[i][1] = f3;
1464     angk1.ang5[i][2] = f4;
1465     return TRUE;
1466     }
1467     }
1468     return FALSE;
1469     }
1470     /* ========================================== */
1471     int check_angle4const(int ia,int ib,int ic,float f1, float f2,float f3,float f4)
1472     {
1473     int i;
1474     char pa[3],pb[3],pc[3],pt[13];
1475    
1476     numeral(ia,pa,3);
1477     numeral(ib,pb,3);
1478     numeral(ic,pc,3);
1479     strcpy(pt,pa);
1480     strcat(pt,pb);
1481     strcat(pt,pc);
1482     for (i=0; i < angk1.nang4; i++)
1483     {
1484     if (strcmp(pt,angk1.ktype4[i]) == 0)
1485     {
1486     angk1.con4[i] = f1;
1487     angk1.ang4[i][0] = f2;
1488     angk1.ang4[i][1] = f3;
1489     angk1.ang4[i][2] = f4;
1490     return TRUE;
1491     }
1492     }
1493     return FALSE;
1494     }
1495     /* ========================================== */
1496     int check_angle3const(int ia,int ib,int ic,float f1, float f2,float f3,float f4)
1497     {
1498     int i;
1499     char pa[3],pb[3],pc[3],pt[13];
1500    
1501     numeral(ia,pa,3);
1502     numeral(ib,pb,3);
1503     numeral(ic,pc,3);
1504     strcpy(pt,pa);
1505     strcat(pt,pb);
1506     strcat(pt,pc);
1507     for (i=0; i < angk1.nang3; i++)
1508     {
1509     if (strcmp(pt,angk1.ktype3[i]) == 0)
1510     {
1511     angk1.con3[i] = f1;
1512     angk1.ang3[i][0] = f2;
1513     angk1.ang3[i][1] = f3;
1514     angk1.ang3[i][2] = f4;
1515     return TRUE;
1516     }
1517     }
1518     return FALSE;
1519     }
1520     /* ========================================== */
1521     int check_bondconst(int ia,int ib,float f1, float f2)
1522     {
1523     int i;
1524     char pa[3],pb[3],pt[13];
1525    
1526     numeral(ia,pa,3);
1527     numeral(ib,pb,3);
1528     strcpy(pt,pa);
1529     strcat(pt,pb);
1530     for (i=0; i < bondk1.nbnd; i++)
1531     {
1532     if (strcmp(pt,bondk1.kb[i]) == 0)
1533     {
1534     bondk1.s[i] = f1;
1535     bondk1.t[i] = f2;
1536     return TRUE;
1537     }
1538     }
1539     return FALSE;
1540     }
1541     /* ========================================== */
1542     int check_bond5const(int ia,int ib,float f1, float f2)
1543     {
1544     int i;
1545     char pa[3],pb[3],pt[13];
1546    
1547     numeral(ia,pa,3);
1548     numeral(ib,pb,3);
1549     strcpy(pt,pa);
1550     strcat(pt,pb);
1551     for (i=0; i < bondk1.nbnd5; i++)
1552     {
1553     if (strcmp(pt,bondk1.kb5[i]) == 0)
1554     {
1555     bondk1.s5[i] = f1;
1556     bondk1.t5[i] = f2;
1557     return TRUE;
1558     }
1559     }
1560     return FALSE;
1561     }
1562     /* ========================================== */
1563     int check_bond4const(int ia,int ib,float f1, float f2)
1564     {
1565     int i;
1566     char pa[3],pb[3],pt[13];
1567    
1568     numeral(ia,pa,3);
1569     numeral(ib,pb,3);
1570     strcpy(pt,pa);
1571     strcat(pt,pb);
1572     for (i=0; i < bondk1.nbnd4; i++)
1573     {
1574     if (strcmp(pt,bondk1.kb4[i]) == 0)
1575     {
1576     bondk1.s4[i] = f1;
1577     bondk1.t4[i] = f2;
1578     return TRUE;
1579     }
1580     }
1581     return FALSE;
1582     }
1583     /* ========================================== */
1584     int check_bond3const(int ia,int ib,float f1, float f2)
1585     {
1586     int i;
1587     char pa[3],pb[3],pt[13];
1588    
1589     numeral(ia,pa,3);
1590     numeral(ib,pb,3);
1591     strcpy(pt,pa);
1592     strcat(pt,pb);
1593     for (i=0; i < bondk1.nbnd3; i++)
1594     {
1595     if (strcmp(pt,bondk1.kb3[i]) == 0)
1596     {
1597     bondk1.s3[i] = f1;
1598     bondk1.t3[i] = f2;
1599     return TRUE;
1600     }
1601     }
1602     return FALSE;
1603     }
1604     /* ============================================= */
1605 tjod 25 void read_datafiles(char *pfile)
1606 tjod 3 {
1607 tjod 25 /*
1608 tjod 3 char filename[256];
1609    
1610     char *ev = NULL;
1611    
1612     if ( (ev=getenv("ENGINE_DIR")) != NULL)
1613     {
1614     strcpy(filename,ev);
1615     strcat(filename,"/");
1616     strcat(filename,datafile);
1617     } else
1618     strcpy(filename,datafile);
1619     read_parameterfile(filename);
1620 tjod 25 */
1621     #include "mmxconst.h"
1622     #include "mmff94.h"
1623     if (!strncmp(pfile, "mmxconst",8)) {
1624     read_parameterfile(mmxconst);
1625     } else if (!strncmp(pfile, "mmff94",6)) {
1626     read_parameterfile(mmff94);
1627     }
1628 tjod 3 return;
1629     }
1630     // =================================================
1631 tjod 25 void read_parameterfile(char **parray)
1632 tjod 3 {
1633 tjod 25 char dumm[21], iptemp[21];
1634 tjod 3 char pa[4],pb[4],pt[7];
1635     char pc[4],pd[4],pang[10],ptor[13];
1636     int iz;
1637     int ihtype,lpde,ihdonor;
1638     float radius,epsilon,alpha,ntmp,atmp,gtmp;
1639     float v1[6];
1640     int se1[6];
1641     int iser, inumber, iligand, iclass, iclass1, iclass2, ivalence;
1642     char symbol[3], descript[25];
1643     float wtt;
1644     float charge, fchrg,tchrg;
1645     int i, iatm1, iatm2, ii, ja,
1646     jb, kk, kp1, kp2, kp3, kp4, kt1, kt2, kt3;
1647    
1648 tjod 25 int irec=0;
1649     char *line;
1650 tjod 3
1651 tjod 25 /*
1652     FILE *datafile;
1653 tjod 3 datafile = fopen(string, "rt");
1654    
1655     if (datafile == NULL)
1656     {
1657     char message[80];
1658     sprintf(message,"Unable to open the data file %s.\nPlease check that this file was installed\n and is in the same directory as PCWIN",string);
1659     message_alert(message,"Read Parameter");
1660 tjod 15 fprintf(stderr,"Unable to open the data file %s.\nPlease check that this file was installed\n and is in the same directory as PCWIN\n",string);
1661 tjod 3 fclose(datafile);
1662     return;
1663     }
1664    
1665     while ( FetchRecord(datafile,line))
1666 tjod 25 */
1667     while (strlen(parray[irec]))
1668 tjod 3 {
1669 tjod 25 line = parray[irec];
1670     ++irec;
1671 tjod 3 sscanf(line,"%s",iptemp);
1672     /* force field descriptors */
1673     if (strcmp(iptemp,"forcefield") == 0)
1674     {
1675     sscanf(line,"%s %s",dumm, field.name);
1676     if (strcmp(field.name,"AMBER95") == 0)
1677     field.type = AMBER;
1678     else if (strcmp(field.name,"MMX") == 0)
1679     field.type = MMX;
1680     else if (strcmp(field.name,"MM2-1991") == 0)
1681     field.type = MM2;
1682     else if (strcmp(field.name,"MM2-PLUS") == 0)
1683     field.type = MM2;
1684     else if (strcmp(field.name,"MM3-1992") == 0 || strcmp(field.name,"MM3-1996") == 0 )
1685     field.type = MM3;
1686     else if (strcmp(field.name,"CHARMM22-PROTEIN") == 0)
1687     field.type = CHARMM;
1688     else if (strcmp(field.name,"MMFF94") == 0)
1689     field.type = MMFF94;
1690     else if (strcmp(field.name,"OPLSAA") == 0)
1691     field.type = OPLSAA;
1692     else
1693     field.type = UNKNOWN;
1694     } else if (strcmp(iptemp,"bondunit") == 0)
1695     {
1696     sscanf(line,"%s %f",dumm, &field.bondunit);
1697     units.bndunit = field.bondunit;
1698     } else if (strcmp(iptemp,"bond-cubic") == 0)
1699     {
1700     sscanf(line,"%s %f",dumm, &field.bond_cubic);
1701     units.cbnd = field.bond_cubic;
1702     } else if (strcmp(iptemp,"bond-quartic") == 0)
1703     {
1704     sscanf(line,"%s %f",dumm, &field.bond_quartic);
1705     units.qbnd = field.bond_quartic;
1706     } else if (strcmp(iptemp,"angleunit") == 0)
1707     {
1708     sscanf(line,"%s %f",dumm, &field.angleunit);
1709     units.angunit = field.angleunit;
1710     } else if (strcmp(iptemp,"angle-cubic") == 0)
1711     {
1712     sscanf(line,"%s %f",dumm, &field.angle_cubic);
1713     units.cang = field.angle_cubic;
1714     } else if (strcmp(iptemp,"angle-quartic") == 0)
1715     {
1716     sscanf(line,"%s %f",dumm, &field.angle_quartic);
1717     units.qang = field.angle_quartic;
1718     } else if (strcmp(iptemp,"angle-pentic") == 0)
1719     {
1720     sscanf(line,"%s %f",dumm, &field.angle_pentic);
1721     units.pang = field.angle_pentic;
1722     } else if (strcmp(iptemp,"angle-sextic") == 0)
1723     {
1724     sscanf(line,"%s %f",dumm, &field.angle_sextic);
1725     units.sang = field.angle_sextic;
1726     } else if (strcmp(iptemp,"str-bndunit") == 0)
1727     {
1728     sscanf(line,"%s %f",dumm, &field.str_bndunit);
1729     units.stbnunit = field.str_bndunit;
1730     } else if (strcmp(iptemp,"ang-angunit") == 0)
1731     {
1732     sscanf(line,"%s %f",dumm, &field.ang_angunit);
1733     units.aaunit = field.ang_angunit;
1734     } else if (strcmp(iptemp,"torsionunit") == 0)
1735     {
1736     sscanf(line,"%s %f",dumm, &field.torsionunit);
1737     units.torsunit = field.torsionunit;
1738     } else if (strcmp(iptemp,"str-torunit") == 0)
1739     {
1740     sscanf(line,"%s %f",dumm, &field.str_torunit);
1741     units.storunit = field.str_torunit;
1742     } else if (strcmp(iptemp,"vdwtype") == 0)
1743     {
1744     sscanf(line,"%s %s",dumm, field.vdwtype);
1745     if (strcmp(field.vdwtype ,"LENNARD-JONES") == 0)
1746     pot.use_lj = TRUE;
1747     else if (strcmp(field.vdwtype ,"BUCKINGHAM") == 0)
1748     pot.use_buck = TRUE;
1749     else if (strcmp(field.vdwtype ,"BUFFERED-14-7") == 0)
1750     pot.use_hal = TRUE;
1751     else if (strcmp(field.vdwtype ,"BGAUSSIAN") == 0)
1752     pot.use_gauss = TRUE;
1753     } else if (strcmp(iptemp,"radiustype") == 0)
1754     {
1755     sscanf(line,"%s %s",dumm, field.radiustype);
1756     } else if (strcmp(iptemp,"radiussize") == 0)
1757     {
1758     sscanf(line,"%s %s",dumm, field.radiussize);
1759     } else if (strcmp(iptemp,"radiusrule") == 0)
1760     {
1761     sscanf(line,"%s %s",dumm, field.radiusrule);
1762     } else if (strcmp(iptemp,"epsilonrule") == 0)
1763     {
1764     sscanf(line,"%s %s",dumm, field.epsrule);
1765     } else if (strcmp(iptemp,"a-expterm") == 0)
1766     {
1767     sscanf(line,"%s %f",dumm, &field.a_expterm);
1768     units.aterm = field.a_expterm;
1769     } else if (strcmp(iptemp,"b-expterm") == 0)
1770     {
1771     sscanf(line,"%s %f",dumm, &field.b_expterm);
1772     units.bterm = field.b_expterm;
1773     } else if (strcmp(iptemp,"c-expterm") == 0)
1774     {
1775     sscanf(line,"%s %f",dumm, &field.c_expterm);
1776     units.cterm = field.c_expterm;
1777     } else if (strcmp(iptemp,"vdw-14-scale") == 0)
1778     {
1779     sscanf(line,"%s %f",dumm, &field.vdw_14scale);
1780     units.v14scale = field.vdw_14scale;
1781     } else if (strcmp(iptemp,"dielectric") == 0)
1782     {
1783     sscanf(line,"%s %f",dumm, &field.dielectric);
1784     if (user.dielec == FALSE)
1785     units.dielec = field.dielectric;
1786     } else if (strcmp(iptemp,"chg-14-scale") == 0)
1787     {
1788     sscanf(line,"%s %f",dumm, &field.chg_14scale);
1789     units.chgscale = field.chg_14scale;
1790     } else if (strcmp(iptemp,"atom") == 0) /* atom type descriptors */
1791     {
1792     if (atom_k.natomtype >= MAXATOMTYPE )
1793     {
1794     message_alert("Error - Too many atom types","Read Parameter");
1795 tjod 15 fprintf(stderr,"Maximum number of atom types exceeded. Others will be ignored\n");
1796 tjod 25 //fclose(datafile);
1797 tjod 3 return;
1798     }
1799    
1800     sscanf(line, "%s %d %s %22c %d %f %d %d %d %d %d",dumm,&iser, symbol, descript, &inumber, &wtt,
1801     &iligand, &ivalence, &iclass, &iclass1, &iclass2);
1802    
1803     atom_k.type[iser] = iser;
1804     strcpy(atom_k.symbol[iser], symbol);
1805     strcpy(atom_k.description[iser], descript);
1806     atom_k.description[iser][20] = '\0';
1807     atom_k.number[iser] = inumber;
1808     atom_k.weight[iser] = wtt;
1809     atom_k.ligands[iser] = iligand;
1810     atom_k.tclass[iser] = iclass;
1811     atom_k.valency[iser] = ivalence;
1812     atom_k.tclass1[iser] = iclass1;
1813     atom_k.tclass2[iser] = iclass2;
1814     atom_k.natomtype++;
1815     } else if ( strcmp(iptemp,"metal") == 0) /* metal data */
1816     {
1817     if (metaldata.nmetal > 200)
1818     {
1819     message_alert("Error-Too many metals in Parameter file","Read Parameter");
1820 tjod 25 //fclose(datafile);
1821 tjod 3 return;
1822     }
1823     sscanf(line,"%s %d %s %f %f",dumm,&metaldata.type[metaldata.nmetal], metaldata.name[metaldata.nmetal],&metaldata.radius[metaldata.nmetal],
1824     &metaldata.eps[metaldata.nmetal]);
1825     metaldata.nmetal++;
1826     } else if( strcmp(iptemp,"bond") == 0 ) /* bond constants */
1827     {
1828     if (bondk1.nbnd >= MAXBONDCONST)
1829     {
1830     message_alert("Error - Too many bond constants","Read Parameter");
1831 tjod 15 fprintf(stderr,"Maximum number of bond constants exceeded. Others will be ignored\n");
1832 tjod 25 //fclose(datafile);
1833 tjod 3 return;
1834     }
1835    
1836     sscanf( line, "%s %d %d %f %f", dumm, &ii, &kk,
1837     &bondk1.s[bondk1.nbnd], &bondk1.t[bondk1.nbnd]);
1838     numeral(ii,pa,3);
1839     numeral(kk,pb,3);
1840     strcpy(pt,pa);
1841     strcat(pt,pb);
1842     strcpy(bondk1.kb[bondk1.nbnd],pt);
1843     bondk1.nbnd++;
1844     } else if( strcmp(iptemp,"bond3") == 0 ) /* bond constants */
1845     {
1846     if (bondk1.nbnd3 >= MAXBOND3CONST)
1847     {
1848     message_alert("Error - Too many bond3 constants","Read Parameter");
1849 tjod 15 fprintf(stderr,"Maximum number of bond3 constants exceeded. Others will be ignored\n");
1850 tjod 25 //fclose(datafile);
1851 tjod 3 return;
1852     }
1853    
1854     sscanf( line, "%s %d %d %f %f", dumm, &ii, &kk,
1855     &bondk1.s3[bondk1.nbnd3], &bondk1.t3[bondk1.nbnd3]);
1856     numeral(ii,pa,3);
1857     numeral(kk,pb,3);
1858     strcpy(pt,pa);
1859     strcat(pt,pb);
1860     strcpy(bondk1.kb3[bondk1.nbnd3],pt);
1861     // bondk1.kb3[bondk1.nbnd3] = ii*100 + kk;
1862     bondk1.nbnd3++;
1863     } else if( strcmp(iptemp,"bond4") == 0 ) /* bond constants */
1864     {
1865     if (bondk1.nbnd4 >= MAXBOND4CONST)
1866     {
1867     message_alert("Error - Too many bond4 constants","Read Parameter");
1868 tjod 15 fprintf(stderr,"Maximum number of bond4 constants exceeded. Others will be ignored\n");
1869 tjod 25 //fclose(datafile);
1870 tjod 3 return;
1871     }
1872    
1873     sscanf( line, "%s %d %d %f %f", dumm, &ii, &kk,
1874     &bondk1.s4[bondk1.nbnd4], &bondk1.t4[bondk1.nbnd4]);
1875     numeral(ii,pa,3);
1876     numeral(kk,pb,3);
1877     strcpy(pt,pa);
1878     strcat(pt,pb);
1879     strcpy(bondk1.kb4[bondk1.nbnd4],pt);
1880     // bondk1.kb4[bondk1.nbnd4] = ii*100 + kk;
1881     bondk1.nbnd4++;
1882     } else if( strcmp(iptemp,"bond5") == 0 ) /* bond constants */
1883     {
1884     if (bondk1.nbnd5>= MAXBOND5CONST)
1885     {
1886     message_alert("Error - Too many bond5 constants","Read Parameter");
1887 tjod 15 fprintf(stderr,"Maximum number of bond5 constants exceeded. Others will be ignored\n");
1888 tjod 25 //fclose(datafile);
1889 tjod 3 return;
1890     }
1891    
1892     sscanf( line, "%s %d %d %f %f", dumm, &ii, &kk,
1893     &bondk1.s5[bondk1.nbnd5], &bondk1.t5[bondk1.nbnd5]);
1894     numeral(ii,pa,3);
1895     numeral(kk,pb,3);
1896     strcpy(pt,pa);
1897     strcat(pt,pb);
1898     strcpy(bondk1.kb5[bondk1.nbnd5],pt);
1899     bondk1.nbnd5++;
1900     } else if( strcmp(iptemp,"bonddel") == 0 ) /* bond constants */
1901     {
1902     if (bondk1.ndeloc>= MAXBONDDELOC)
1903     {
1904     message_alert("Error - Too many delocalized bond constants","Read Parameter");
1905 tjod 15 fprintf(stderr,"Maximum number of delocalized bond constants exceeded. Others will be ignored\n");
1906 tjod 25 //fclose(datafile);
1907 tjod 3 return;
1908     }
1909    
1910     sscanf( line, "%s %d %d %f %f", dumm, &ii, &kk,
1911     &bondk1.sdel[bondk1.ndeloc], &bondk1.tdel[bondk1.ndeloc]);
1912     numeral(ii,pa,3);
1913     numeral(kk,pb,3);
1914     strcpy(pt,pa);
1915     strcat(pt,pb);
1916     strcpy(bondk1.kbdel[bondk1.ndeloc],pt);
1917     bondk1.ndeloc++;
1918     } else if (strcmp(iptemp, "electroi") == 0)
1919     {
1920     sscanf(line,"%s %d %d %d %f",dumm,&electroneg.itype[electroneg.nelecti],
1921     &electroneg.ibond[electroneg.nelecti],
1922     &electroneg.iattach[electroneg.nelecti],
1923     &electroneg.icorr[electroneg.nelecti]);
1924     electroneg.nelecti++;
1925     } else if (strcmp(iptemp, "electroj") == 0)
1926     {
1927     sscanf(line,"%s %d %d %d %f",dumm,&electroneg.jtype[electroneg.nelectj],
1928     &electroneg.jbond[electroneg.nelectj],
1929     &electroneg.jattach[electroneg.nelectj],
1930     &electroneg.jcorr[electroneg.nelectj]);
1931     electroneg.nelectj++;
1932     } else if (strcmp(iptemp, "piatom") == 0)
1933     {
1934     if (piatomk.npiatom > MAXPIATOM)
1935     {
1936     message_alert("Error - Too many piatom constants","Read Parameter");
1937 tjod 15 fprintf(stderr,"Maximum number of piatom constants exceeded. Others will be ignored\n");
1938 tjod 25 //fclose(datafile);
1939 tjod 3 return;
1940     }
1941     sscanf(line,"%s %d %f %f %f %f %f %f %f %f %d", dumm, &piatomk.kat[piatomk.npiatom],
1942     &piatomk.q[piatomk.npiatom], &piatomk.ion[piatomk.npiatom], &piatomk.emz[piatomk.npiatom],
1943     &piatomk.zor[piatomk.npiatom],&piatomk.zz[piatomk.npiatom], &piatomk.w1[piatomk.npiatom],
1944     &piatomk.w2[piatomk.npiatom], &piatomk.w3[piatomk.npiatom],&piatomk.qp[piatomk.npiatom]);
1945     piatomk.npiatom++;
1946     } else if (strcmp(iptemp, "pibond40") == 0)
1947     {
1948     if (bondpk.npibond40 > 20 )
1949     {
1950     }
1951     bondpk.use_pibond40 = TRUE;
1952     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
1953     &bondpk.bk40[bondpk.npibond40], &bondpk.bl40[bondpk.npibond40], &bondpk.bmom40[bondpk.npibond40],
1954     &bondpk.sslop40[bondpk.npibond40], &bondpk.tslop40[bondpk.npibond40], &bondpk.tslop240[bondpk.npibond40]);
1955     numeral(ii,pa,3);
1956     numeral(kk,pb,3);
1957     strcpy(pt,pa);
1958     strcat(pt,pb);
1959     strcpy(bondpk.kb40[bondpk.npibond40],pt);
1960     bondpk.npibond40++;
1961     } else if (strcmp(iptemp, "pibond50") == 0)
1962     {
1963     if (bondpk.npibond50 > 20 )
1964     {
1965     }
1966     bondpk.use_pibond50 = TRUE;
1967     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
1968     &bondpk.bk50[bondpk.npibond50], &bondpk.bl50[bondpk.npibond50], &bondpk.bmom50[bondpk.npibond50],
1969     &bondpk.sslop50[bondpk.npibond50], &bondpk.tslop50[bondpk.npibond50], &bondpk.tslop250[bondpk.npibond50]);
1970     numeral(ii,pa,3);
1971     numeral(kk,pb,3);
1972     strcpy(pt,pa);
1973     strcat(pt,pb);
1974     strcpy(bondpk.kb50[bondpk.npibond50],pt);
1975     bondpk.npibond50++;
1976     } else if (strcmp(iptemp, "pibond60") == 0)
1977     {
1978     if (bondpk.npibond60 > 20 )
1979     {
1980     }
1981     bondpk.use_pibond60 = TRUE;
1982     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
1983     &bondpk.bk60[bondpk.npibond60], &bondpk.bl60[bondpk.npibond60], &bondpk.bmom60[bondpk.npibond60],
1984     &bondpk.sslop60[bondpk.npibond60], &bondpk.tslop60[bondpk.npibond60], &bondpk.tslop260[bondpk.npibond60]);
1985     numeral(ii,pa,3);
1986     numeral(kk,pb,3);
1987     strcpy(pt,pa);
1988     strcat(pt,pb);
1989     strcpy(bondpk.kb60[bondpk.npibond60],pt);
1990     bondpk.npibond60++;
1991     } else if (strcmp(iptemp, "pibond44i") == 0)
1992     {
1993     if (bondpk.npibond44i > 20 )
1994     {
1995     }
1996     bondpk.use_pibond44i = TRUE;
1997     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
1998     &bondpk.bk44i[bondpk.npibond44i], &bondpk.bl44i[bondpk.npibond44i], &bondpk.bmom44i[bondpk.npibond44i],
1999     &bondpk.sslop44i[bondpk.npibond44i], &bondpk.tslop44i[bondpk.npibond44i], &bondpk.tslop244i[bondpk.npibond44i]);
2000     numeral(ii,pa,3);
2001     numeral(kk,pb,3);
2002     strcpy(pt,pa);
2003     strcat(pt,pb);
2004     strcpy(bondpk.kb44i[bondpk.npibond44i],pt);
2005     bondpk.npibond44i++;
2006     } else if (strcmp(iptemp, "pibond45i") == 0)
2007     {
2008     if (bondpk.npibond45i > 20 )
2009     {
2010     }
2011     bondpk.use_pibond45i = TRUE;
2012     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2013     &bondpk.bk45i[bondpk.npibond45i], &bondpk.bl45i[bondpk.npibond45i], &bondpk.bmom45i[bondpk.npibond45i],
2014     &bondpk.sslop45i[bondpk.npibond45i], &bondpk.tslop45i[bondpk.npibond45i], &bondpk.tslop245i[bondpk.npibond45i]);
2015     numeral(ii,pa,3);
2016     numeral(kk,pb,3);
2017     strcpy(pt,pa);
2018     strcat(pt,pb);
2019     strcpy(bondpk.kb45i[bondpk.npibond45i],pt);
2020     bondpk.npibond45i++;
2021     } else if (strcmp(iptemp, "pibond46i") == 0)
2022     {
2023     if (bondpk.npibond46i > 20 )
2024     {
2025     }
2026     bondpk.use_pibond46i = TRUE;
2027     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2028     &bondpk.bk46i[bondpk.npibond46i], &bondpk.bl46i[bondpk.npibond46i], &bondpk.bmom46i[bondpk.npibond46i],
2029     &bondpk.sslop46i[bondpk.npibond46i], &bondpk.tslop46i[bondpk.npibond46i], &bondpk.tslop246i[bondpk.npibond46i]);
2030     numeral(ii,pa,3);
2031     numeral(kk,pb,3);
2032     strcpy(pt,pa);
2033     strcat(pt,pb);
2034     strcpy(bondpk.kb46i[bondpk.npibond46i],pt);
2035     bondpk.npibond46i++;
2036     } else if (strcmp(iptemp, "pibond55i") == 0)
2037     {
2038     if (bondpk.npibond55i > 20 )
2039     {
2040     }
2041     bondpk.use_pibond55i = TRUE;
2042     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2043     &bondpk.bk55i[bondpk.npibond55i], &bondpk.bl55i[bondpk.npibond55i], &bondpk.bmom55i[bondpk.npibond55i],
2044     &bondpk.sslop55i[bondpk.npibond55i], &bondpk.tslop55i[bondpk.npibond55i], &bondpk.tslop255i[bondpk.npibond55i]);
2045     numeral(ii,pa,3);
2046     numeral(kk,pb,3);
2047     strcpy(pt,pa);
2048     strcat(pt,pb);
2049     strcpy(bondpk.kb55i[bondpk.npibond55i],pt);
2050     bondpk.npibond55i++;
2051     } else if (strcmp(iptemp, "pibond56i") == 0)
2052     {
2053     if (bondpk.npibond56i > 20 )
2054     {
2055     }
2056     bondpk.use_pibond56i = TRUE;
2057     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2058     &bondpk.bk56i[bondpk.npibond56i], &bondpk.bl56i[bondpk.npibond56i], &bondpk.bmom56i[bondpk.npibond56i],
2059     &bondpk.sslop56i[bondpk.npibond56i], &bondpk.tslop56i[bondpk.npibond56i], &bondpk.tslop256i[bondpk.npibond56i]);
2060     numeral(ii,pa,3);
2061     numeral(kk,pb,3);
2062     strcpy(pt,pa);
2063     strcat(pt,pb);
2064     strcpy(bondpk.kb56i[bondpk.npibond56i],pt);
2065     bondpk.npibond56i++;
2066     } else if (strcmp(iptemp, "pibond66i") == 0)
2067     {
2068     if (bondpk.npibond66i > 20 )
2069     {
2070     }
2071     bondpk.use_pibond66i = TRUE;
2072     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2073     &bondpk.bk66i[bondpk.npibond66i], &bondpk.bl66i[bondpk.npibond66i], &bondpk.bmom66i[bondpk.npibond66i],
2074     &bondpk.sslop66i[bondpk.npibond66i], &bondpk.tslop66i[bondpk.npibond66i], &bondpk.tslop266i[bondpk.npibond66i]);
2075     numeral(ii,pa,3);
2076     numeral(kk,pb,3);
2077     strcpy(pt,pa);
2078     strcat(pt,pb);
2079     strcpy(bondpk.kb66i[bondpk.npibond66i],pt);
2080     bondpk.npibond66i++;
2081     } else if (strcmp(iptemp, "pibond44o") == 0)
2082     {
2083     if (bondpk.npibond44o > 20 )
2084     {
2085     }
2086     bondpk.use_pibond44o = TRUE;
2087     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2088     &bondpk.bk44o[bondpk.npibond44o], &bondpk.bl44o[bondpk.npibond44o], &bondpk.bmom44o[bondpk.npibond44o],
2089     &bondpk.sslop44o[bondpk.npibond44o], &bondpk.tslop44o[bondpk.npibond44o], &bondpk.tslop244o[bondpk.npibond44o]);
2090     numeral(ii,pa,3);
2091     numeral(kk,pb,3);
2092     strcpy(pt,pa);
2093     strcat(pt,pb);
2094     strcpy(bondpk.kb44o[bondpk.npibond44o],pt);
2095     bondpk.npibond44o++;
2096     } else if (strcmp(iptemp, "pibond45o") == 0)
2097     {
2098     if (bondpk.npibond45o > 20 )
2099     {
2100     }
2101     bondpk.use_pibond45o = TRUE;
2102     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2103     &bondpk.bk45o[bondpk.npibond45o], &bondpk.bl45o[bondpk.npibond45o], &bondpk.bmom45o[bondpk.npibond45o],
2104     &bondpk.sslop45o[bondpk.npibond45o], &bondpk.tslop45o[bondpk.npibond45o], &bondpk.tslop245o[bondpk.npibond45o]);
2105     numeral(ii,pa,3);
2106     numeral(kk,pb,3);
2107     strcpy(pt,pa);
2108     strcat(pt,pb);
2109     strcpy(bondpk.kb45o[bondpk.npibond45o],pt);
2110     bondpk.npibond45o++;
2111     } else if (strcmp(iptemp, "pibond46o") == 0)
2112     {
2113     if (bondpk.npibond46o > 20 )
2114     {
2115     }
2116     bondpk.use_pibond46o = TRUE;
2117     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2118     &bondpk.bk46o[bondpk.npibond46o], &bondpk.bl46o[bondpk.npibond46o], &bondpk.bmom46o[bondpk.npibond46o],
2119     &bondpk.sslop46o[bondpk.npibond46o], &bondpk.tslop46o[bondpk.npibond46o], &bondpk.tslop246o[bondpk.npibond46o]);
2120     numeral(ii,pa,3);
2121     numeral(kk,pb,3);
2122     strcpy(pt,pa);
2123     strcat(pt,pb);
2124     strcpy(bondpk.kb46o[bondpk.npibond46o],pt);
2125     bondpk.npibond46o++;
2126     } else if (strcmp(iptemp, "pibond55o") == 0)
2127     {
2128     if (bondpk.npibond55o > 20 )
2129     {
2130     }
2131     bondpk.use_pibond55o = TRUE;
2132     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2133     &bondpk.bk55o[bondpk.npibond55o], &bondpk.bl55o[bondpk.npibond55o], &bondpk.bmom55o[bondpk.npibond55o],
2134     &bondpk.sslop55o[bondpk.npibond55o], &bondpk.tslop55o[bondpk.npibond55o], &bondpk.tslop255o[bondpk.npibond55o]);
2135     numeral(ii,pa,3);
2136     numeral(kk,pb,3);
2137     strcpy(pt,pa);
2138     strcat(pt,pb);
2139     strcpy(bondpk.kb55o[bondpk.npibond55o],pt);
2140     bondpk.npibond55o++;
2141     } else if (strcmp(iptemp, "pibond56o") == 0)
2142     {
2143     if (bondpk.npibond56o > 20 )
2144     {
2145     }
2146     bondpk.use_pibond56o = TRUE;
2147     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2148     &bondpk.bk56o[bondpk.npibond56o], &bondpk.bl56o[bondpk.npibond56o], &bondpk.bmom56o[bondpk.npibond56o],
2149     &bondpk.sslop56o[bondpk.npibond56o], &bondpk.tslop56o[bondpk.npibond56o], &bondpk.tslop256o[bondpk.npibond56o]);
2150     numeral(ii,pa,3);
2151     numeral(kk,pb,3);
2152     strcpy(pt,pa);
2153     strcat(pt,pb);
2154     strcpy(bondpk.kb56o[bondpk.npibond56o],pt);
2155     bondpk.npibond56o++;
2156     } else if (strcmp(iptemp, "pibond66o") == 0)
2157     {
2158     if (bondpk.npibond66o > 20 )
2159     {
2160     }
2161     bondpk.use_pibond66o = TRUE;
2162     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2163     &bondpk.bk66o[bondpk.npibond66o], &bondpk.bl66o[bondpk.npibond66o], &bondpk.bmom66o[bondpk.npibond66o],
2164     &bondpk.sslop66o[bondpk.npibond66o], &bondpk.tslop66o[bondpk.npibond66o], &bondpk.tslop266o[bondpk.npibond66o]);
2165     numeral(ii,pa,3);
2166     numeral(kk,pb,3);
2167     strcpy(pt,pa);
2168     strcat(pt,pb);
2169     strcpy(bondpk.kb66o[bondpk.npibond66o],pt);
2170     bondpk.npibond66o++;
2171     } else if (strcmp(iptemp, "pibond") == 0)
2172     {
2173     if (bondpk.npibond > MAXPIBONDCONST )
2174     {
2175     }
2176     bondpk.use_pibond = TRUE;
2177     sscanf(line,"%s %d %d %f %f %f %f %f %f",dumm, &ii, &kk,
2178     &bondpk.bk[bondpk.npibond], &bondpk.bl[bondpk.npibond], &bondpk.bmom[bondpk.npibond],
2179     &bondpk.sslop[bondpk.npibond], &bondpk.tslop[bondpk.npibond], &bondpk.tslop2[bondpk.npibond]);
2180     numeral(ii,pa,3);
2181     numeral(kk,pb,3);
2182     strcpy(pt,pa);
2183     strcat(pt,pb);
2184     strcpy(bondpk.kb[bondpk.npibond],pt);
2185     bondpk.npibond++;
2186     } else if (strcmp(iptemp, "picorrect") == 0)
2187     {
2188     sscanf(line,"%s %f %f %f %f",dumm,&picorrect.fudge1,&picorrect.fudge2,&picorrect.fudge3,
2189     &picorrect.fudge4);
2190     } else if (strcmp(iptemp,"piangle") == 0)
2191     {
2192     if (angkp1.npiang > MAXPIANGCONST)
2193     {
2194     }
2195     sscanf(line,"%s %d %d %d %f %f %f %f",dumm,&kt1,&kt2,&kt3,&angkp1.pacon[angkp1.npiang],
2196     &angkp1.panat[angkp1.npiang][0],&angkp1.panat[angkp1.npiang][1],&angkp1.panat[angkp1.npiang][1]);
2197     numeral(kt1,pa,3);
2198     numeral(kt2,pb,3);
2199     numeral(kt3,pc,3);
2200     strcpy(pt,pa);
2201     strcat(pt,pb); strcat(pt,pc);
2202     strcpy(angkp1.kpa[angkp1.npiang],pt);
2203     angkp1.npiang++;
2204     } else if (strcmp(iptemp,"dipole") == 0)
2205     {
2206     if (dipole_k.ndipole > MAXBONDCONST)
2207     {
2208     message_alert("Error - Too many bond dipoles","Read Parameter");
2209 tjod 15 fprintf(stderr,"Maximum number of bond dipole constants exceeded. Others will be ignored\n");
2210 tjod 25 //fclose(datafile);
2211 tjod 3 return;
2212     }
2213     sscanf( line,"%s %d %d %f",dumm, &ii, &kk, &radius);
2214     numeral(ii,pa,3);
2215     numeral(kk,pb,3);
2216     strcpy(pt,pa); strcat(pt,pb);
2217     strcpy(dipole_k.kb[dipole_k.ndipole],pt);
2218     dipole_k.bmom[dipole_k.ndipole] = radius;
2219     dipole_k.ndipole++;
2220     } else if (strcmp(iptemp, "charge") == 0)
2221     {
2222     if (charge_k.ncharge > MAXATOMTYPE)
2223     {
2224     message_alert("Error - Too many atom charges","Read Parameter");
2225 tjod 15 fprintf(stderr,"Maximum number of atom charge constants exceeded. Others will be ignored\n");
2226 tjod 25 //fclose(datafile);
2227 tjod 3 return;
2228     }
2229     sscanf( line,"%s %d %f", dumm, &charge_k.type[charge_k.ncharge], &charge_k.charge[charge_k.ncharge]);
2230     charge_k.ncharge++;
2231     } else if (strcmp(iptemp, "oplschrg") == 0)
2232     {
2233     if (oplschrg_k.ncharge > MAXOPLS)
2234     {
2235     message_alert("Error - Too many opls charges","Read Parameter");
2236 tjod 15 fprintf(stderr,"Maximum number of opls charge constants exceeded. Others will be ignored\n");
2237 tjod 25 //fclose(datafile);
2238 tjod 3 return;
2239     }
2240     sscanf(line,"%s %d %f %s\n",dumm,&oplschrg_k.type[oplschrg_k.ncharge],&oplschrg_k.chrg[oplschrg_k.ncharge],
2241     oplschrg_k.chrgstring[oplschrg_k.ncharge]);
2242     oplschrg_k.ncharge++;
2243     } else if (strcmp(iptemp, "amberchrg") == 0)
2244     {
2245     if (amberchrg_k.ncharge > MAXATOMTYPE)
2246     {
2247     message_alert("Error - Too many amber charges","Read Parameter");
2248 tjod 15 fprintf(stderr,"Maximum number of amber charge constants exceeded. Others will be ignored\n");
2249 tjod 25 //fclose(datafile);
2250 tjod 3 return;
2251     }
2252     for (iz=0;iz < 25; iz++) descript[iz] = ' ';
2253    
2254     sscanf(line,"%s %d %f %s %25c",dumm,&amberchrg_k.type[amberchrg_k.ncharge],
2255     &amberchrg_k.chrg[amberchrg_k.ncharge],amberchrg_k.symbol[amberchrg_k.ncharge], descript);
2256     for (iz=24; iz > 0; iz--)
2257     {
2258     if ( isalpha(descript[iz]) != 0)
2259     {
2260     descript[iz+1] = '\0';
2261     break;
2262     }
2263     }
2264     if (strcmp(descript,"Glycine") == 0)
2265     amberchrg_k.res_type[amberchrg_k.ncharge] = GLY;
2266     else if (strcmp(descript,"Alanine") == 0)
2267     amberchrg_k.res_type[amberchrg_k.ncharge] = ALA;
2268     else if (strcmp(descript,"Valine") == 0)
2269     amberchrg_k.res_type[amberchrg_k.ncharge] = VAL;
2270     else if (strcmp(descript,"Leucine") == 0)
2271     amberchrg_k.res_type[amberchrg_k.ncharge] = LEU;
2272     else if (strcmp(descript,"Isoleucine") == 0)
2273     amberchrg_k.res_type[amberchrg_k.ncharge] = ILE;
2274     else if (strcmp(descript,"Serine") == 0)
2275     amberchrg_k.res_type[amberchrg_k.ncharge] = SER;
2276     else if (strcmp(descript,"Threonine") == 0)
2277     amberchrg_k.res_type[amberchrg_k.ncharge] = THR;
2278     else if (strcmp(descript,"Cysteine (-SH") == 0)
2279     amberchrg_k.res_type[amberchrg_k.ncharge] = CYSH;
2280     else if (strcmp(descript,"Cystine (-SS") == 0)
2281     amberchrg_k.res_type[amberchrg_k.ncharge] = CYS;
2282     else if (strcmp(descript,"Proline") == 0)
2283     amberchrg_k.res_type[amberchrg_k.ncharge] = PRO;
2284     else if (strcmp(descript,"Phenylalanine") == 0)
2285     amberchrg_k.res_type[amberchrg_k.ncharge] = PHE;
2286     else if (strcmp(descript,"Tyrosine") == 0)
2287     amberchrg_k.res_type[amberchrg_k.ncharge] = TYR;
2288     else if (strcmp(descript,"Tryptophan") == 0)
2289     amberchrg_k.res_type[amberchrg_k.ncharge] = TRP;
2290     else if (strcmp(descript,"Histidine") == 0)
2291     amberchrg_k.res_type[amberchrg_k.ncharge] = HISP;
2292     else if (strcmp(descript,"Histidine (HD") == 0)
2293     amberchrg_k.res_type[amberchrg_k.ncharge] = HIS;
2294     else if (strcmp(descript,"Histidine (HE") == 0)
2295     amberchrg_k.res_type[amberchrg_k.ncharge] = HISE;
2296     else if (strcmp(descript,"Aspartic Acid") == 0)
2297     amberchrg_k.res_type[amberchrg_k.ncharge] = ASP;
2298     else if (strcmp(descript,"Asparagine") == 0)
2299     amberchrg_k.res_type[amberchrg_k.ncharge] = ASN;
2300     else if (strcmp(descript,"Glutamic Acid") == 0)
2301     amberchrg_k.res_type[amberchrg_k.ncharge] = GLU;
2302     else if (strcmp(descript,"Glutamine") == 0)
2303     amberchrg_k.res_type[amberchrg_k.ncharge] = GLN;
2304     else if (strcmp(descript,"Methionine") == 0)
2305     amberchrg_k.res_type[amberchrg_k.ncharge] = MET;
2306     else if (strcmp(descript,"Lysine") == 0)
2307     amberchrg_k.res_type[amberchrg_k.ncharge] = LYS;
2308     else if (strcmp(descript,"Arginine") == 0)
2309     amberchrg_k.res_type[amberchrg_k.ncharge] = ARG;
2310     else
2311     amberchrg_k.res_type[amberchrg_k.ncharge] = UNK;
2312     amberchrg_k.ncharge++;
2313     } else if (strcmp(iptemp, "mmffchrg") == 0)
2314     {
2315     if (charge_k.ncharge > MAXATOMTYPE)
2316     {
2317     message_alert("Error - Too many atom charges","Read Parameter");
2318 tjod 15 fprintf(stderr,"Maximum number of atom charge constants exceeded. Others will be ignored\n");
2319 tjod 25 //fclose(datafile);
2320 tjod 3 return;
2321     }
2322     sscanf( line,"%s %d %f %f %f", dumm, &i, &charge, &fchrg, &tchrg);
2323     charge_k.type[i] = i;
2324     charge_k.charge[i] = charge;
2325     charge_k.formchrg[i] = fchrg;
2326     charge_k.typechrg[i] = tchrg;
2327     } else if (strcmp(iptemp,"bndchrgdel") == 0)
2328     {
2329     if (charge_k.nbndchrgdel > MAXBONDCONST)
2330     {
2331     message_alert("Error - Too many bond charges","Read Parameter");
2332 tjod 15 fprintf(stderr,"Maximum number of bond charge constants exceeded. Others will be ignored\n");
2333 tjod 25 //fclose(datafile);
2334 tjod 3 return;
2335     }
2336     sscanf(line,"%s %d %d %f",dumm,&kp1,&kp2,&charge_k.bchargedel[charge_k.nbndchrgdel]);
2337     charge_k.btypedel[charge_k.nbndchrgdel] = kp1*100+kp2;
2338     charge_k.nbndchrgdel++;
2339     } else if (strcmp(iptemp,"bndchrg") == 0)
2340     {
2341     if (charge_k.nbndchrg > MAXBONDCONST)
2342     {
2343     message_alert("Error - Too many bond charges","Read Parameter");
2344 tjod 15 fprintf(stderr,"Maximum number of bond charge constants exceeded. Others will be ignored\n");
2345 tjod 25 //fclose(datafile);
2346 tjod 3 return;
2347     }
2348     sscanf(line,"%s %d %d %f",dumm,&kp1,&kp2,&charge_k.bcharge[charge_k.nbndchrg]);
2349     charge_k.btype[charge_k.nbndchrg] = kp1*100+kp2;
2350     charge_k.nbndchrg++;
2351     } else if (strcmp(iptemp,"vdwmmff") == 0)
2352     {
2353     if (vdw1.nvdw > MAXVDWCONST)
2354     {
2355     message_alert("Error - Too many vdw constants","Read Parameter");
2356 tjod 15 fprintf(stderr,"Maximum number of vdw constants exceeded. Others will be ignored\n");
2357 tjod 25 //fclose(datafile);
2358 tjod 3 return;
2359     }
2360     sscanf( line, "%s %d %f %f %f %f %s", dumm, &i, &alpha, &ntmp, &atmp, &gtmp, symbol);
2361     vdw1.alpha[i] = alpha;
2362     vdw1.n[i] = ntmp;
2363     vdw1.a[i] = atmp;
2364     vdw1.g[i] = gtmp;
2365     strcpy(vdw1.da[i],symbol);
2366     vdw1.nvdw++;
2367     } else if( strcmp(iptemp,"vdw") == 0 )
2368     {
2369     if (vdw1.nvdw > MAXVDWCONST)
2370     {
2371     message_alert("Error - Too many vdw constants","Read Parameter");
2372 tjod 15 fprintf(stderr,"Maximum number of vdw constants exceeded. Others will be ignored\n");
2373 tjod 25 //fclose(datafile);
2374 tjod 3 return;
2375     }
2376     sscanf( line, "%s %d %f %f %d %d %d", dumm, &i, &radius, &epsilon, &lpde, &ihtype, &ihdonor);
2377    
2378     vdw1.rad[i] = radius;
2379     vdw1.eps[i] = epsilon;
2380     vdw1.lpd[i] = lpde;
2381     vdw1.ihdon[i] = ihdonor;
2382     vdw1.ihtyp[i] = ihtype;
2383     vdw1.nvdw++;
2384     }
2385     /* torsion read */
2386     else if( strcmp(iptemp,"torsion4") == 0 )
2387     {
2388     if (torkn1.ntor4 >= MAXTOR4CONST)
2389     {
2390     message_alert("Error - Too many torsion 4 constants","Read Parameter");
2391 tjod 15 fprintf(stderr,"Maximum number of torsion 4 constants exceeded. Others will be ignored\n");
2392 tjod 25 //fclose(datafile);
2393 tjod 3 return;
2394     }
2395     for(i=0; i < 6; i++)
2396     v1[i] = se1[i] = 0.0;
2397     iz = sscanf(line, "%s %d %d %d %d %f %d %f %d %f %d %f %d %f %d %f %d", dumm, &kp1,
2398     &kp2, &kp3, &kp4 ,
2399     &v1[0], &se1[0],&v1[1], &se1[1],&v1[2], &se1[2],
2400     &v1[3], &se1[3],&v1[4], &se1[4],&v1[5], &se1[5]);
2401     torphase(6, v1, se1);
2402     torkn1.tv41[torkn1.ntor4] = v1[0];
2403     torkn1.tv42[torkn1.ntor4] = v1[1];
2404     torkn1.tv43[torkn1.ntor4] = v1[2];
2405     torkn1.phase41[torkn1.ntor4] = se1[0];
2406     torkn1.phase42[torkn1.ntor4] = se1[1];
2407     torkn1.phase43[torkn1.ntor4] = se1[2];
2408    
2409     numeral(kp1,pa,3);
2410     numeral(kp2,pb,3);
2411     numeral(kp3,pc,3);
2412     numeral(kp4,pd,3);
2413     strcpy(ptor,pa);
2414     strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2415     strcpy(torkn1.kv4[torkn1.ntor4],ptor);
2416     torkn1.ntor4++;
2417     }else if( strcmp(iptemp,"torsion5") == 0 )
2418     {
2419     if (torkn1.ntor5 >= MAXTOR5CONST)
2420     {
2421     message_alert("Error - Too many torsion 5 constants","Read Parameter");
2422 tjod 15 fprintf(stderr,"Maximum number of torsion 5 constants exceeded. Others will be ignored\n");
2423 tjod 25 //fclose(datafile);
2424 tjod 3 return;
2425     }
2426     for(i=0; i < 6; i++)
2427     v1[i] = se1[i] = 0.0;
2428     iz = sscanf(line, "%s %d %d %d %d %f %d %f %d %f %d %f %d %f %d %f %d", dumm, &kp1,
2429     &kp2, &kp3, &kp4 ,
2430     &v1[0], &se1[0],&v1[1], &se1[1],&v1[2], &se1[2],
2431     &v1[3], &se1[3],&v1[4], &se1[4],&v1[5], &se1[5]);
2432     torphase(6, v1, se1);
2433     torkn1.tv51[torkn1.ntor5] = v1[0];
2434     torkn1.tv52[torkn1.ntor5] = v1[1];
2435     torkn1.tv53[torkn1.ntor5] = v1[2];
2436     torkn1.phase51[torkn1.ntor5] = se1[0];
2437     torkn1.phase52[torkn1.ntor5] = se1[1];
2438     torkn1.phase53[torkn1.ntor5] = se1[2];
2439     numeral(kp1,pa,3);
2440     numeral(kp2,pb,3);
2441     numeral(kp3,pc,3);
2442     numeral(kp4,pd,3);
2443     strcpy(ptor,pa);
2444     strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2445     strcpy(torkn1.kv5[torkn1.ntor5],ptor);
2446     torkn1.ntor5++;
2447     }else if( strcmp(iptemp,"torsiondel") == 0 )
2448     {
2449     if (torkn1.ntordel >= MAXTORDEL)
2450     {
2451     message_alert("Error - Too many delocalized torsion constants","Read Parameter");
2452 tjod 15 fprintf(stderr,"Maximum number of delocalized torsion constants exceeded. Others will be ignored\n");
2453 tjod 25 //fclose(datafile);
2454 tjod 3 return;
2455     }
2456     for(i=0; i < 6; i++)
2457     v1[i] = se1[i] = 0.0;
2458     iz = sscanf(line, "%s %d %d %d %d %d %f %d %f %d %f %d %f %d %f %d %f %d", dumm, &torkn1.torindex[torkn1.ntordel],
2459     &kp1, &kp2, &kp3, &kp4 ,
2460     &v1[0], &se1[0],&v1[1], &se1[1],&v1[2], &se1[2],
2461     &v1[3], &se1[3],&v1[4], &se1[4],&v1[5], &se1[5]);
2462     torphase(6, v1, se1);
2463     torkn1.tvdel1[torkn1.ntordel] = v1[0];
2464     torkn1.tvdel2[torkn1.ntordel] = v1[1];
2465     torkn1.tvdel3[torkn1.ntordel] = v1[2];
2466     torkn1.phasedel1[torkn1.ntordel] = se1[0];
2467     torkn1.phasedel2[torkn1.ntordel] = se1[1];
2468     torkn1.phasedel3[torkn1.ntordel] = se1[2];
2469     numeral(kp1,pa,3);
2470     numeral(kp2,pb,3);
2471     numeral(kp3,pc,3);
2472     numeral(kp4,pd,3);
2473     strcpy(ptor,pa);
2474     strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2475     strcpy(torkn1.kvdel[torkn1.ntordel],ptor);
2476     torkn1.ntordel++;
2477     }else if( strcmp(iptemp, "torsion") == 0 )
2478     {
2479     if (torkn1.ntor >= MAXTORCONST)
2480     {
2481     message_alert("Error - Too many torsion constants","Read Parameter");
2482 tjod 15 fprintf(stderr,"Maximum number of torsion constants exceeded. Others will be ignored\n");
2483 tjod 25 //fclose(datafile);
2484 tjod 3 return;
2485     }
2486     for(i=0; i < 6; i++)
2487     v1[i] = se1[i] = 0.0;
2488    
2489     sscanf( line, "%s %d %d %d %d %f %d %f %d %f %d %f %d %f %d %f %d", dumm,
2490     &kp1, &kp2, &kp3, &kp4,
2491     &v1[0], &se1[0],&v1[1], &se1[1],&v1[2], &se1[2],
2492     &v1[3], &se1[3],&v1[4], &se1[4],&v1[5], &se1[5]);
2493     torphase(6, v1, se1);
2494     torkn1.tv1[torkn1.ntor] = v1[0];
2495     torkn1.phase1[torkn1.ntor] = se1[0];
2496     torkn1.tv2[torkn1.ntor] = v1[1];
2497     torkn1.phase2[torkn1.ntor] = se1[1];
2498     torkn1.tv3[torkn1.ntor] = v1[2];
2499     torkn1.phase3[torkn1.ntor] = se1[2];
2500     torkn1.tv4[torkn1.ntor] = v1[3];
2501     torkn1.phase4[torkn1.ntor] = se1[3];
2502     torkn1.tv5[torkn1.ntor] = v1[4];
2503     torkn1.phase5[torkn1.ntor] = se1[4];
2504     torkn1.tv6[torkn1.ntor] = v1[5];
2505     torkn1.phase6[torkn1.ntor] = se1[5];
2506     numeral(kp1,pa,3);
2507     numeral(kp2,pb,3);
2508     numeral(kp3,pc,3);
2509     numeral(kp4,pd,3);
2510     strcpy(ptor,pa);
2511     strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2512     strcpy(torkn1.kv[torkn1.ntor],ptor);
2513    
2514     // torkn1.kv[torkn1.ntor] = ((kp1*100 + kp2)*100 + kp3)*100 + kp4;
2515     torkn1.ntor++;
2516     }else if (strcmp(iptemp, "pitorsion") == 0)
2517     {
2518     if (torknp.npitor > MAXPITORCONST)
2519     {
2520     }
2521     sscanf( line,"%s %d %d %d %d %f %d %f %d %f %d", dumm, &kp1, &kp2, &kp3, &kp4,
2522     &v1[0], &se1[0], &v1[1], &se1[1], &v1[2], &se1[2]);
2523    
2524     torphase(3, v1,se1);
2525     torknp.tv1[torknp.npitor] = v1[0];
2526     torknp.ph1[torknp.npitor] = se1[0];
2527     torknp.tv2[torknp.npitor] = v1[1];
2528     torknp.ph2[torknp.npitor] = se1[1];
2529     torknp.tv3[torknp.npitor] = v1[2];
2530     torknp.ph3[torknp.npitor] = se1[2];
2531     numeral(kp1,pa,3);
2532     numeral(kp2,pb,3);
2533     numeral(kp3,pc,3);
2534     numeral(kp4,pd,3);
2535     strcpy(ptor,pa);
2536     strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2537     strcpy(torknp.kv[torknp.npitor],ptor);
2538     // torknp.kv[torknp.npitor] = ((kp1*100 + kp2)*100 + kp3)*100 + kp4;
2539     torknp.npitor++;
2540     } else if( strcmp(iptemp, "imptors") == 0)
2541     {
2542     if (improptor_k.nimptor >= MAXIMP)
2543     {
2544     message_alert("Error - Too many improper torsion constants","Read Parameter");
2545 tjod 15 fprintf(stderr,"Maximum number of improper torsion constants exceeded. Others will be ignored\n");
2546 tjod 25 //fclose(datafile);
2547 tjod 3 return;
2548     }
2549     sscanf( line, "%s %d %d %d %d %f %d %f %d %f %d", dumm,
2550     &kp1, &kp2, &kp3, &kp4,
2551     &v1[0], &se1[0], &v1[1], &se1[1], &v1[2], &se1[2]);
2552     torphase(3,v1,se1);
2553     improptor_k.v1[improptor_k.nimptor] = v1[0];
2554     improptor_k.v2[improptor_k.nimptor] = v1[1];
2555     improptor_k.v3[improptor_k.nimptor] = v1[2];
2556     improptor_k.ph1[improptor_k.nimptor] = se1[0];
2557     improptor_k.ph2[improptor_k.nimptor] = se1[1];
2558     improptor_k.ph3[improptor_k.nimptor] = se1[2];
2559     numeral(kp1,pa,3);
2560     numeral(kp2,pb,3);
2561     numeral(kp3,pc,3);
2562     numeral(kp4,pd,3);
2563     strcpy(ptor,pa);
2564     strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2565     strcpy(improptor_k.kv[improptor_k.nimptor],ptor);
2566     // improptor_k.kv[improptor_k.nimptor] = ((kp1*100 + kp2)*100 + kp3)*100 + kp4;
2567     improptor_k.nimptor++;
2568    
2569     } else if( strcmp(iptemp, "improper") == 0)
2570     {
2571     if (improptor_k.nimptor >= MAXIMP)
2572     {
2573     message_alert("Error - Too many improper torsion constants","Read Parameter");
2574 tjod 15 fprintf(stderr,"Maximum number of improper torsion constants exceeded. Others will be ignored\n");
2575 tjod 25 //fclose(datafile);
2576 tjod 3 return;
2577     }
2578     sscanf( line, "%s %d %d %d %d %f %f", dumm,
2579     &kp1, &kp2, &kp3, &kp4,
2580     &improptor_k.cimptor[improptor_k.nimptor], &improptor_k.tdi[improptor_k.nimptor]);
2581     numeral(kp1,pa,3);
2582     numeral(kp2,pb,3);
2583     numeral(kp3,pc,3);
2584     numeral(kp4,pd,3);
2585     strcpy(ptor,pa);
2586     strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2587     strcpy(improptor_k.kv[improptor_k.nimptor],ptor);
2588     improptor_k.nimptor++;
2589     } else if( strcmp(iptemp, "ureybrad") == 0)
2590     {
2591     if (ureybrad_k.nurey_brad > MAXUREY)
2592     {
2593     message_alert("Error - Too many UreyBrad constants","Read Parameter");
2594 tjod 15 fprintf(stderr,"Maximum number of ureybradley constants exceeded. Others will be ignored\n");
2595 tjod 25 //fclose(datafile);
2596 tjod 3 return;
2597     }
2598     sscanf(line, "%s %d %d %d %f %f",dumm,
2599     &kp1, &kp2, &kp3, &ureybrad_k.ubconst[ureybrad_k.nurey_brad],
2600     &ureybrad_k.ubdist[ureybrad_k.nurey_brad]);
2601     numeral(kp1,pa,3);
2602     numeral(kp2,pb,3);
2603     numeral(kp3,pc,3);
2604     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2605     strcpy(ureybrad_k.kubang[ureybrad_k.nurey_brad],pang);
2606     ureybrad_k.nurey_brad++;
2607     } else if (strcmp(iptemp, "anglef") == 0)
2608     {
2609     sscanf(line,"%s %d %d %d %f %f %f %f %f %f %f %f",dumm, &kt1,&kt2,&kt3, &angf.fcon[angf.nfang],
2610     &angf.fc0[angf.nfang],&angf.fc1[angf.nfang],&angf.fc2[angf.nfang],&angf.fc3[angf.nfang],
2611     &angf.fc4[angf.nfang],&angf.fc5[angf.nfang],&angf.fc6[angf.nfang]);
2612     numeral(kt1,pa,3);
2613     numeral(kt2,pb,3);
2614     numeral(kt3,pc,3);
2615     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2616     strcpy(angf.kftype[angf.nfang],pang);
2617     angf.nfang++;
2618     } else if( strcmp(iptemp, "angle5") == 0 )
2619     {
2620     if (angk1.nang5 >= MAXANG5CONST)
2621     {
2622     message_alert("Error - Too many angle 5 constants","Read Parameter");
2623 tjod 15 fprintf(stderr,"Maximum number of angle 5 constants exceeded. Others will be ignored\n");
2624 tjod 25 //fclose(datafile);
2625 tjod 3 return;
2626     }
2627     sscanf( line, "%s %d %d %d %f %f %f %f", dumm,
2628     &kt1, &kt2, &kt3, &angk1.con5[angk1.nang5], &angk1.ang5[angk1.nang5][0],
2629     &angk1.ang5[angk1.nang5][1], &angk1.ang5[angk1.nang5][2] );
2630     numeral(kt1,pa,3);
2631     numeral(kt2,pb,3);
2632     numeral(kt3,pc,3);
2633     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2634     strcpy(angk1.ktype5[angk1.nang5],pang);
2635     angk1.nang5++;
2636     } else if( strcmp(iptemp, "angle4") == 0 )
2637     {
2638     if (angk1.nang4 >= MAXANG4CONST)
2639     {
2640     message_alert("Error - Too many angle 4 constants","Read Parameter");
2641 tjod 15 fprintf(stderr,"Maximum number of angle 4 constants exceeded. Others will be ignored\n");
2642 tjod 25 //fclose(datafile);
2643 tjod 3 return;
2644     }
2645     sscanf( line, "%s %d %d %d %f %f %f %f", dumm,
2646     &kt1, &kt2, &kt3, &angk1.con4[angk1.nang4], &angk1.ang4[angk1.nang4][0],
2647     &angk1.ang4[angk1.nang4][1], &angk1.ang4[angk1.nang4][2] );
2648     numeral(kt1,pa,3);
2649     numeral(kt2,pb,3);
2650     numeral(kt3,pc,3);
2651     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2652     strcpy(angk1.ktype4[angk1.nang4],pang);
2653     angk1.nang4++;
2654     } else if( strcmp(iptemp, "angle3") == 0 )
2655     {
2656     if (angk1.nang3 >= MAXANG3CONST )
2657     {
2658     message_alert("Error - Too many angle 3 constants","Read Parameter");
2659 tjod 15 fprintf(stderr,"Maximum number of angle 3 constants exceeded. Others will be ignored\n");
2660 tjod 25 //fclose(datafile);
2661 tjod 3 return;
2662     }
2663     sscanf( line, "%s %d %d %d %f %f %f %f", dumm,
2664     &kt1, &kt2, &kt3, &angk1.con3[angk1.nang3], &angk1.ang3[angk1.nang3][0],
2665     &angk1.ang3[angk1.nang3][1], &angk1.ang3[angk1.nang3][2] );
2666     numeral(kt1,pa,3);
2667     numeral(kt2,pb,3);
2668     numeral(kt3,pc,3);
2669     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2670     strcpy(angk1.ktype3[angk1.nang3],pang);
2671     angk1.nang3++;
2672     } else if( strcmp(iptemp, "angle") == 0 )
2673     {
2674     if (angk1.nang >= MAXANGCONST)
2675     {
2676     message_alert("Error - Too many angle constants","Read Parameter");
2677 tjod 15 fprintf(stderr,"Maximum number of angle constants exceeded. Others will be ignored\n");
2678 tjod 25 //fclose(datafile);
2679 tjod 3 return;
2680     }
2681     sscanf( line, "%s %d %d %d %f %f %f %f", dumm,
2682     &kt1, &kt2, &kt3, &angk1.con[angk1.nang], &angk1.ang[angk1.nang][0],
2683     &angk1.ang[angk1.nang][1], &angk1.ang[angk1.nang][2] );
2684     numeral(kt1,pa,3);
2685     numeral(kt2,pb,3);
2686     numeral(kt3,pc,3);
2687     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2688     strcpy(angk1.ktype[angk1.nang],pang);
2689     angk1.nang++;
2690     } else if( strcmp(iptemp, "angdel") == 0 )
2691     {
2692     if (angk1.ndel >= MAXANGDEL)
2693     {
2694     message_alert("Error - Too many delocalized angle constants","Read Parameter");
2695 tjod 15 fprintf(stderr,"Maximum number of delocalized angle constants exceeded. Others will be ignored\n");
2696 tjod 25 //fclose(datafile);
2697 tjod 3 return;
2698     }
2699     sscanf( line, "%s %d %d %d %f %f %f %f", dumm,
2700     &kt1, &kt2, &kt3, &angk1.condel[angk1.ndel], &angk1.angdel[angk1.ndel][0],
2701     &angk1.angdel[angk1.ndel][1], &angk1.angdel[angk1.ndel][2] );
2702     numeral(kt1,pa,3);
2703     numeral(kt2,pb,3);
2704     numeral(kt3,pc,3);
2705     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2706     strcpy(angk1.kdel[angk1.ndel],pang);
2707     angk1.ndel++;
2708     } else if( strcmp(iptemp, "ang3del") == 0 )
2709     {
2710     if (angk1.ndel3 >= MAXANG3DEL)
2711     {
2712     message_alert("Error - Too many delocalized angle constants","Read Parameter");
2713 tjod 15 fprintf(stderr,"Maximum number of delocalized angle constants exceeded. Others will be ignored\n");
2714 tjod 25 //fclose(datafile);
2715 tjod 3 return;
2716     }
2717     sscanf( line, "%s %d %d %d %f %f %f %f", dumm,
2718     &kt1, &kt2, &kt3, &angk1.condel3[angk1.ndel3], &angk1.angdel3[angk1.ndel3][0],
2719     &angk1.angdel3[angk1.ndel3][1], &angk1.angdel3[angk1.ndel3][2] );
2720     numeral(kt1,pa,3);
2721     numeral(kt2,pb,3);
2722     numeral(kt3,pc,3);
2723     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2724     strcpy(angk1.kdel3[angk1.ndel3],pang);
2725     angk1.ndel3++;
2726     } else if( strcmp(iptemp, "ang4del") == 0 )
2727     {
2728     if (angk1.ndel4 >= MAXANG4DEL)
2729     {
2730     message_alert("Error - Too many delocalized angle constants","Read Parameter");
2731 tjod 15 fprintf(stderr,"Maximum number of delocalized angle constants exceeded. Others will be ignored\n");
2732 tjod 25 //fclose(datafile);
2733 tjod 3 return;
2734     }
2735     sscanf( line, "%s %d %d %d %f %f %f %f", dumm,
2736     &kt1, &kt2, &kt3, &angk1.condel4[angk1.ndel4], &angk1.angdel4[angk1.ndel4][0],
2737     &angk1.angdel4[angk1.ndel4][1], &angk1.angdel4[angk1.ndel4][2] );
2738     numeral(kt1,pa,3);
2739     numeral(kt2,pb,3);
2740     numeral(kt3,pc,3);
2741     strcpy(pang,pa); strcat(pang,pb); strcat(pang,pc);
2742     strcpy(angk1.kdel4[angk1.ndel4],pang);
2743     angk1.ndel4++;
2744     } else if( strcmp(iptemp, "opbend") == 0 )
2745     {
2746     if (ooplane_k.nopbend >= MAXOOP)
2747     {
2748     message_alert("Error - Too many out of plane bending constants","Read Parameter");
2749 tjod 15 fprintf(stderr,"Maximum number of out of plane bend constants exceeded. Others will be ignored\n");
2750 tjod 25 //fclose(datafile);
2751 tjod 3 return;
2752     }
2753     sscanf( line, "%s %d %d %d %d %f", dumm, &kp1, &kp2, &kp3,&kp4,
2754     &ooplane_k.copb[ooplane_k.nopbend] );
2755     numeral(kp1,pa,3);
2756     numeral(kp2,pb,3);
2757     numeral(kp3,pc,3);
2758     numeral(kp4,pd,3);
2759     strcpy(ptor,pa); strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2760     strcpy(ooplane_k.iopb[ooplane_k.nopbend],ptor);
2761     ooplane_k.nopbend++;
2762     } else if (strcmp(iptemp, "strbnd") == 0)
2763     {
2764     if (crossterm_k.nstrbnd > MAXSTBN)
2765     {
2766     message_alert("Error - Too many str bend constants","Read Parameter");
2767 tjod 15 fprintf(stderr,"Maximum number of str bend constants exceeded. Others will be ignored\n");
2768 tjod 25 //fclose(datafile);
2769 tjod 3 return;
2770     }
2771     kp4 = 0;
2772     sscanf(line,"%s %d %d %d %f %f %f %d",dumm, &kp1, &kp2, &kp3,
2773     &crossterm_k.stbncon[crossterm_k.nstrbnd][0],
2774     &crossterm_k.stbncon[crossterm_k.nstrbnd][1],
2775     &crossterm_k.stbncon[crossterm_k.nstrbnd][2], &kp4);
2776     numeral(kp1,pa,3);
2777     numeral(kp2,pb,3);
2778     numeral(kp3,pc,3);
2779     strcpy(ptor,pa); strcat(ptor,pb); strcat(ptor,pc);
2780     strcpy(crossterm_k.stbn[crossterm_k.nstrbnd],ptor);
2781     crossterm_k.stbnindex[crossterm_k.nstrbnd] = kp4;
2782     crossterm_k.nstrbnd++;
2783     } else if (strcmp(iptemp, "angang") == 0)
2784     {
2785     if (crossterm_k.nangang > MAXAA)
2786     {
2787     message_alert("Error - Too many ang ang constants","Read Parameter");
2788 tjod 15 fprintf(stderr,"Maximum number of ang ang constants exceeded. Others will be ignored\n");
2789 tjod 25 //fclose(datafile);
2790 tjod 3 return;
2791     }
2792     sscanf(line,"%s %d %f %f %f",dumm, &crossterm_k.ang_ang[crossterm_k.nangang],
2793     &crossterm_k.aacon[crossterm_k.nangang][0],
2794     &crossterm_k.aacon[crossterm_k.nangang][1],
2795     &crossterm_k.aacon[crossterm_k.nangang][2]);
2796     crossterm_k.nangang++;
2797     } else if (strcmp(iptemp, "strtors") == 0)
2798     {
2799     if (crossterm_k.nstrtor > MAXSTRTOR)
2800     {
2801     message_alert("Error - Too many str tor constants","Read Parameter");
2802 tjod 15 fprintf(stderr,"Maximum number of str tor constants exceeded. Others will be ignored\n");
2803 tjod 25 //fclose(datafile);
2804 tjod 3 return;
2805     }
2806     sscanf(line,"%s %d %d %f",dumm, &kt1, &kt2, &crossterm_k.str_torcon[crossterm_k.nstrtor]);
2807     numeral(kt1,pa,3);
2808     numeral(kt2,pb,3);
2809     strcpy(pt,pa); strcat(pt,pb);
2810     strcpy(crossterm_k.str_tor[crossterm_k.nstrtor],pt);
2811     crossterm_k.nstrtor++;
2812     } else if (strcmp(iptemp, "vdwpr") == 0)
2813     {
2814     if (vdwpr_k .nvdwpr > MAXBONDCONST)
2815     {
2816     message_alert("Error - Too many vdwpr constants","Read Parameter");
2817 tjod 15 fprintf(stderr,"Maximum number of vdwpr constants exceeded. Others will be ignored\n");
2818 tjod 25 //fclose(datafile);
2819 tjod 3 return;
2820     }
2821     sscanf(line,"%s %d %d %f %f",dumm, &ii, &kk, &vdwpr_k.radius[vdwpr_k .nvdwpr],
2822     &vdwpr_k.eps[vdwpr_k .nvdwpr]);
2823     vdwpr_k.ia1[vdwpr_k .nvdwpr] = ii;
2824     vdwpr_k.ia2[vdwpr_k .nvdwpr] = kk;
2825     numeral(ii,pa,3);
2826     numeral(kk,pb,3);
2827     strcpy(pt,pa); strcat(pt,pb);
2828     strcpy(vdwpr_k.kv[vdwpr_k .nvdwpr],pt);
2829     vdwpr_k .nvdwpr++;
2830     } else if( strcmp(iptemp, "BDHF") == 0 )
2831     {
2832     if (ehpara.neheat >= 220)
2833     {
2834 tjod 15 fprintf(stderr,"Maximum number of BDHF constants exceeded. Others will be ignored\n");
2835 tjod 25 //fclose(datafile);
2836 tjod 3 return;
2837     }
2838     sscanf( line, "%s %d %d %20c %f %f %f", dumm, &ja, &jb,
2839     ehpara.cc[ehpara.neheat], &ehpara.ee[ehpara.neheat], &ehpara.estr[ehpara.neheat],
2840     &ehpara.et[ehpara.neheat] );
2841     ehpara.nn[ehpara.neheat] = ja*100 + jb;
2842     ehpara.neheat++;
2843     } else if( strcmp(iptemp, "EVHF") == 0 )
2844     {
2845     if (ehpara.nevhf >= 155)
2846     {
2847 tjod 15 fprintf(stderr,"Maximum number of EVHF constants exceeded. Others will be ignored\n");
2848 tjod 25 //fclose(datafile);
2849 tjod 3 return;
2850     }
2851     sscanf( line, "%s %d %d %d %d %d %20c %f %f", dumm, &ehpara.cent[ehpara.nevhf],
2852     &kp1,&kp2,&kp3,&kp4, ehpara.mm[ehpara.nevhf],
2853     &ehpara.ss[ehpara.nevhf], &ehpara.ssless[ehpara.nevhf]);
2854     numeral(kp1,pa,3);
2855     numeral(kp2,pb,3);
2856     numeral(kp3,pc,3);
2857     numeral(kp4,pd,3);
2858     strcpy(ptor,pa); strcat(ptor,pb); strcat(ptor,pc); strcat(ptor,pd);
2859     strcpy(ehpara.kheat[ehpara.nevhf],ptor);
2860     ehpara.nevhf++;
2861     } else if( strcmp(iptemp, "PIHF") == 0 )
2862     {
2863     if (epiheat.npihf >= 125)
2864     {
2865 tjod 15 fprintf(stderr,"Maximum number of PIHF constants exceeded. Others will be ignored\n");
2866 tjod 25 //fclose(datafile);
2867 tjod 3 return;
2868     }
2869     sscanf( line, "%s %d %d %20c %f %f %f %f ", dumm,
2870     &iatm1, &iatm2, epiheat.ccp[epiheat.npihf], &epiheat.eep[epiheat.npihf],
2871     &epiheat.aa[epiheat.npihf], &epiheat.bb[epiheat.npihf], &epiheat.ccc[epiheat.npihf] );
2872     epiheat.nnp[epiheat.npihf] = iatm1*100 + iatm2;
2873     epiheat.npihf++;
2874     }
2875     }
2876    
2877 tjod 25 //fclose(datafile);
2878 tjod 3
2879 tjod 15 /*
2880     fprintf(stderr," field : %s\n",field.name);
2881     fprintf(stderr," Atom Types: %d\n",atom_k.natomtype);
2882     fprintf(stderr," Bonds: %d Bond3: %d Bond4: %d Bond5: %d\n",
2883 tjod 3 bondk1.nbnd,bondk1.nbnd3, bondk1.nbnd4, bondk1.nbnd5);
2884 tjod 15 fprintf(stderr," Angle: %d Angle3: %d Angle4: %d Angle5: %d\n",
2885 tjod 3 angk1.nang,angk1.nang3, angk1.nang4, angk1.nang5);
2886 tjod 15 fprintf(stderr," Torsion: %d Torsion4: %d Torsion5: %d\n",
2887 tjod 3 torkn1.ntor, torkn1.ntor4, torkn1.ntor5);
2888 tjod 15 fprintf(stderr," Vdw: %d OOP: %d Dipole: %d Charge: %d Improper: %d\n",
2889 tjod 3 vdw1.nvdw,ooplane_k.nopbend, dipole_k.ndipole, charge_k.ncharge, improptor_k.nimptor);
2890 tjod 15 fprintf(stderr," STBN: %d ANGANG: %d STRTOR: %d VDWPR: %d\n",
2891 tjod 3 crossterm_k.nstrbnd, crossterm_k.nangang, crossterm_k.nstrtor, vdwpr_k.nvdwpr);
2892 tjod 15 */
2893 tjod 3
2894     return;
2895     }
2896    
2897    
2898     void torphase(int icount, float v1[6], int se[6])
2899     {
2900     int i;
2901     float amp[6], phase[6];
2902     int fold[6];
2903    
2904     for (i=0; i < 6; i++)
2905     {
2906     fold[i] = 0;
2907     amp[i] = phase[i] = 0.0;
2908     }
2909    
2910     for (i=0; i < icount; i++)
2911     {
2912     amp[i] = v1[i];
2913     fold[i] = abs( (int)(2.0*se[i]) - (int) se[i]);
2914     if (se[i] > 0.0)
2915     phase[i] = 1.0;
2916     else
2917     phase[i] = -1.0;
2918     v1[i] = 0.0;
2919     se[i] = 0;
2920     }
2921    
2922     for (i=0; i < icount; i++)
2923     {
2924     if (fold[i] != 0 && fold[i] <= icount)
2925     {
2926     v1[fold[i]-1] = amp[i];
2927     se[fold[i]-1] = phase[i];
2928     }
2929     }
2930     }
2931    
2932