3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Green Red Orange Magenta Azure Cyan Skyblue
44 static void my_calc_xcm(int nbb, atom_id bbind[], rvec x[], rvec xcm)
49 for (i = 0; (i < nbb); i++)
54 for (m = 0; (m < DIM); m++)
60 static void my_sub_xcm(int nbb, atom_id bbind[], rvec x[], rvec xcm)
64 for (i = 0; (i < nbb); i++)
71 real fit_ahx(int nres, t_bb bb[], int natoms, int nall, atom_id allindex[],
73 atom_id caindex[], matrix box, gmx_bool bFit)
75 static rvec *xref = NULL;
76 static real *mass = NULL;
77 const real d = 0.15; /* Rise per residue (nm) */
78 const real tw = 1.745; /* Twist per residue (rad) */
79 const real rad = 0.23; /* Radius of the helix (nm) */
86 gmx_fatal(FARGS, "Need at least 3 Calphas to fit to, (now %d)...\n", nca);
95 for (i = 0; (i < nca); i++)
98 xref[ai][XX] = rad*cos(phi0);
99 xref[ai][YY] = rad*sin(phi0);
102 /* Set the mass to select that this Calpha contributes to fitting */
106 fprintf(stderr, "%5d %8.3f %8.3f %8.3f %8.3f %8.3f %8.3f\n", ai,
107 x[ai][XX], x[ai][YY], x[ai][ZZ],
108 xref[ai][XX], xref[ai][YY], xref[ai][ZZ]);
113 /* Center the referece around the origin */
114 my_calc_xcm(nca, caindex, xref, xcm);
115 my_sub_xcm(nca, caindex, xref, xcm);
119 /* Now center the to-be-fitted coords around the origin */
120 my_calc_xcm(nca, caindex, x, xcm);
121 my_sub_xcm(nall, allindex, x, xcm);
124 dump_ahx(nres, bb, xref, box, 0);
128 for (i = 0; (i < natoms); i++)
137 gmx_fatal(FARGS, "nmass=%d, nca=%d\n", nmass, nca);
140 /* Now call the fitting routine */
143 do_fit(natoms, mass, xref, x);
146 /* Reset masses and calc rms */
148 for (i = 0; (i < nres); i++)
154 rvec_sub(x[ai], xref[ai], dx);
156 bb[i].rmsa += sqrt(rms);
162 return sqrt(trms/nca);