4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
11 * Copyright (c) 1991-2001, University of Groningen, The Netherlands
12 * This program is free software; you can redistribute it and/or
13 * modify it under the terms of the GNU General Public License
14 * as published by the Free Software Foundation; either version 2
15 * of the License, or (at your option) any later version.
17 * If you want to redistribute modifications, please consider that
18 * scientific software is very special. Version control is crucial -
19 * bugs must be traceable. We will be happy to consider code for
20 * inclusion in the official distribution, but derived work must not
21 * be called official GROMACS. Details are found in the README & COPYING
22 * files - if they are missing, get the official version at www.gromacs.org.
24 * To help us fund GROMACS development, we humbly ask that you cite
25 * the papers on the package - you can find them in the top README file.
27 * For more info, check our website at http://www.gromacs.org
30 * Gromacs Runs One Microsecond At Cannonball Speeds
32 static char *SRCID_fitahx_c
= "$Id$";
38 static void my_calc_xcm(int nbb
,atom_id bbind
[],rvec x
[],rvec xcm
)
43 for(i
=0; (i
<nbb
); i
++) {
47 for(m
=0; (m
<DIM
); m
++)
51 static void my_sub_xcm(int nbb
,atom_id bbind
[],rvec x
[],rvec xcm
)
55 for(i
=0; (i
<nbb
); i
++) {
61 real
fit_ahx(int nres
,t_bb bb
[],int natoms
,int nall
,atom_id allindex
[],
63 atom_id caindex
[],matrix box
,bool bFit
)
65 static rvec
*xref
=NULL
;
66 static real
*mass
=NULL
;
67 const real d
=0.15; /* Rise per residue (nm) */
68 const real tw
=1.745; /* Twist per residue (rad) */
69 const real rad
=0.23; /* Radius of the helix (nm) */
75 fatal_error(0,"Need at least 3 Calphas to fit to, (now %d)...\n",nca
);
82 for(i
=0; (i
<nca
); i
++) {
84 xref
[ai
][XX
]=rad
*cos(phi0
);
85 xref
[ai
][YY
]=rad
*sin(phi0
);
88 /* Set the mass to select that this Calpha contributes to fitting */
92 fprintf(stderr
,"%5d %8.3f %8.3f %8.3f %8.3f %8.3f %8.3f\n",ai
,
93 x
[ai
][XX
],x
[ai
][YY
],x
[ai
][ZZ
],
94 xref
[ai
][XX
],xref
[ai
][YY
],xref
[ai
][ZZ
]);
99 /* Center the referece around the origin */
100 my_calc_xcm(nca
,caindex
,xref
,xcm
);
101 my_sub_xcm(nca
,caindex
,xref
,xcm
);
104 /* Now center the to-be-fitted coords around the origin */
105 my_calc_xcm(nca
,caindex
,x
,xcm
);
106 my_sub_xcm(nall
,allindex
,x
,xcm
);
109 dump_ahx(nres
,bb
,xref
,box
,0);
113 for(i
=0; (i
<natoms
); i
++)
117 fatal_error(0,"nmass=%d, nca=%d\n",nmass
,nca
);
119 /* Now call the fitting routine */
121 do_fit(natoms
,mass
,xref
,x
);
123 /* Reset masses and calc rms */
125 for(i
=0; (i
<nres
); i
++) {
128 if (mass
[ai
] > 0.0) {
129 rvec_sub(x
[ai
],xref
[ai
],dx
);
131 bb
[i
].rmsa
+=sqrt(rms
);
137 return sqrt(trms
/nca
);