FUCK

diff --git a/src/dssp.cpp b/src/dssp.cpp
index 91c562e4d225c25480d55eefc4c0548dfe4506ab..51074df087404d3a7430f64c30921f05a725a5db 100644 (file)
--- a/src/dssp.cpp
+++ b/src/dssp.cpp
@@ -112,13 +112,12 @@ void Dssp::optionsFinished(TrajectoryAnalysisSettings* /* settings */) {}
 
 void Dssp::initAnalysis(const TrajectoryAnalysisSettings &settings, const TopologyInformation& top)
 {
-    DT.setInitValues(initParams);
-    DT.initAnalysis(top);
+    DT.initAnalysis(top, initParams);
 }
 
 void Dssp::analyzeFrame(int frnr, const t_trxframe& fr, t_pbc* pbc, TrajectoryAnalysisModuleData* /* pdata*/)
 {
-    DT.analyzeFrame(frnr, fr, pbc, alternateSearch);
+    DT.analyzeFrame(frnr, fr, pbc);
 }
 
 void Dssp::finishAnalysis(int /*nframes*/) {

diff --git a/src/dssptools.cpp b/src/dssptools.cpp
index d59aea9fbc5ce7946941acb66b0c0669ee13e75a..f5762465ddfa9dab695a3d480afe24ee25864007 100644 (file)
--- a/src/dssptools.cpp
+++ b/src/dssptools.cpp
@@ -165,11 +165,11 @@ bridgeTypes secondaryStructures::calculateBridge(std::size_t i, std::size_t j) c
     }
 
     if(NoChainBreaksBetween(i - 1, i + 1) && NoChainBreaksBetween(j - 1, j + 1)){
-        if((hasHBondBetween(i + 1, j) && hasHBondBetween(j, i - 1)) || (hasHBondBetween(j + 1, i) && hasHBondBetween(i, j - 1)) ){
-            return bridgeTypes::Bridge;
+        if((hasHBondBetween(i + 1, j) && hasHBondBetween(j, i - 1)) || (hasHBondBetween(j + 1, i) && hasHBondBetween(i, j - 1)) ){ //possibly swap
+            return bridgeTypes::ParallelBridge;
         }
-        else if((hasHBondBetween(i + 1, j - 1) && hasHBondBetween(j + 1, i - 1)) || (hasHBondBetween(j, i) && hasHBondBetween(i, j)) ){
-            return bridgeTypes::AntiBridge;
+        else if((hasHBondBetween(i + 1, j - 1) && hasHBondBetween(j + 1, i - 1)) || (hasHBondBetween(j, i) && hasHBondBetween(i, j)) ){ //possibly swap
+            return bridgeTypes::AntiParallelBridge;
         }
     }
 
@@ -181,7 +181,11 @@ void secondaryStructures::analyzeBridgesAndLaddersPatterns(){
     for(std::size_t i {1}; i + 4 < SecondaryStructuresStatusMap.size(); ++i){
 
         for(std::size_t j {i + 3}; j + 1 < SecondaryStructuresStatusMap.size(); ++j ){
-            bridgeTypes type = calculateBridge(i, j);
+            bridgeTypes type {calculateBridge(i, j)};
+            if (type == bridgeTypes::None){
+                continue;
+            }
+            bool found {false};
         }
     }
 
@@ -295,6 +299,8 @@ void secondaryStructures::analyzeTurnsAndHelicesPatterns(){
                     }
                     Helix = secondaryStructureTypes::Helix_4;
                     break;
+                default:
+                    break;
                 }
                 if ( empty ){
                     for(std::size_t k {0}; k < (static_cast<std::size_t>(i) + strideFactor); ++k ){
@@ -325,7 +331,7 @@ void secondaryStructures::analyzeTurnsAndHelicesPatterns(){
 
 }
 
-std::string secondaryStructures::patternSearch(const std::vector<std::vector<bool>> &HBondsMap){
+std::string secondaryStructures::patternSearch(){
 
     analyzeBridgesAndLaddersPatterns();
     analyzeTurnsAndHelicesPatterns();
@@ -502,14 +508,6 @@ DsspTool::DsspStorage DsspTool::Storage;
 DsspTool::DsspTool(){
 }
 
-void DsspTool::setInitValues(const initParameters &initParams){
-    cutoff_ = initParams.cutoff_;
-    sel_ = initParams.sel_;
-    modeType_ = initParams.mode;
-    NBSM_ = initParams.NBS;
-    verbose_ = initParams.verbose;
-}
-
 void DsspTool::calculateBends(const t_trxframe &fr, const t_pbc *pbc)
 {
    const float benddegree{ 70.0 }, maxdist{ 2.5 };
@@ -595,7 +593,7 @@ bool DsspTool::isHbondExist(const ResInfo& resA,
    distanceCN = CalculateAtomicDistances(
            resA.getIndex(backboneAtomTypes::AtomC), resB.getIndex(backboneAtomTypes::AtomN), fr, pbc);
 
-   switch(modeType_){
+   switch(initParams.mode){
    case modeType::Classique:
        distanceOH = distanceON;
        distanceCH = distanceCN;
@@ -613,7 +611,7 @@ bool DsspTool::isHbondExist(const ResInfo& resA,
                     resA.getIndex(backboneAtomTypes::AtomCA), resB.getIndex(backboneAtomTypes::AtomCA), fr, pbc) << std::endl;
 
    if (resA.getIndex(backboneAtomTypes::AtomC) && resA.getIndex(backboneAtomTypes::AtomO)
-       && resB.getIndex(backboneAtomTypes::AtomN) && ( resB.getIndex(backboneAtomTypes::AtomH) || (modeType_ == modeType::Classique) ) )
+       && resB.getIndex(backboneAtomTypes::AtomN) && ( resB.getIndex(backboneAtomTypes::AtomH) || (initParams.mode == modeType::Classique) ) )
    {
        if (CalculateAtomicDistances(
                    resA.getIndex(backboneAtomTypes::AtomCA), resB.getIndex(backboneAtomTypes::AtomCA), fr, pbc)
@@ -644,11 +642,13 @@ float DsspTool::CalculateAtomicDistances(const int &A, const int &B, const t_trx
    return r.norm(); // TODO * gmx::c_nm2A; if not PDB, i guess????
 }
 
-void DsspTool::initAnalysis(/*const TrajectoryAnalysisSettings &settings,*/const TopologyInformation& top)
+void DsspTool::initAnalysis(/*const TrajectoryAnalysisSettings &settings,*/const TopologyInformation& top, const initParameters &initParamz)
 {
+   initParams = initParamz;
+
    ResInfo _backboneAtoms;
    std::size_t                 i{ 0 };
-   int resicompare{ top.atoms()->atom[static_cast<std::size_t>(*(sel_.atomIndices().begin()))].resind };
+   int resicompare{ top.atoms()->atom[static_cast<std::size_t>(*(initParams.sel_.atomIndices().begin()))].resind };
    IndexMap.resize(0);
    IndexMap.push_back(_backboneAtoms);
    IndexMap[i].setInfo(top.atoms()->resinfo[i]);
@@ -658,7 +658,7 @@ void DsspTool::initAnalysis(/*const TrajectoryAnalysisSettings &settings,*/const
 //       IndexMap[i].setIndex(backboneAtomTypes::PRO, 1);
 //   }
 
-   for (gmx::ArrayRef<const int>::iterator ai{ sel_.atomIndices().begin() }; (ai != sel_.atomIndices().end()); ++ai){
+   for (gmx::ArrayRef<const int>::iterator ai{ initParams.sel_.atomIndices().begin() }; (ai != initParams.sel_.atomIndices().end()); ++ai){
        if (resicompare != top.atoms()->atom[static_cast<std::size_t>(*ai)].resind)
        {
            ++i;
@@ -725,50 +725,66 @@ void DsspTool::initAnalysis(/*const TrajectoryAnalysisSettings &settings,*/const
    }
 }
 
-void DsspTool::analyzeFrame(int frnr, const t_trxframe &fr, t_pbc *pbc, const bool &alternateSearch)
+void DsspTool::analyzeFrame(int frnr, const t_trxframe &fr, t_pbc *pbc)
 {
 
-   if(alternateSearch){
-       alternateNeighborhoodSearch as_;
+    switch(initParams.NBS){
+    case (NBSearchMethod::Classique): {
 
-       as_.setCutoff(cutoff_);
+        // store positions of CA atoms to use them for nbSearch
+        std::vector<gmx::RVec> positionsCA_;
+        for (std::size_t i{ 0 }; i < IndexMap.size(); ++i)
+        {
+            positionsCA_.emplace_back(fr.x[IndexMap[i].getIndex(backboneAtomTypes::AtomCA)]);
+        }
 
-       as_.AltPairSearch(fr, IndexMap);
+        AnalysisNeighborhood nb_;
+        nb_.setCutoff(initParams.cutoff_);
+        AnalysisNeighborhoodPositions       nbPos_(positionsCA_);
+        gmx::AnalysisNeighborhoodSearch     start      = nb_.initSearch(pbc, nbPos_);
+        gmx::AnalysisNeighborhoodPairSearch pairSearch = start.startPairSearch(nbPos_);
+        gmx::AnalysisNeighborhoodPair       pair;
+        while (pairSearch.findNextPair(&pair))
+        {
+            if(isHbondExist(IndexMap[pair.refIndex()], IndexMap[pair.testIndex()], fr, pbc)){
+                IndexMap[pair.refIndex()].setHBondedPair(*(IndexMap[pair.testIndex()].getInfo()));
+            }
+        }
 
-       while (as_.findNextPair()){
-           if(isHbondExist(IndexMap[as_.getResiI()], IndexMap[as_.getResiJ()], fr, pbc)){
-               IndexMap[as_.getResiI()].setHBondedPair(*(IndexMap[as_.getResiJ()].getInfo()));
-           }
-       }
-   }
-   else{
+        break;
+    }
+    case (NBSearchMethod::Experimental): {
 
-       // store positions of CA atoms to use them for nbSearch
-       std::vector<gmx::RVec> positionsCA_;
-       for (std::size_t i{ 0 }; i < IndexMap.size(); ++i)
-       {
-           positionsCA_.emplace_back(fr.x[IndexMap[i].getIndex(backboneAtomTypes::AtomCA)]);
-       }
+        alternateNeighborhoodSearch as_;
 
-       AnalysisNeighborhood nb_;
-       nb_.setCutoff(cutoff_);
-       AnalysisNeighborhoodPositions       nbPos_(positionsCA_);
-       gmx::AnalysisNeighborhoodSearch     start      = nb_.initSearch(pbc, nbPos_);
-       gmx::AnalysisNeighborhoodPairSearch pairSearch = start.startPairSearch(nbPos_);
-       gmx::AnalysisNeighborhoodPair       pair;
-       while (pairSearch.findNextPair(&pair))
-       {
-           if(isHbondExist(IndexMap[pair.refIndex()], IndexMap[pair.testIndex()], fr, pbc)){
-               IndexMap[pair.refIndex()].setHBondedPair(*(IndexMap[pair.testIndex()].getInfo()));
-           }
-       }
+        as_.setCutoff(initParams.cutoff_);
 
-   }
+        as_.AltPairSearch(fr, IndexMap);
+
+        while (as_.findNextPair()){
+            if(isHbondExist(IndexMap[as_.getResiI()], IndexMap[as_.getResiJ()], fr, pbc)){
+                IndexMap[as_.getResiI()].setHBondedPair(*(IndexMap[as_.getResiJ()].getInfo()));
+            }
+        }
+
+        break;
+    }
+    default: {
+        for(std::size_t i {0}; i < IndexMap.size(); ++i){
+            for(std::size_t j {0}; j < IndexMap.size() && j != i; ++j){
+                if(isHbondExist(IndexMap[i], IndexMap[j], fr, pbc)){
+                    IndexMap[i].setHBondedPair(*(IndexMap[j].getInfo()));
+                }
+            }
+        }
+        break;
+    }
+    }
 
    PatternSearch.initiateSearch(IndexMap, true); // TODO
    calculateBends(fr, pbc);
 
-   Storage.storageData(frnr, PatternSearch.calculateLinuxSucks(HBondsMap));
+   Storage.storageData(frnr, PatternSearch.patternSearch());
 
 }
 

diff --git a/src/dssptools.h b/src/dssptools.h
index f52688e2f8859cb51e5e5beb2de5366c6611cefd..d3e984040dde1371b247e722b846dbf8df174681 100644 (file)
--- a/src/dssptools.h
+++ b/src/dssptools.h
@@ -156,8 +156,8 @@ enum class HelixPositions : std::size_t {
 
 enum class bridgeTypes : std::size_t {
     None = 0,
-    AntiBridge,
-    Bridge,
+    AntiParallelBridge,
+    ParallelBridge,
     Count
 };
 
@@ -165,9 +165,7 @@ class secondaryStructures{
 public:
     secondaryStructures();
     void                    initiateSearch(const std::vector<ResInfo> &ResInfoMatrix, const bool PiHelicesPreferencez = true);
-    std::string             patternSearch(const std::vector<std::vector<bool>> &HBondsMap);
-//    void                    setStatus(const std::size_t resiIndex, const secondaryStructureTypes secondaryStructureTypeName, const bool secondaryStructureStatus);
-//    bool                    getStatus(const std::size_t resiIndex, const secondaryStructureTypes secondaryStructureTypeName) const;
+    std::string             patternSearch();
     ~secondaryStructures();
 
     class secondaryStructuresData{
@@ -231,18 +229,13 @@ class DsspTool
 {
 public:
    DsspTool();
-   void setInitValues(const initParameters &initParams);
-   void initAnalysis(/*const TrajectoryAnalysisSettings &settings,*/const TopologyInformation& top);
+   void initAnalysis(/*const TrajectoryAnalysisSettings &settings,*/const TopologyInformation& top, const initParameters &initParamz);
 
-   void analyzeFrame(int frnr, const t_trxframe& fr, t_pbc* pbc, const bool &alternateSearch);
+   void analyzeFrame(int frnr, const t_trxframe& fr, t_pbc* pbc);
 
    std::vector<std::pair<int, std::string>> getData();
 private:
-   real                                  cutoff_ = 4.0;
-   Selection                             sel_;
-   modeType                              modeType_;
-   NBSearchMethod                        NBSM_;
-   bool                                  verbose_;
+   initParameters                        initParams;
    AtomsDataPtr                          atoms_;
    std::string                           filename_;
    void                                  calculateBends(const t_trxframe& fr, const t_pbc* pbc);