added sel_list and masses

diff --git a/src/spirals.cpp b/src/spirals.cpp
index 3fb3a6b5e9e0e2c136b56c13cfc59eaec7e498b6..94ba8d7a247f59b1b0fd38473df59c4ee6c3a122 100644 (file)
--- a/src/spirals.cpp
+++ b/src/spirals.cpp
@@ -134,7 +134,9 @@ class Spirals : public TrajectoryAnalysisModule
 
 
         Selection                                           selec;
+        SelectionList                                       sel_;
         std::vector< std::vector < RVec > >                 trajectory;
+        std::vector< std::vector< int > >                   index_all;
         std::vector< int >                                  index;
         int                                                 frames       = 0;
 };
@@ -155,14 +157,14 @@ Spirals::initOptions(IOptionsContainer          *options,
 
     // Add the descriptive text (program help text) to the options
     settings->setHelpText(desc);
-    // Add option for selecting a subset of atoms
-    options->addOption(SelectionOption("select")
-                           .store(&selec).required()
-                           .description("Atoms that are considered as part of the excluded volume"));
     // Add option for output file name
     options->addOption(FileNameOption("on").filetype(eftIndex).outputFile()
                             .store(&fnNdx_).defaultBasename("rcore")
                             .description("Index file from the rcore"));
+    // Add option for selection list
+    options->addOption(SelectionOption("Select groups").storeVector(&sel_)
+                           .required().dynamicMask().multiValue()
+                           .description("Domains to form rigid skeleton"));
 }
 
 // -s '/home/toluk/Develop/samples/reca_rd/reca_rd.mono.tpr' -f '/home/toluk/Develop/samples/reca_rd/reca_rd.mono.xtc' -n '/home/toluk/Develop/samples/reca_rd/test.ndx' -on '/home/toluk/Develop/samples/reca_rd/core.ndx'
@@ -171,11 +173,17 @@ void
 Spirals::initAnalysis(const TrajectoryAnalysisSettings &settings,
                                const TopologyInformation         & /*top*/)
 {
+    std::vector< int > a;
+    a.resize(0);
+    index_all.resize(0);
     index.resize(0);
-    ConstArrayRef< int > atomind = selec.atomIndices();
-    for (ConstArrayRef< int >::iterator ai = atomind.begin(); (ai < atomind.end()); ai++)
-    {
-        index.push_back(*ai);
+    for (int i = 0; i < sel_.size(); i++) {
+        index_all.push_back(a);
+        ConstArrayRef< int > atomind = sel_[i].atomIndices();
+        for (ConstArrayRef< int >::iterator ai = atomind.begin(); (ai < atomind.end()); ai++) {
+            index_all.back().push_back(*ai);
+            index.push_back(*ai);
+        }
     }
 }
 
@@ -224,11 +232,33 @@ Spirals::finishAnalysis(int /*nframes*/)
         }
     }
 
+    std::vector< std::vector< long double > > masses;
+    masses.resize(0);
+    int count = 0;
+    for (int i = 0; i < index_all.size(); i++) {
+        std::vector< long double > center;
+        center.resize(0, 3);
+        for (int j = 0; j < index_all[i].size(); j++) {
+            center[0] += trajectory[0][count][0];
+            center[1] += trajectory[0][count][1];
+            center[2] += trajectory[0][count][2];
+            count++;
+        }
+        center[0] /= index_all[i].size();
+        center[1] /= index_all[i].size();
+        center[2] /= index_all[i].size();
+        masses.push_back(center);
+    }
+
+
     for (int i = 0; i < 1; i++) {
-        for (int j = 0; j < index.size() - 2; j++) {
-            test = return_crcl( trajectory[i][j][0], trajectory[i][j][1], trajectory[i][j][2],
+        for (int j = 0; j < masses.size() - 2; j++) {
+            /*test = return_crcl( trajectory[i][j][0], trajectory[i][j][1], trajectory[i][j][2],
                                 trajectory[i][j + 1][0], trajectory[i][j + 1][1], trajectory[i][j + 1][2],
-                                trajectory[i][j + 2][0], trajectory[i][j + 2][1], trajectory[i][j + 2][2]);
+                                trajectory[i][j + 2][0], trajectory[i][j + 2][1], trajectory[i][j + 2][2]);*/
+            test = return_crcl( masses[j][0], masses[j][1], masses[j][2],
+                                masses[j + 1][0], masses[j + 1][1], masses[j + 1][2],
+                                masses[j + 2][0], masses[j + 2][1], masses[j + 2][2]);
             std::cout << test.x << " " << test.y << " " << test.z << " " << test.r << "\n";
         }
     }