src/fitng.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2013,2014,2015, by the GROMACS development team, led by
   5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   6  * and including many others, as listed in the AUTHORS file in the
   7  * top-level source directory and at http://www.gromacs.org.
   8  *
   9  * GROMACS is free software; you can redistribute it and/or
  10  * modify it under the terms of the GNU Lesser General Public License
  11  * as published by the Free Software Foundation; either version 2.1
  12  * of the License, or (at your option) any later version.
  13  *
  14  * GROMACS is distributed in the hope that it will be useful,
  15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17  * Lesser General Public License for more details.
  18  *
  19  * You should have received a copy of the GNU Lesser General Public
  20  * License along with GROMACS; if not, see
  21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  23  *
  24  * If you want to redistribute modifications to GROMACS, please
  25  * consider that scientific software is very special. Version
  26  * control is crucial - bugs must be traceable. We will be happy to
  27  * consider code for inclusion in the official distribution, but
  28  * derived work must not be called official GROMACS. Details are found
  29  * in the README & COPYING files - if they are missing, get the
  30  * official version at http://www.gromacs.org.
  31  *
  32  * To help us fund GROMACS development, we humbly ask that you cite
  33  * the research papers on the package. Check out http://www.gromacs.org.
  34  */
  35 /*! \internal \file
  36  * \brief
  37  * Implements gmx::analysismodules::Freevolume.
  38  *
  39  * \author Titov Anatoly <Wapuk-cobaka@yandex.ru>
  40  * \ingroup module_trajectoryanalysis
  41  */
  42
  43 #include <iostream>
  44 #include <chrono>
  45 #include <string.h>
  46 #include <algorithm>
  47 #include <cstdio>
  48
  49 #include <omp.h>
  50
  51 #include <gromacs/fileio/trxio.h>
  52 #include <gromacs/trajectoryanalysis/topologyinformation.h>
  53
  54 #include <newfit.h>
  55
  56 class Fitng : public TrajectoryAnalysisModule
  57 {
  58     public:
  59
  60         Fitng();
  61         virtual ~Fitng();
  62
  63         //! Set the options and setting
  64         virtual void initOptions(IOptionsContainer          *options,
  65                                  TrajectoryAnalysisSettings *settings);
  66
  67         //! First routine called by the analysis framework
  68         // virtual void initAnalysis(const t_trxframe &fr, t_pbc *pbc);
  69         virtual void initAnalysis(const TrajectoryAnalysisSettings &settings,
  70                                   const TopologyInformation        &top);
  71
  72         //! Call for each frame of the trajectory
  73         // virtual void analyzeFrame(const t_trxframe &fr, t_pbc *pbc);
  74         virtual void analyzeFrame(int frnr, const t_trxframe &fr, t_pbc *pbc,
  75                                   TrajectoryAnalysisModuleData *pdata);
  76
  77         //! Last routine called by the analysis framework
  78         // virtual void finishAnalysis(t_pbc *pbc);
  79         virtual void finishAnalysis(int nframes);
  80
  81         //! Routine to write output, that is additional over the built-in
  82         virtual void writeOutput();
  83
  84     private:
  85         Selection                                   sel_;
  86         std::vector< int >                          index;
  87         std::vector < RVec >                        trajectoryFrame;
  88         std::vector< RVec >                         reference;
  89         std::vector< std::pair< size_t, size_t > >  fitPairs;
  90         double                                      fitConst {0.000'001};
  91         std::string                                 outputTrjName;
  92         t_trxframe                                  tempFrame;
  93         t_trxstatus                                 *op;
  94 };
  95
  96 Fitng::Fitng(): TrajectoryAnalysisModule()
  97 {
  98 }
  99
 100 Fitng::~Fitng()
 101 {
 102 }
 103
 104 void
 105 Fitng::initOptions( IOptionsContainer          *options,
 106                     TrajectoryAnalysisSettings *settings)
 107 {
 108     static const char *const desc[] = {
 109         "[THISMODULE] to be done"
 110     };
 111     // Add the descriptive text (program help text) to the options
 112     settings->setHelpText(desc);
 113     // Add option for selection list
 114     options->addOption(SelectionOption("select")
 115                         .store(&sel_).required()
 116                         .description("Atoms that are considered as part of the excluded volume"));
 117     // Add option for Fit constant
 118     options->addOption(DoubleOption("fitConst")
 119                         .store(&fitConst)
 120                         .description("Fitting untill diff <= FitConst, by default == 0.000'001"));
 121     // Add option for output pdb traj (meh edition)
 122     options->addOption(StringOption("out")
 123                         .store(&outputTrjName)
 124                         .description("transformed trajectory"));
 125     // Control input settings
 126     settings->setFlags(TrajectoryAnalysisSettings::efNoUserPBC);
 127     settings->setFlag(TrajectoryAnalysisSettings::efUseTopX);
 128     settings->setPBC(true);
 129 }
 130
 131 void
 132 Fitng::initAnalysis(const TrajectoryAnalysisSettings    &settings,
 133                     const TopologyInformation           &top)
 134 {
 135     index.resize(0);
 136     for (ArrayRef< const int >::iterator ai {sel_.atomIndices().begin()}; ai < sel_.atomIndices().end(); ++ai) {
 137         index.push_back(*ai);
 138     }
 139
 140     reference.resize(0);
 141     if (top.hasFullTopology()) {
 142         for (size_t i {0}; i < index.size(); ++i) {
 143             reference.push_back(top.x().at(index[i]));
 144         }
 145     }
 146
 147     op = open_trx(outputTrjName.c_str(), "w+");
 148     trajectoryFrame.reserve(index.size());
 149
 150     fitPairs.resize(index.size());
 151     for (size_t i {0}; i < index.size(); ++i) {
 152         fitPairs[i] = std::make_pair(i, i);
 153     }
 154 }
 155
 156 void
 157 Fitng::analyzeFrame(int                             frnr,
 158                     const t_trxframe               &fr,
 159                     t_pbc                          *pbc,
 160                     TrajectoryAnalysisModuleData   *pdata)
 161 {
 162     trajectoryFrame.resize(0);
 163     for (size_t i {0}; i < index.size(); ++i) {
 164         trajectoryFrame.push_back(fr.x[index[i]]);
 165     }
 166
 167     MyFitNew(reference, trajectoryFrame, fitPairs, fitConst);
 168
 169     tempFrame = fr;
 170     for (size_t i {0}; i < index.size(); ++i) {
 171         copy_rvec(trajectoryFrame[i].as_vec(), tempFrame.x[index[i]]);
 172     }
 173
 174     if (write_trxframe_indexed(op, static_cast<const t_trxframe *>(&tempFrame), index.size(), static_cast<const int *>(&index.front()), NULL) != 0) {
 175         std::cout << "\nFileOutPut error | frame number = " << frnr << "\n";
 176     }
 177 }
 178
 179 void
 180 Fitng::finishAnalysis(int nframes)
 181 {
 182     std::cout << "\nGAME OVER\n";
 183 }
 184
 185 void
 186 Fitng::writeOutput()
 187 {
 188
 189 }
 190
 191
 192 /*! \brief
 193  * The main function for the analysis template.
 194  */
 195 int
 196 main(int argc, char *argv[])
 197 {
 198     return gmx::TrajectoryAnalysisCommandLineRunner::runAsMain<Fitng>(argc, argv);
 199 }