Changeset 375dcf for molecuilder/src/analyzer.cpp

Timestamp:

Jun 27, 2008, 9:14:58 PM (17 years ago)

Author:

Frederik Heber <heber@…>

Children:

74baec

Parents:

d7b1bf

git-author:

Frederik Heber <heber@…> (06/27/08 21:12:52)

git-committer:

Frederik Heber <heber@…> (06/27/08 21:14:58)

Message:

introduced shieldings to analyzer and joiner

both now handle pcp.sigma_all...csv files just as pcp.forces.all. Therefore the data format in pcp/perturbed.c was adapted a bit, as we need a header.
periodentafel.hpp got periodentafel and element class from molecules.hpp

File:

• molecuilder/src/analyzer.cpp (modified) (7 diffs)

molecuilder/src/analyzer.cpp

@@ -8 +8 @@
 //============================ INCLUDES ===========================
 
+#include "datacreator.hpp"
 #include "helpers.hpp" 
 #include "parser.hpp"
-#include "datacreator.hpp"
+#include "periodentafel.hpp" 
+
+// include config.h
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
 
 
@@ -17 +23 @@
 int main(int argc, char **argv)
 {
+  periodentafel *periode = NULL; // and a period table of all elements
   EnergyMatrix Energy;
   ForceMatrix Force;
+  ForceMatrix Shielding;
+  ForceMatrix ShieldingPAS;
   MatrixContainer Time;
   EnergyMatrix EnergyFragments;
   ForceMatrix ForceFragments;
-  KeySetsContainer KeySet;  
+  ForceMatrix ShieldingFragments;
+  ForceMatrix ShieldingPASFragments;
+  KeySetsContainer KeySet;
   ofstream output;
   ofstream output2;
   ofstream output3;
   ofstream output4;
+  ifstream input;
   stringstream filename;
   time_t t = time(NULL);
@@ -39 +51 @@
   // Get the command line options
   if (argc < 4) {
-    cout << "Usage: " << argv[0] << " <inputdir> <prefix> <outputdir>" << endl;
+    cout << "Usage: " << argv[0] << " <inputdir> <prefix> <outputdir> [elementsdb]" << endl;
     cout << "<inputdir>\ttherein the output of a molecuilder fragmentation is expected, each fragment with a subdir containing an energy.all and a forces.all file." << endl;
     cout << "<prefix>\tprefix of energy and forces file." << endl;
     cout << "<outputdir>\tcreated plotfiles and datafiles are placed into this directory " << endl;
+    cout << "[elementsdb]\tpath to elements database, needed for shieldings." << endl;
     return 1;
+  }
+  if (argc > 4) {
+    cout << "Loading periodentafel." << endl;
+    periode = new periodentafel;
+    periode->LoadPeriodentafel(argv[4]);
   }
   
@@ -56 +74 @@
   if (!Force.ParseMatrix(argv[1], argv[2], ForcesSuffix,0,0)) return 1;
   if (!Time.ParseMatrix(argv[1], argv[2], TimeSuffix, 10,1)) return 1;
+  if (periode != NULL) { // also look for PAS values
+    if (!Shielding.ParseMatrix(argv[1], argv[2], ShieldingSuffix, 1, 0)) return 1;
+    if (!ShieldingPAS.ParseMatrix(argv[1], argv[2], ShieldingPASSuffix, 1, 0)) return 1;
+  }
 
   // ---------- Parse the TE Factors into an array -----------------
@@ -62 +84 @@
   // ---------- Parse the Force indices into an array ---------------
   if (!Force.ParseIndices(argv[1])) return 1;
+
+  // ---------- Parse the shielding indices into an array ---------------
+  if (periode != NULL) { // also look for PAS values
+    if(!Shielding.ParseIndices(argv[1])) return 1;
+    if(!ShieldingPAS.ParseIndices(argv[1])) return 1;
+  }
 
   // ---------- Parse the KeySets into an array ---------------
   if (!KeySet.ParseKeySets(argv[1], Force.RowCounter, Force.MatrixCounter)) return 1;
   if (!KeySet.ParseManyBodyTerms()) return 1;
+  
+  // ---------- Parse fragment files created by 'joiner' into an array -------------
   if (!EnergyFragments.ParseMatrix(argv[1], argv[2], EnergyFragmentSuffix,0,0)) return 1;
   if (!ForceFragments.ParseMatrix(argv[1], argv[2], ForceFragmentSuffix,0,0)) return 1;
+  if (periode != NULL) { // also look for PAS values
+    if (!ShieldingFragments.ParseMatrix(argv[1], argv[2], ShieldingSuffix, 1, 0)) return 1;
+    if (!ShieldingPASFragments.ParseMatrix(argv[1], argv[2], ShieldingPASSuffix, 1, 0)) return 1;
+  }
 
   // +++++++++++++++ TESTING ++++++++++++++++++++++++++++++
@@ -91 +125 @@
   output << endl;
 
+  if (periode != NULL) { // also look for PAS values
+    output << endl << "Total Shieldings" << endl << "===============" << endl << Shielding.Header << endl;
+    for(int j=0;j<Shielding.RowCounter[Shielding.MatrixCounter];j++) {
+      for(int k=0;k<Shielding.ColumnCounter;k++)
+        output << scientific << Shielding.Matrix[ Shielding.MatrixCounter ][j][k] << "\t";
+      output << endl;
+    }
+    output << endl;
+  
+    output << endl << "Total Shieldings PAS" << endl << "===============" << endl << ShieldingPAS.Header << endl;
+    for(int j=0;j<ShieldingPAS.RowCounter[ShieldingPAS.MatrixCounter];j++) {
+      for(int k=0;k<ShieldingPAS.ColumnCounter;k++)
+        output << scientific << ShieldingPAS.Matrix[ ShieldingPAS.MatrixCounter ][j][k] << "\t";
+      output << endl;
+    }
+    output << endl;
+  }
+  
   output << endl << "Total Times" << endl << "===============" << endl << Time.Header << endl;
   Time.SetLastMatrix(0., 0);
@@ -380 +432 @@
   output.close();
  
-  // ++++++++++++++++ exit ++++++++++++++++++++++++++++++++++  
+  // ++++++++++++++++ exit ++++++++++++++++++++++++++++++++++
+  delete(periode);
   cout << "done." << endl;
   return 0;