src/gromacs/trajectoryanalysis/modules/distance.cpp

   1 /*
   2  *
   3  *                This source code is part of
   4  *
   5  *                 G   R   O   M   A   C   S
   6  *
   7  *          GROningen MAchine for Chemical Simulations
   8  *
   9  * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
  10  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  11  * Copyright (c) 2001-2009, The GROMACS development team,
  12  * check out http://www.gromacs.org for more information.
  13
  14  * This program is free software; you can redistribute it and/or
  15  * modify it under the terms of the GNU General Public License
  16  * as published by the Free Software Foundation; either version 2
  17  * of the License, or (at your option) any later version.
  18  *
  19  * If you want to redistribute modifications, please consider that
  20  * scientific software is very special. Version control is crucial -
  21  * bugs must be traceable. We will be happy to consider code for
  22  * inclusion in the official distribution, but derived work must not
  23  * be called official GROMACS. Details are found in the README & COPYING
  24  * files - if they are missing, get the official version at www.gromacs.org.
  25  *
  26  * To help us fund GROMACS development, we humbly ask that you cite
  27  * the papers on the package - you can find them in the top README file.
  28  *
  29  * For more info, check our website at http://www.gromacs.org
  30  */
  31 /*! \internal \file
  32  * \brief
  33  * Implements gmx::analysismodules::Distance.
  34  *
  35  * \author Teemu Murtola <teemu.murtola@cbr.su.se>
  36  * \ingroup module_trajectoryanalysis
  37  */
  38 #include "distance.h"
  39
  40 #include "gromacs/legacyheaders/pbc.h"
  41 #include "gromacs/legacyheaders/vec.h"
  42
  43 #include "gromacs/analysisdata/analysisdata.h"
  44 #include "gromacs/analysisdata/modules/plot.h"
  45 #include "gromacs/options/basicoptions.h"
  46 #include "gromacs/options/filenameoption.h"
  47 #include "gromacs/options/options.h"
  48 #include "gromacs/selection/selection.h"
  49 #include "gromacs/selection/selectionoption.h"
  50 #include "gromacs/trajectoryanalysis/analysissettings.h"
  51 #include "gromacs/utility/exceptions.h"
  52 #include "gromacs/utility/stringutil.h"
  53
  54 namespace gmx
  55 {
  56
  57 namespace analysismodules
  58 {
  59
  60 const char Distance::name[] = "distance";
  61 const char Distance::shortDescription[] =
  62     "Calculate distances";
  63
  64 Distance::Distance()
  65     : TrajectoryAnalysisModule(name, shortDescription),
  66       avem_(new AnalysisDataAverageModule())
  67 {
  68     data_.setColumnCount(4);
  69     registerAnalysisDataset(&data_, "distance");
  70 }
  71
  72
  73 Distance::~Distance()
  74 {
  75 }
  76
  77
  78 void
  79 Distance::initOptions(Options *options, TrajectoryAnalysisSettings * /*settings*/)
  80 {
  81     static const char *const desc[] = {
  82         "g_dist can calculate the distance between two positions as",
  83         "a function of time. The total distance and its",
  84         "x, y and z components are plotted."
  85     };
  86
  87     options->setDescription(concatenateStrings(desc));
  88
  89     options->addOption(FileNameOption("o").filetype(eftPlot).outputFile()
  90                            .store(&fnDist_).defaultBasename("dist")
  91                            .description("Computed distances"));
  92     options->addOption(SelectionOption("select").required().valueCount(2)
  93                            .store(sel_));
  94 }
  95
  96
  97 void
  98 Distance::initAnalysis(const TrajectoryAnalysisSettings &settings,
  99                        const TopologyInformation & /*top*/)
 100 {
 101     if (sel_[0].posCount() != 1)
 102     {
 103         GMX_THROW(InvalidInputError("The first selection does not define a single position"));
 104     }
 105     if (sel_[1].posCount() != 1)
 106     {
 107         GMX_THROW(InvalidInputError("The second selection does not define a single position"));
 108     }
 109
 110     data_.addModule(avem_);
 111     AnalysisDataPlotModulePointer plotm_(new AnalysisDataPlotModule());
 112     plotm_->setSettings(settings.plotSettings());
 113     plotm_->setFileName(fnDist_);
 114     plotm_->setTitle("Distance");
 115     plotm_->setXAxisIsTime();
 116     plotm_->setYLabel("Distance (nm)");
 117     data_.addModule(plotm_);
 118 }
 119
 120
 121 void
 122 Distance::analyzeFrame(int frnr, const t_trxframe &fr, t_pbc *pbc,
 123                        TrajectoryAnalysisModuleData *pdata)
 124 {
 125     AnalysisDataHandle  dh = pdata->dataHandle(data_);
 126     const Selection    &sel1 = pdata->parallelSelection(sel_[0]);
 127     const Selection    &sel2 = pdata->parallelSelection(sel_[1]);
 128     rvec                dx;
 129     real                r;
 130     const SelectionPosition &p1 = sel1.position(0);
 131     const SelectionPosition &p2 = sel2.position(0);
 132
 133     if (pbc != NULL)
 134     {
 135         pbc_dx(pbc, p1.x(), p2.x(), dx);
 136     }
 137     else
 138     {
 139         rvec_sub(p1.x(), p2.x(), dx);
 140     }
 141     r = norm(dx);
 142     dh.startFrame(frnr, fr.time);
 143     dh.setPoint(0, r);
 144     dh.setPoints(1, 3, dx);
 145     dh.finishFrame();
 146 }
 147
 148
 149 void
 150 Distance::finishAnalysis(int /*nframes*/)
 151 {
 152 }
 153
 154
 155 void
 156 Distance::writeOutput()
 157 {
 158     fprintf(stderr, "Average distance: %f\n", avem_->average(0));
 159     fprintf(stderr, "Std. deviation:   %f\n", avem_->stddev(0));
 160 }
 161
 162 } // namespace analysismodules
 163
 164 } // namespace gmx