2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2021, 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.
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.
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.
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.
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.
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.
37 * Declares trajectory analysis module for secondary structure asignment.
39 * \author Sergey Gorelov <gorelov_sv@pnpi.nrcki.ru>
40 * \author Anatoly Titov <titov_ai@pnpi.nrcki.ru>
41 * \author Alexey Shvetsov <alexxyum@gmail.com>
42 * \ingroup module_trajectoryanalysis
45 #ifndef GMX_TRAJECTORYANALYSIS_MODULES_DSSPTOOLS_H
46 #define GMX_TRAJECTORYANALYSIS_MODULES_DSSPTOOLS_H
49 #include "gromacs/math/units.h"
50 #include "gromacs/pbcutil/pbc.h"
51 #include <gromacs/trajectoryanalysis.h>
52 #include "gromacs/trajectoryanalysis/topologyinformation.h"
60 namespace analysismodules
63 enum class NBSearchMethod : std::size_t
71 const gmx::EnumerationArray<NBSearchMethod, const char*> NBSearchMethodNames = {
72 { "Classique", "Experimental", "DSSP" }
75 struct initParameters {
77 real cutoff_; // = 4.0; ???
79 bool verbose, PPHelices, addHydrogens;
82 enum class backboneAtomTypes : std::size_t
91 //! String values corresponding to backbone atom types
92 const gmx::EnumerationArray<backboneAtomTypes, const char*> backboneAtomTypeNames = {
93 { "CA", "C", "O", "N", "H" }
97 std::array<std::size_t, static_cast<std::size_t>(backboneAtomTypes::Count)> _backboneIndices{ 0, 0, 0, 0, 0 }; // TODO something with zeroes
98 std::size_t getIndex(backboneAtomTypes atomTypeName) const;
99 t_resinfo *info{nullptr}, *donor[2]{nullptr, nullptr}, *acceptor[2]{nullptr, nullptr};
100 ResInfo *prevResi{nullptr}, *nextResi{nullptr};
101 float donorEnergy[2]{}, acceptorEnergy[2]{};
102 bool is_proline{false};
109 // void setIndex(backboneAtomTypes atomTypeName, std::size_t atomIndex);
110 // std::size_t getIndex(backboneAtomTypes atomTypeName) const;
112 // ResInfo *donor[2]{nullptr, nullptr}, *acceptor[2]{nullptr, nullptr};
113 // float donorEnergy[2]{}, acceptorEnergy[2]{};
114 // std::string resiName;
117 // std::array<std::size_t, static_cast<std::size_t>(backboneAtomTypes::Count)> _ResInfo{ 0, 0, 0, 0, 0 };
120 enum class secondaryStructureTypes : std::size_t { // TODO
135 enum class turnsTypes : std::size_t {
142 enum class HelixPositions : std::size_t {
151 enum class bridgeTypes : std::size_t {
158 class secondaryStructures{ // PatterSearch Wrapper
160 secondaryStructures();
161 void initiateSearch(const std::vector<ResInfo> &ResInfoMatrix, const bool PiHelicesPreferencez = true);
162 std::string patternSearch();
163 ~secondaryStructures();
165 class secondaryStructuresData{ // PatternSearch Tool
167 void setStatus(const secondaryStructureTypes secondaryStructureTypeName);
168 void setStatus(const HelixPositions helixPosition, const turnsTypes turn);
169 bool getStatus(const secondaryStructureTypes secondaryStructureTypeName) const, isBreakPartnerWith(const secondaryStructuresData *partner) const;
170 HelixPositions getStatus(const turnsTypes turn) const;
171 secondaryStructureTypes getStatus() const;
172 void setBreak(secondaryStructuresData *breakPartner), setBridge(secondaryStructuresData *bridgePartner, std::size_t bridgePartnerIndex, bridgeTypes bridgeType);
173 bool hasBridges() const, hasBridges(bridgeTypes bridgeType) const, isBridgePartnerWith(secondaryStructuresData *bridgePartner, bridgeTypes bridgeType) const;
174 std::size_t getBridgePartnerIndex(bridgeTypes bridgeType) const;
175 secondaryStructuresData getBridgePartner(bridgeTypes bridgeType) const;
177 std::array<bool, static_cast<std::size_t>(secondaryStructureTypes::Count)> SecondaryStructuresStatusArray{ true, 0, 0, 0, 0, 0, 0 };
178 secondaryStructuresData *breakPartners[2]{nullptr, nullptr};
179 secondaryStructuresData *bridgePartners[2]{nullptr, nullptr};
180 std::size_t bridgePartnersIndexes[2]{};
181 secondaryStructureTypes SecondaryStructuresStatus {secondaryStructureTypes::Loop};
182 std::array<HelixPositions, static_cast<std::size_t>(turnsTypes::Count)> TurnsStatusArray {HelixPositions::None, HelixPositions::None, HelixPositions::None};
185 std::vector<secondaryStructuresData> SecondaryStructuresStatusMap;
188 const std::vector<ResInfo> *ResInfoMap;
190 const gmx::EnumerationArray<secondaryStructureTypes, const char> secondaryStructureTypeNames = {
191 { '~', '=', 'S', 'T', 'P', 'I', 'G', 'E', 'B', 'H'}
193 std::string SecondaryStructuresStringLine;
195 bool hasHBondBetween(std::size_t resi1, std::size_t resi2) const;
197 bool NoChainBreaksBetween(std::size_t ResiStart, std::size_t ResiEnd) const, isLoop(const std::size_t resiIndex) const, PiHelixPreference;
198 bridgeTypes calculateBridge(std::size_t i, std::size_t j) const;
199 void analyzeBridgesAndStrandsPatterns(), analyzeTurnsAndHelicesPatterns(), analyzePPHelicesPatterns();
201 const float HBondEnergyCutOff{ -0.5 };
204 class alternateNeighborhoodSearch{
206 alternateNeighborhoodSearch(){}
207 void setCutoff(const real &cutoff_init);
208 void AltPairSearch(const t_trxframe &fr, const std::vector<ResInfo> &IndexMap);
210 std::size_t getResiI() const, getResiJ() const;
213 std::size_t resiIpos{0}, resiJpos{0};
214 std::vector<std::vector<bool>> PairMap;
215 std::vector<std::size_t> MiniBoxMap, partners, num_of_miniboxes;
216 std::vector<std::vector<std::size_t>> MiniBoxesReverseMap;
217 std::vector<std::vector<std::vector<std::vector<std::size_t>>>> MiniBoxesMap;
218 std::array<int, 3> MiniBoxSize{0, 0, 0}, offset{-1, -1, -1}, fixBox{0,0,0};
219 void GetMiniBoxesMap(const t_trxframe &fr, const std::vector<ResInfo> &IndexMap);
220 void FixAtomCoordinates(real &coordinate, const real vector_length);
221 void ReCalculatePBC(int &x, const int &x_max);
222 std::vector<std::vector<bool>>::iterator ResiI;
223 std::vector<bool>::iterator ResiJ;
231 void initAnalysis(/*const TrajectoryAnalysisSettings &settings,*/const TopologyInformation& top, const initParameters &initParamz);
233 void analyzeFrame(int frnr, const t_trxframe& fr, t_pbc* pbc);
237 std::vector<std::pair<int, std::string>> getData();
239 initParameters initParams;
241 std::string filename_;
242 void calculateBends(const t_trxframe& fr, const t_pbc* pbc);
243 std::vector<std::string> ResiNames;
244 std::vector<std::size_t> AtomResi;
245 std::vector<ResInfo> IndexMap;
246 std::vector<std::vector<std::size_t>> nturnsmap, Bridges, AntiBridges;
247 secondaryStructures PatternSearch;
248 const float HBondEnergyCutOff{ -0.5 },
249 minimalCAdistance{ 9.0 };
250 void calculateHBondEnergy(ResInfo& Donor,
252 const t_trxframe& fr,
254 float CalculateAtomicDistances(const int &A, const int &B, const t_trxframe &fr, const t_pbc *pbc);
255 float CalculateAtomicDistances(const rvec &A, const int &B, const t_trxframe &fr, const t_pbc *pbc);
262 /*Storages our dirty data, duh*/
263 void storageData(int frnr, std::string data);
265 /*Clear data after usage for the next pdb file*/
268 /*Perform several checks and fucks and shits and returns sorted sexy data*/
269 std::vector<std::pair<int, std::string>> returnData();
271 std::vector<std::pair<int, std::string>> storaged_data;
272 static std::mutex mx;
275 static DsspStorage Storage;
279 //end of namespace analysismodules
282 //end of namespace gmx