/*
+ * This file is part of the GROMACS molecular simulation package.
*
- * This source code is part of
+ * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
*
- * G R O M A C S
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
*
- * GROningen MAchine for Chemical Simulations
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
*
- * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
- * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2009, The GROMACS development team,
- * check out http://www.gromacs.org for more information.
-
- * This program is free software; you can redistribute it and/or
- * modify it under the terms of the GNU General Public License
- * as published by the Free Software Foundation; either version 2
- * of the License, or (at your option) any later version.
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
*
- * If you want to redistribute modifications, please consider that
- * scientific software is very special. Version control is crucial -
- * bugs must be traceable. We will be happy to consider code for
- * inclusion in the official distribution, but derived work must not
- * be called official GROMACS. Details are found in the README & COPYING
- * files - if they are missing, get the official version at www.gromacs.org.
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
*
* To help us fund GROMACS development, we humbly ask that you cite
- * the papers on the package - you can find them in the top README file.
- *
- * For more info, check our website at http://www.gromacs.org
+ * the research papers on the package. Check out http://www.gromacs.org.
*/
-/*! \page page_module_selection_insolidangle Selection method: insolidangle
+/*! \internal
+ * \page page_module_selection_insolidangle Selection method: insolidangle
*
* This method selects a subset of particles that are located in a solid
* angle defined by a center and a set of points.
* point is in the solid angle if it lies within any of these cones.
* The width of the cones can be adjusted.
*
- * \internal
- *
* The method is implemented by partitioning the surface of the unit sphere
* into bins using the polar coordinates \f$(\theta, \phi)\f$.
* The partitioning is always uniform in the zenith angle \f$\theta\f$,
* (along the lines it is handled in selection.h and trajana.h in the old C
* API).
*
- * \author Teemu Murtola <teemu.murtola@cbr.su.se>
+ * \author Teemu Murtola <teemu.murtola@gmail.com>
* \ingroup module_selection
*/
-#include <algorithm>
+#include "gmxpre.h"
#include <math.h>
-#include "gromacs/legacyheaders/macros.h"
-#include "gromacs/legacyheaders/maths.h"
-#include "gromacs/legacyheaders/pbc.h"
-#include "gromacs/legacyheaders/physics.h"
-#include "gromacs/legacyheaders/smalloc.h"
-#include "gromacs/legacyheaders/vec.h"
+#include <algorithm>
+#include "gromacs/legacyheaders/macros.h"
+#include "gromacs/math/units.h"
+#include "gromacs/math/utilities.h"
+#include "gromacs/math/vec.h"
+#include "gromacs/pbcutil/pbc.h"
#include "gromacs/selection/indexutil.h"
#include "gromacs/selection/position.h"
#include "gromacs/selection/selection.h"
-#include "gromacs/selection/selmethod.h"
#include "gromacs/utility/exceptions.h"
+#include "gromacs/utility/smalloc.h"
#include "selelem.h"
+#include "selmethod.h"
using std::min;
using std::max;
-/*! \internal \brief
+/*! \internal
+ * \brief
* Internal data structure for the \p insolidangle selection method.
*
* \see \c t_partition
+ *
+ * \ingroup module_selection
*/
typedef struct
{
int bin;
} t_partition_item;
-/*! \internal \brief
+/*! \internal
+ * \brief
* Internal data structure for the \p insolidangle selection method.
*
* Describes the surface partitioning within one slice along the zenith angle.
* The slice from azimuthal angle \p p[i].left to \p p[i+1].left belongs to
* bin \p p[i].bin.
+ *
+ * \ingroup module_selection
*/
typedef struct
{
t_partition_item *p;
} t_partition;
-/*! \internal \brief
+/*! \internal
+ * \brief
* Internal data structure for the \p insolidangle selection method.
*
* Contains the reference points that partially cover a certain region on the
* surface of the unit sphere.
* If \p n is -1, the whole region described by the bin is covered.
+ *
+ * \ingroup module_selection
*/
typedef struct
{
rvec *x;
} t_spheresurfacebin;
-/*! \internal \brief
+/*! \internal
+ * \brief
* Data structure for the \p insolidangle selection method.
*
* All angle values are in the units of radians.
+ *
+ * \ingroup module_selection
*/
typedef struct
{
t_spheresurfacebin *bin;
} t_methoddata_insolidangle;
-/** Allocates data for the \p insolidangle selection method. */
+/*! \brief
+ * Allocates data for the \p insolidangle selection method.
+ *
+ * \param[in] npar Not used (should be 3).
+ * \param[in,out] param Method parameters (should point to
+ * \ref smparams_insolidangle).
+ * \returns Pointer to the allocated data (\ref t_methoddata_insolidangle).
+ *
+ * Allocates memory for a \ref t_methoddata_insolidangle structure and
+ * initializes the parameter as follows:
+ * - \p center defines the value for t_methoddata_insolidangle::center.
+ * - \p span defines the value for t_methoddata_insolidangle::span.
+ * - \p cutoff defines the value for t_methoddata_insolidangle::angcut.
+ */
static void *
init_data_insolidangle(int npar, gmx_ana_selparam_t *param);
-/** Initializes the \p insolidangle selection method. */
+/*! \brief
+ * Initializes the \p insolidangle selection method.
+ *
+ * \param top Not used.
+ * \param npar Not used.
+ * \param param Not used.
+ * \param data Pointer to \ref t_methoddata_insolidangle to initialize.
+ * \returns 0 on success, -1 on failure.
+ *
+ * Converts t_methoddata_insolidangle::angcut to radians and allocates
+ * and allocates memory for the bins used during the evaluation.
+ */
static void
-init_insolidangle(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data);
+init_insolidangle(t_topology * top, int npar, gmx_ana_selparam_t * param, void *data);
/** Frees the data allocated for the \p insolidangle selection method. */
static void
free_data_insolidangle(void *data);
-/** Initializes the evaluation of the \p insolidangle selection method for a frame. */
+/*! \brief
+ * Initializes the evaluation of the \p insolidangle selection method for a frame.
+ *
+ * \param[in] top Not used.
+ * \param[in] fr Not used.
+ * \param[in] pbc PBC structure.
+ * \param data Should point to a \ref t_methoddata_insolidangle.
+ *
+ * Creates a lookup structure that enables fast queries of whether a point
+ * is within the solid angle or not.
+ */
static void
-init_frame_insolidangle(t_topology *top, t_trxframe *fr, t_pbc *pbc, void *data);
+init_frame_insolidangle(t_topology * top, t_trxframe * fr, t_pbc *pbc, void *data);
/** Internal helper function for evaluate_insolidangle(). */
static bool
accept_insolidangle(rvec x, t_pbc *pbc, void *data);
/** Evaluates the \p insolidangle selection method. */
static void
-evaluate_insolidangle(t_topology *top, t_trxframe *fr, t_pbc *pbc,
+evaluate_insolidangle(t_topology * /* top */, t_trxframe * /* fr */, t_pbc *pbc,
gmx_ana_pos_t *pos, gmx_ana_selvalue_t *out, void *data);
/** Calculates the distance between unit vectors. */
/** Adds a single reference point and updates the surface bins. */
static void
store_surface_point(t_methoddata_insolidangle *surf, rvec x);
-/** Optimizes the surface bins for faster searching. */
+/*! \brief
+ * Optimizes the surface bins for faster searching.
+ *
+ * \param[in,out] surf Surface data structure.
+ *
+ * Currently, this function does nothing.
+ */
static void
optimize_surface_points(t_methoddata_insolidangle *surf);
/** Estimates the area covered by the reference cones. */
"of these cones. The cutoff determines the width of the cones.",
};
-/** \internal Selection method data for the \p insolidangle method. */
+/** Selection method data for the \p insolidangle method. */
gmx_ana_selmethod_t sm_insolidangle = {
"insolidangle", GROUP_VALUE, SMETH_DYNAMIC,
asize(smparams_insolidangle), smparams_insolidangle,
asize(help_insolidangle), help_insolidangle},
};
-/*!
- * \param[in] npar Not used (should be 3).
- * \param[in,out] param Method parameters (should point to
- * \ref smparams_insolidangle).
- * \returns Pointer to the allocated data (\ref t_methoddata_insolidangle).
- *
- * Allocates memory for a \ref t_methoddata_insolidangle structure and
- * initializes the parameter as follows:
- * - \p center defines the value for t_methoddata_insolidangle::center.
- * - \p span defines the value for t_methoddata_insolidangle::span.
- * - \p cutoff defines the value for t_methoddata_insolidangle::angcut.
- */
static void *
-init_data_insolidangle(int npar, gmx_ana_selparam_t *param)
+init_data_insolidangle(int /* npar */, gmx_ana_selparam_t *param)
{
- t_methoddata_insolidangle *data;
+ t_methoddata_insolidangle *data = new t_methoddata_insolidangle();
+ data->angcut = 5.0;
+ data->cfrac = 0.0;
+
+ data->distccut = 0.0;
+ data->targetbinsize = 0.0;
+
+ data->ntbins = 0;
+ data->tbinsize = 0.0;
+ data->tbin = NULL;
+ data->maxbins = 0;
+ data->nbins = 0;
+ data->bin = NULL;
- snew(data, 1);
- data->angcut = 5.0;
param[0].val.u.p = &data->center;
param[1].val.u.p = &data->span;
param[2].val.u.r = &data->angcut;
return data;
}
-/*!
- * \param top Not used.
- * \param npar Not used.
- * \param param Not used.
- * \param data Pointer to \ref t_methoddata_insolidangle to initialize.
- * \returns 0 on success, -1 on failure.
- *
- * Converts t_methoddata_insolidangle::angcut to radians and allocates
- * and allocates memory for the bins used during the evaluation.
- */
static void
-init_insolidangle(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data)
+init_insolidangle(t_topology * /* top */, int /* npar */, gmx_ana_selparam_t * /* param */, void *data)
{
t_methoddata_insolidangle *surf = (t_methoddata_insolidangle *)data;
int i, c;
surf->angcut *= DEG2RAD;
- surf->distccut = -cos(surf->angcut);
+ surf->distccut = -cos(surf->angcut);
surf->targetbinsize = surf->angcut / 2;
- surf->ntbins = static_cast<int>(M_PI / surf->targetbinsize);
- surf->tbinsize = (180.0 / surf->ntbins)*DEG2RAD;
+ surf->ntbins = static_cast<int>(M_PI / surf->targetbinsize);
+ surf->tbinsize = (180.0 / surf->ntbins)*DEG2RAD;
snew(surf->tbin, static_cast<int>(M_PI / surf->tbinsize) + 1);
surf->maxbins = 0;
}
free_surface_points(d);
sfree(d->bin);
- sfree(d);
+ delete d;
}
-/*!
- * \param[in] top Not used.
- * \param[in] fr Current frame.
- * \param[in] pbc PBC structure.
- * \param data Should point to a \ref t_methoddata_insolidangle.
- *
- * Creates a lookup structure that enables fast queries of whether a point
- * is within the solid angle or not.
- */
static void
-init_frame_insolidangle(t_topology *top, t_trxframe *fr, t_pbc *pbc, void *data)
+init_frame_insolidangle(t_topology * /* top */, t_trxframe * /* fr */, t_pbc *pbc, void *data)
{
t_methoddata_insolidangle *d = (t_methoddata_insolidangle *)data;
rvec dx;
free_surface_points(d);
clear_surface_points(d);
- for (i = 0; i < d->span.nr; ++i)
+ for (i = 0; i < d->span.count(); ++i)
{
if (pbc)
{
* \c t_methoddata_insolidangle::center, and stores the result in \p out->u.g.
*/
static void
-evaluate_insolidangle(t_topology *top, t_trxframe *fr, t_pbc *pbc,
+evaluate_insolidangle(t_topology * /* top */, t_trxframe * /* fr */, t_pbc *pbc,
gmx_ana_pos_t *pos, gmx_ana_selvalue_t *out, void *data)
{
- int b;
-
out->u.g->isize = 0;
- for (b = 0; b < pos->nr; ++b)
+ for (int b = 0; b < pos->count(); ++b)
{
if (accept_insolidangle(pos->x[b], pbc, data))
{
- gmx_ana_pos_append(NULL, out->u.g, pos, b, 0);
+ gmx_ana_pos_add_to_group(out->u.g, pos, b);
}
}
}
{
real theta, phi;
int tbin, pbin;
-
+
theta = acos(x[ZZ]);
- phi = atan2(x[YY], x[XX]);
- tbin = static_cast<int>(floor(theta / surf->tbinsize));
+ phi = atan2(x[YY], x[XX]);
+ tbin = static_cast<int>(floor(theta / surf->tbinsize));
if (tbin >= surf->ntbins)
{
tbin = surf->ntbins - 1;
surf->nbins = 0;
for (i = 0; i < surf->ntbins; ++i)
{
- c = static_cast<int>(min(sin(surf->tbinsize*i),
+ c = static_cast<int>(min(sin(surf->tbinsize*i),
sin(surf->tbinsize*(i+1)))
* M_2PI / surf->targetbinsize) + 1;
if (c <= 0)
surf->tbin[i].n = c;
for (j = 0; j < c; ++j)
{
- surf->tbin[i].p[j].left = -M_PI + j*M_2PI/c - 0.0001;
- surf->tbin[i].p[j].bin = surf->nbins;
+ surf->tbin[i].p[j].left = -M_PI + j*M_2PI/c - 0.0001;
+ surf->tbin[i].p[j].bin = surf->nbins;
surf->bin[surf->nbins].n = 0;
surf->nbins++;
}
surf->tbin[i].p[c].left = M_PI + 0.0001;
- surf->tbin[i].p[c].bin = -1;
+ surf->tbin[i].p[c].bin = -1;
}
}
sfree(surf->bin[i].x);
}
surf->bin[i].n_alloc = 0;
- surf->bin[i].x = NULL;
+ surf->bin[i].x = NULL;
}
}
bin = surf->tbin[tbin].p[pbin].bin;
/* Return if bin is already completely covered */
if (surf->bin[bin].n == -1)
+ {
return;
+ }
/* Allocate more space if necessary */
- if (surf->bin[bin].n == surf->bin[bin].n_alloc) {
+ if (surf->bin[bin].n == surf->bin[bin].n_alloc)
+ {
surf->bin[bin].n_alloc += 10;
srenew(surf->bin[bin].x, surf->bin[bin].n_alloc);
}
{
int bin;
- bin = surf->tbin[tbin].p[pbin].bin;
+ bin = surf->tbin[tbin].p[pbin].bin;
surf->bin[bin].n = -1;
}
/* Find the edges of the bins affected */
pdelta = max(max(pdelta1, pdelta2), pdeltamax);
- phi1 = phi - pdelta;
+ phi1 = phi - pdelta;
if (phi1 >= -M_PI)
{
- pbin = find_partition_bin(&surf->tbin[tbin], phi1);
+ pbin = find_partition_bin(&surf->tbin[tbin], phi1);
pbin1 = pbin;
}
else
{
- pbin = find_partition_bin(&surf->tbin[tbin], phi1 + M_2PI);
+ pbin = find_partition_bin(&surf->tbin[tbin], phi1 + M_2PI);
pbin1 = pbin - surf->tbin[tbin].n;
}
phi2 = phi + pdelta;
}
else
{
- pbin2 = find_partition_bin(&surf->tbin[tbin], phi2 - M_2PI);
+ pbin2 = find_partition_bin(&surf->tbin[tbin], phi2 - M_2PI);
pbin2 += surf->tbin[tbin].n;
}
++pbin2;
}
/* Find the edges of completely covered region */
pdelta = min(pdelta1, pdelta2);
- phi1 = phi - pdelta;
+ phi1 = phi - pdelta;
if (phi1 < -M_PI)
{
phi1 += M_2PI;
/* Wrap bin around if end reached */
if (pbin == surf->tbin[tbin].n)
{
- pbin = 0;
+ pbin = 0;
phi1 -= M_2PI;
phi2 -= M_2PI;
}
int tbin;
theta = acos(x[ZZ]);
- phi = atan2(x[YY], x[XX]);
+ phi = atan2(x[YY], x[XX]);
/* Find the maximum extent in the phi direction */
if (theta <= surf->angcut)
{
pdeltamax = M_PI;
- tmax = 0;
+ tmax = 0;
}
else if (theta >= M_PI - surf->angcut)
{
pdeltamax = M_PI;
- tmax = M_PI;
+ tmax = M_PI;
}
else
{
pdeltamax = asin(sin(surf->angcut) / sin(theta));
- tmax = acos(cos(theta) / cos(surf->angcut));
+ tmax = acos(cos(theta) / cos(surf->angcut));
}
/* Find the first affected bin */
- tbin = max(static_cast<int>(floor((theta - surf->angcut) / surf->tbinsize)), 0);
+ tbin = max(static_cast<int>(floor((theta - surf->angcut) / surf->tbinsize)), 0);
theta1 = tbin * surf->tbinsize;
if (theta1 < theta - surf->angcut)
{
* such that the case above catches this instead of falling through
* here. */
pdelta2 = 2*asin(sqrt(
- (sqr(sin(surf->angcut/2)) - sqr(sin((theta2-theta)/2))) /
- (sin(theta) * sin(theta2))));
+ (sqr(sin(surf->angcut/2)) - sqr(sin((theta2-theta)/2))) /
+ (sin(theta) * sin(theta2))));
}
/* Update the bin */
if (tmax >= theta1 && tmax <= theta2)
update_surface_bin(surf, tbin, phi, pdelta1, pdelta2, 0, x);
}
/* Next bin */
- theta1 = theta2;
+ theta1 = theta2;
pdelta1 = pdelta2;
++tbin;
}
}
-/*!
- * \param[in,out] surf Surface data structure.
- *
- * Currently, this function does nothing.
- */
static void
-optimize_surface_points(t_methoddata_insolidangle *surf)
+optimize_surface_points(t_methoddata_insolidangle * /* surf */)
{
/* TODO: Implement */
}
for (p = 0; p < surf->tbin[t].n; ++p)
{
pfrac = surf->tbin[t].p[p+1].left - surf->tbin[t].p[p].left;
- n = surf->bin[surf->tbin[t].p[p].bin].n;
+ n = surf->bin[surf->tbin[t].p[p].bin].n;
if (n == -1) /* Bin completely covered */
{
cfrac += tfrac * pfrac;
}
- else if (n > 0) /* Bin partially covered */
+ else if (n > 0) /* Bin partially covered */
{
cfrac += tfrac * pfrac / 2; /* A rough estimate */
}
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff