Changeset 6ac7ee for src/datacreator.cpp


Ignore:
Timestamp:
Feb 9, 2009, 5:24:10 PM (17 years ago)
Author:
Frederik Heber <heber@…>
Branches:
Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, Candidate_v1.7.0, Candidate_v1.7.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable
Children:
d8b94a
Parents:
124df1
git-author:
Frederik Heber <heber@…> (02/09/09 15:55:37)
git-committer:
Frederik Heber <heber@…> (02/09/09 17:24:10)
Message:

Merge branch 'ConcaveHull' of ../espack2 into ConcaveHull

Conflicts:

molecuilder/src/boundary.cpp
molecuilder/src/boundary.hpp
molecuilder/src/builder.cpp
molecuilder/src/linkedcell.cpp
molecuilder/src/linkedcell.hpp
molecuilder/src/vector.cpp
molecuilder/src/vector.hpp
util/src/NanoCreator.c

Basically, this resulted from a lot of conversions two from spaces to one tab, which is my standard indentation. The mess was caused by eclipse auto-indenting. And in espack2:ConcaveHull was the new stuff, so all from ConcaveHull was replaced in case of doubt.
Additionally, vector had ofstream << operator instead ostream << ...

File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/datacreator.cpp

    • Property mode changed from 100644 to 100755
    r124df1 r6ac7ee  
    11/** \file datacreator.cpp
    22 *
    3  * Declarations of assisting functions in creating data and plot files. 
    4  *   
     3 * Declarations of assisting functions in creating data and plot files.
     4 *
    55 */
    66
     
    1919bool OpenOutputFile(ofstream &output, const char *dir, const char *filename)
    2020{
    21   stringstream name;
    22   name << dir << "/" << filename;
    23   output.open(name.str().c_str(), ios::out);
    24   if (output == NULL) {
    25     cout << "Unable to open " << name.str() << " for writing, is directory correct?" << endl;
    26     return false;
    27   }
    28   return true;
    29 }; 
     21        stringstream name;
     22        name << dir << "/" << filename;
     23        output.open(name.str().c_str(), ios::out);
     24        if (output == NULL) {
     25                cout << "Unable to open " << name.str() << " for writing, is directory correct?" << endl;
     26                return false;
     27        }
     28        return true;
     29};
    3030
    3131/** Opens a file for appending with \a *filename in \a *dir.
     
    3737bool AppendOutputFile(ofstream &output, const char *dir, const char *filename)
    3838{
    39   stringstream name;
    40   name << dir << "/" << filename;
    41   output.open(name.str().c_str(), ios::app);
    42   if (output == NULL) {
    43     cout << "Unable to open " << name.str() << " for writing, is directory correct?" << endl;
    44     return false;
    45   }
    46   return true;
    47 }; 
     39        stringstream name;
     40        name << dir << "/" << filename;
     41        output.open(name.str().c_str(), ios::app);
     42        if (output == NULL) {
     43                cout << "Unable to open " << name.str() << " for writing, is directory correct?" << endl;
     44                return false;
     45        }
     46        return true;
     47};
    4848
    4949/** Plots an energy vs. order.
     
    5454 * \return true if file was written successfully
    5555 */
    56 bool CreateDataEnergyOrder(class EnergyMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum) 
    57 {
    58   stringstream filename;
    59   ofstream output;
    60 
    61   filename << prefix << ".dat";
    62   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    63   cout << msg << endl;
    64   output << "# " << msg << ", created on " << datum;
    65   output << "#Order\tFrag.No.\t" << Fragments.Header << endl;
    66   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
    67     for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
    68       for(int j=Fragments.RowCounter[ KeySet.OrderSet[BondOrder][i] ];j--;)
    69         for(int k=Fragments.ColumnCounter;k--;)
    70           Fragments.Matrix[Fragments.MatrixCounter][j][k] += Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    71     }
    72     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    73     for (int l=0;l<Fragments.ColumnCounter;l++)
    74       output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l];
    75     output << endl;
    76   }
    77   output.close();
    78   return true;
     56bool CreateDataEnergyOrder(class EnergyMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum)
     57{
     58        stringstream filename;
     59        ofstream output;
     60
     61        filename << prefix << ".dat";
     62        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     63        cout << msg << endl;
     64        output << "# " << msg << ", created on " << datum;
     65        output << "#Order\tFrag.No.\t" << Fragments.Header << endl;
     66        for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     67                for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
     68                        for(int j=Fragments.RowCounter[ KeySet.OrderSet[BondOrder][i] ];j--;)
     69                                for(int k=Fragments.ColumnCounter;k--;)
     70                                        Fragments.Matrix[Fragments.MatrixCounter][j][k] += Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
     71                }
     72                output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
     73                for (int l=0;l<Fragments.ColumnCounter;l++)
     74                        output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l];
     75                output << endl;
     76        }
     77        output.close();
     78        return true;
    7979};
    8080
     
    8787 * \return true if file was written successfully
    8888 */
    89 bool CreateDataDeltaEnergyOrder(class EnergyMatrix &Energy, class EnergyMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum) 
    90 {
    91   stringstream filename;
    92   ofstream output;
    93 
    94   filename << prefix << ".dat";
    95   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    96   cout << msg << endl;
    97   output << "# " << msg << ", created on " << datum;
    98   output << "#Order\tFrag.No.\t" << Fragments.Header << endl;
    99   Fragments.SetLastMatrix(Energy.Matrix[Energy.MatrixCounter],0);
    100   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
    101     for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
    102       for(int j=Fragments.RowCounter[ KeySet.OrderSet[BondOrder][i] ];j--;)
    103         for(int k=Fragments.ColumnCounter;k--;)
    104           Fragments.Matrix[Fragments.MatrixCounter][j][k] -= Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    105     }
    106     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    107     for (int l=0;l<Fragments.ColumnCounter;l++)
    108       if (fabs(Energy.Matrix[Energy.MatrixCounter][ Energy.RowCounter[Energy.MatrixCounter]-1 ][l]) < MYEPSILON)
    109         output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l];
    110       else
    111         output << scientific << "\t" << (Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l] / Energy.Matrix[Energy.MatrixCounter][ Energy.RowCounter[Energy.MatrixCounter]-1 ][l]);
    112     output << endl;
    113   }
    114   output.close();
    115   return true;
     89bool CreateDataDeltaEnergyOrder(class EnergyMatrix &Energy, class EnergyMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum)
     90{
     91        stringstream filename;
     92        ofstream output;
     93
     94        filename << prefix << ".dat";
     95        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     96        cout << msg << endl;
     97        output << "# " << msg << ", created on " << datum;
     98        output << "#Order\tFrag.No.\t" << Fragments.Header << endl;
     99        Fragments.SetLastMatrix(Energy.Matrix[Energy.MatrixCounter],0);
     100        for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     101                for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
     102                        for(int j=Fragments.RowCounter[ KeySet.OrderSet[BondOrder][i] ];j--;)
     103                                for(int k=Fragments.ColumnCounter;k--;)
     104                                        Fragments.Matrix[Fragments.MatrixCounter][j][k] -= Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
     105                }
     106                output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
     107                for (int l=0;l<Fragments.ColumnCounter;l++)
     108                        if (fabs(Energy.Matrix[Energy.MatrixCounter][ Energy.RowCounter[Energy.MatrixCounter]-1 ][l]) < MYEPSILON)
     109                                output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l];
     110                        else
     111                                output << scientific << "\t" << (Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter]-1 ][l] / Energy.Matrix[Energy.MatrixCounter][ Energy.RowCounter[Energy.MatrixCounter]-1 ][l]);
     112                output << endl;
     113        }
     114        output.close();
     115        return true;
    116116};
    117117
     
    126126bool CreateDataForcesOrder(class ForceMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum, void (*CreateForce)(class MatrixContainer &, int))
    127127{
    128   stringstream filename;
    129   ofstream output;
    130 
    131   filename << prefix << ".dat";
    132   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    133   cout << msg << endl;
    134   output << "# " << msg << ", created on " << datum;
    135   output << "# Order\tFrag.No.\t" << Fragments.Header << endl;
    136   Fragments.SetLastMatrix(0.,0);
    137   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
    138     Fragments.SumSubForces(Fragments, KeySet, BondOrder, 1.);
    139     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    140     CreateForce(Fragments, Fragments.MatrixCounter);
    141     for (int l=0;l<Fragments.ColumnCounter;l++)
    142       output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter] ][l];
    143     output << endl;
    144   }
    145   output.close();
    146   return true;
     128        stringstream filename;
     129        ofstream output;
     130
     131        filename << prefix << ".dat";
     132        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     133        cout << msg << endl;
     134        output << "# " << msg << ", created on " << datum;
     135        output << "# Order\tFrag.No.\t" << Fragments.Header << endl;
     136        Fragments.SetLastMatrix(0.,0);
     137        for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     138                Fragments.SumSubForces(Fragments, KeySet, BondOrder, 1.);
     139                output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
     140                CreateForce(Fragments, Fragments.MatrixCounter);
     141                for (int l=0;l<Fragments.ColumnCounter;l++)
     142                        output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter] ][l];
     143                output << endl;
     144        }
     145        output.close();
     146        return true;
    147147};
    148148
     
    158158bool CreateDataDeltaForcesOrder(class ForceMatrix &Force, class ForceMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum, void (*CreateForce)(class MatrixContainer &, int))
    159159{
    160   stringstream filename;
    161   ofstream output;
    162 
    163   filename << prefix << ".dat";
    164   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    165   cout << msg << endl;
    166   output << "# " << msg << ", created on " << datum;
    167   output << "# Order\tFrag.No.\t" << Fragments.Header << endl;
    168   Fragments.SetLastMatrix(Force.Matrix[Force.MatrixCounter],0);
    169   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
    170     Fragments.SumSubForces(Fragments, KeySet, BondOrder, -1.);
    171     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    172     CreateForce(Fragments, Fragments.MatrixCounter);
    173     for (int l=0;l<Fragments.ColumnCounter;l++)
    174       output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter] ][l];
    175     output << endl;
    176   }
    177   output.close();
    178   return true;
     160        stringstream filename;
     161        ofstream output;
     162
     163        filename << prefix << ".dat";
     164        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     165        cout << msg << endl;
     166        output << "# " << msg << ", created on " << datum;
     167        output << "# Order\tFrag.No.\t" << Fragments.Header << endl;
     168        Fragments.SetLastMatrix(Force.Matrix[Force.MatrixCounter],0);
     169        for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     170                Fragments.SumSubForces(Fragments, KeySet, BondOrder, -1.);
     171                output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
     172                CreateForce(Fragments, Fragments.MatrixCounter);
     173                for (int l=0;l<Fragments.ColumnCounter;l++)
     174                        output << scientific << "\t" << Fragments.Matrix[Fragments.MatrixCounter][ Fragments.RowCounter[Fragments.MatrixCounter] ][l];
     175                output << endl;
     176        }
     177        output.close();
     178        return true;
    179179};
    180180
     
    188188 * \return true if file was written successfully
    189189 */
    190 bool CreateDataDeltaForcesOrderPerAtom(class ForceMatrix &Force, class ForceMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum) 
    191 {
    192   stringstream filename;
    193   ofstream output;
    194   double norm = 0.;
    195 
    196   filename << prefix << ".dat";
    197   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    198   cout << msg << endl;
    199   output << "# " << msg << ", created on " << datum;
    200   output << "# AtomNo\t" << Fragments.Header << endl;
    201   Fragments.SetLastMatrix(Force.Matrix[Force.MatrixCounter], 0);
    202   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
    203     //cout << "Current order is " << BondOrder << "." << endl;
    204     Fragments.SumSubForces(Fragments, KeySet, BondOrder, -1.);
    205     // errors per atom
    206     output << endl << "#Order\t" << BondOrder+1 << endl;
    207     for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
    208       output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
    209       for (int l=0;l<Fragments.ColumnCounter;l++) {
    210         if (((l+1) % 3) == 0) {
    211           norm = 0.;
    212           for (int m=0;m<NDIM;m++)
    213             norm += Force.Matrix[Force.MatrixCounter][ j ][l+m]*Force.Matrix[Force.MatrixCounter][ j ][l+m];
    214           norm = sqrt(norm);
    215         }                                                                                                           
    216 //        if (norm < MYEPSILON)
    217           output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
    218 //        else
    219 //          output << scientific << (Fragments.Matrix[Fragments.MatrixCounter][ j ][l] / norm) << "\t";
    220       }
    221       output << endl;
    222     }
    223     output << endl;
    224   }
    225   output.close();
    226   return true;
     190bool CreateDataDeltaForcesOrderPerAtom(class ForceMatrix &Force, class ForceMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum)
     191{
     192        stringstream filename;
     193        ofstream output;
     194        double norm = 0.;
     195
     196        filename << prefix << ".dat";
     197        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     198        cout << msg << endl;
     199        output << "# " << msg << ", created on " << datum;
     200        output << "# AtomNo\t" << Fragments.Header << endl;
     201        Fragments.SetLastMatrix(Force.Matrix[Force.MatrixCounter], 0);
     202        for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     203                //cout << "Current order is " << BondOrder << "." << endl;
     204                Fragments.SumSubForces(Fragments, KeySet, BondOrder, -1.);
     205                // errors per atom
     206                output << endl << "#Order\t" << BondOrder+1 << endl;
     207                for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
     208                        output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
     209                        for (int l=0;l<Fragments.ColumnCounter;l++) {
     210                                if (((l+1) % 3) == 0) {
     211                                        norm = 0.;
     212                                        for (int m=0;m<NDIM;m++)
     213                                                norm += Force.Matrix[Force.MatrixCounter][ j ][l+m]*Force.Matrix[Force.MatrixCounter][ j ][l+m];
     214                                        norm = sqrt(norm);
     215                                }                                                                                                                                                                                                                                                                                                                                                                                                               
     216//                              if (norm < MYEPSILON)
     217                                        output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
     218//                              else
     219//                                      output << scientific << (Fragments.Matrix[Fragments.MatrixCounter][ j ][l] / norm) << "\t";
     220                        }
     221                        output << endl;
     222                }
     223                output << endl;
     224        }
     225        output.close();
     226        return true;
    227227};
    228228
     
    235235 * \return true if file was written successfully
    236236 */
    237 bool CreateDataForcesOrderPerAtom(class ForceMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum) 
    238 {
    239   stringstream filename;
    240   ofstream output;
    241 
    242   filename << prefix << ".dat";
    243   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    244   cout << msg << endl;
    245   output << "# " << msg << ", created on " << datum;
    246   output << "# AtomNo\t" << Fragments.Header << endl;
    247   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
    248     //cout << "Current order is " << BondOrder << "." << endl;
    249     Fragments.SumSubForces(Fragments, KeySet, BondOrder, 1.);
    250     // errors per atom
    251     output << endl << "#Order\t" << BondOrder+1 << endl;
    252     for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
    253       output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
    254       for (int l=0;l<Fragments.ColumnCounter;l++)
    255         output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
    256       output << endl;
    257     }
    258     output << endl;
    259   }
    260   output.close();
    261   return true;
     237bool CreateDataForcesOrderPerAtom(class ForceMatrix &Fragments, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum)
     238{
     239        stringstream filename;
     240        ofstream output;
     241
     242        filename << prefix << ".dat";
     243        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     244        cout << msg << endl;
     245        output << "# " << msg << ", created on " << datum;
     246        output << "# AtomNo\t" << Fragments.Header << endl;
     247        for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     248                //cout << "Current order is " << BondOrder << "." << endl;
     249                Fragments.SumSubForces(Fragments, KeySet, BondOrder, 1.);
     250                // errors per atom
     251                output << endl << "#Order\t" << BondOrder+1 << endl;
     252                for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
     253                        output << Fragments.Indices[Fragments.MatrixCounter][j] << "\t";
     254                        for (int l=0;l<Fragments.ColumnCounter;l++)
     255                                output << scientific << Fragments.Matrix[Fragments.MatrixCounter][ j ][l] << "\t";
     256                        output << endl;
     257                }
     258                output << endl;
     259        }
     260        output.close();
     261        return true;
    262262};
    263263
     
    266266bool CreateDataFragment(class MatrixContainer &Fragment, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum, void (*CreateFragment)(class MatrixContainer &, int))
    267267{
    268   stringstream filename;
    269   ofstream output;
    270 
    271   filename << prefix << ".dat";
    272   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    273   cout << msg << endl;
    274   output << "# " << msg << ", created on " << datum << endl;
    275   output << "#Order\tFrag.No.\t" << Fragment.Header << endl;
    276   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
    277     for(int i=0;i<KeySet.FragmentsPerOrder[BondOrder];i++) {
    278       output << BondOrder+1 << "\t" << KeySet.OrderSet[BondOrder][i]+1;
    279       CreateFragment(Fragment, KeySet.OrderSet[BondOrder][i]);
    280       for (int l=0;l<Fragment.ColumnCounter;l++)
    281         output << scientific << "\t" << Fragment.Matrix[ KeySet.OrderSet[BondOrder][i] ][ Fragment.RowCounter[ KeySet.OrderSet[BondOrder][i] ] ][l];
    282       output << endl;
    283     }
    284   }
    285   output.close();
    286   return true;
     268        stringstream filename;
     269        ofstream output;
     270
     271        filename << prefix << ".dat";
     272        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     273        cout << msg << endl;
     274        output << "# " << msg << ", created on " << datum << endl;
     275        output << "#Order\tFrag.No.\t" << Fragment.Header << endl;
     276        for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     277                for(int i=0;i<KeySet.FragmentsPerOrder[BondOrder];i++) {
     278                        output << BondOrder+1 << "\t" << KeySet.OrderSet[BondOrder][i]+1;
     279                        CreateFragment(Fragment, KeySet.OrderSet[BondOrder][i]);
     280                        for (int l=0;l<Fragment.ColumnCounter;l++)
     281                                output << scientific << "\t" << Fragment.Matrix[ KeySet.OrderSet[BondOrder][i] ][ Fragment.RowCounter[ KeySet.OrderSet[BondOrder][i] ] ][l];
     282                        output << endl;
     283                }
     284        }
     285        output.close();
     286        return true;
    287287};
    288288
     
    291291 * \param &KeySet KeySetsContainer with associations of each fragment to a bond order
    292292 * \param BondOrder current bond order
    293  */ 
     293 */
    294294void CreateMaxFragmentOrder(class MatrixContainer &Fragments, class KeySetsContainer &KeySet, int BondOrder)
    295295{
    296   for(int j=Fragments.RowCounter[ Fragments.MatrixCounter ];j--;) {
    297     for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
    298       if (fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) < fabs(Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][1])) {
    299         for (int k=Fragments.ColumnCounter;k--;)
    300           Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    301       }
    302     }     
    303   }
     296        for(int j=Fragments.RowCounter[ Fragments.MatrixCounter ];j--;) {
     297                for(int i=KeySet.FragmentsPerOrder[BondOrder];i--;) {
     298                        if (fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) < fabs(Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][1])) {
     299                                for (int k=Fragments.ColumnCounter;k--;)
     300                                        Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
     301                        }
     302                }
     303        }
    304304};
    305305
     
    308308 * \param &KeySet KeySetsContainer with associations of each fragment to a bond order
    309309 * \param BondOrder current bond order
    310  */ 
     310 */
    311311void CreateMinFragmentOrder(class MatrixContainer &Fragments, class KeySetsContainer &KeySet, int BondOrder)
    312312{
    313   for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
    314     int i=0;
    315     do {  // first get a minimum value unequal to 0
    316       for (int k=Fragments.ColumnCounter;k--;)
    317         Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    318       i++;
    319     } while ((fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) < MYEPSILON) && (i<KeySet.FragmentsPerOrder[BondOrder]));
    320     for(;i<KeySet.FragmentsPerOrder[BondOrder];i++) { // then find lowest
    321       if (fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) > fabs(Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][1])) {
    322         for (int k=Fragments.ColumnCounter;k--;)
    323           Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
    324       }
    325     }     
    326   }
     313        for(int j=0;j<Fragments.RowCounter[ Fragments.MatrixCounter ];j++) {
     314                int i=0;
     315                do {    // first get a minimum value unequal to 0
     316                        for (int k=Fragments.ColumnCounter;k--;)
     317                                Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
     318                        i++;
     319                } while ((fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) < MYEPSILON) && (i<KeySet.FragmentsPerOrder[BondOrder]));
     320                for(;i<KeySet.FragmentsPerOrder[BondOrder];i++) { // then find lowest
     321                        if (fabs(Fragments.Matrix[ Fragments.MatrixCounter ][j][1]) > fabs(Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][1])) {
     322                                for (int k=Fragments.ColumnCounter;k--;)
     323                                        Fragments.Matrix[ Fragments.MatrixCounter ][j][k] = Fragments.Matrix[ KeySet.OrderSet[BondOrder][i] ][j][k];
     324                        }
     325                }
     326        }
    327327};
    328328
     
    331331bool CreateDataFragmentOrder(class MatrixContainer &Fragment, class KeySetsContainer &KeySet, char *dir, char *prefix, char *msg, char *datum, void (*CreateFragmentOrder)(class MatrixContainer &, class KeySetsContainer &, int))
    332332{
    333   stringstream filename;
    334   ofstream output;
    335 
    336   filename << prefix << ".dat";
    337   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    338   cout << msg << endl;
    339   output << "# " << msg << ", created on " << datum;
    340   output << "#Order\tFrag.No.\t" << Fragment.Header << endl;
    341   // max
    342   for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
    343     Fragment.SetLastMatrix(0.,0);
    344     CreateFragmentOrder(Fragment, KeySet, BondOrder);
    345     output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
    346     for (int l=0;l<Fragment.ColumnCounter;l++)
    347       output << scientific << "\t" << Fragment.Matrix[ Fragment.MatrixCounter ][ Fragment.RowCounter[ Fragment.MatrixCounter ]-1 ][l];
    348     output << endl;
    349   }
    350   output.close();
    351   return true;
     333        stringstream filename;
     334        ofstream output;
     335
     336        filename << prefix << ".dat";
     337        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     338        cout << msg << endl;
     339        output << "# " << msg << ", created on " << datum;
     340        output << "#Order\tFrag.No.\t" << Fragment.Header << endl;
     341        // max
     342        for (int BondOrder=0;BondOrder<KeySet.Order;BondOrder++) {
     343                Fragment.SetLastMatrix(0.,0);
     344                CreateFragmentOrder(Fragment, KeySet, BondOrder);
     345                output << BondOrder+1 << "\t" << KeySet.FragmentsPerOrder[BondOrder];
     346                for (int l=0;l<Fragment.ColumnCounter;l++)
     347                        output << scientific << "\t" << Fragment.Matrix[ Fragment.MatrixCounter ][ Fragment.RowCounter[ Fragment.MatrixCounter ]-1 ][l];
     348                output << endl;
     349        }
     350        output.close();
     351        return true;
    352352};
    353353
     
    358358void CreateEnergy(class MatrixContainer &Energy, int MatrixNumber)
    359359{
    360   for(int k=0;k<Energy.ColumnCounter;k++)
    361     Energy.Matrix[MatrixNumber][ Energy.RowCounter[MatrixNumber] ] [k] =  Energy.Matrix[MatrixNumber][ Energy.RowCounter[MatrixNumber]-1 ] [k];
     360        for(int k=0;k<Energy.ColumnCounter;k++)
     361                Energy.Matrix[MatrixNumber][ Energy.RowCounter[MatrixNumber] ] [k] =    Energy.Matrix[MatrixNumber][ Energy.RowCounter[MatrixNumber]-1 ] [k];
    362362};
    363363
     
    369369void CreateMinimumForce(class MatrixContainer &Force, int MatrixNumber)
    370370{
    371   for (int l=Force.ColumnCounter;l--;)
    372     Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
    373   for (int l=5;l<Force.ColumnCounter;l+=3) {
    374     double stored = 0;
    375     int k=0;
    376     do {
    377       for (int m=NDIM;m--;) {
    378         stored += Force.Matrix[MatrixNumber][ k ][l+m]
    379               * Force.Matrix[MatrixNumber][ k ][l+m];
    380         Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l+m]  = Force.Matrix[MatrixNumber][ k ][l+m];
    381       }
    382       k++;
    383     } while ((fabs(stored) < MYEPSILON) && (k<Force.RowCounter[MatrixNumber]));
    384     for (;k<Force.RowCounter[MatrixNumber];k++) {
    385       double tmp = 0;
    386       for (int m=NDIM;m--;)
    387         tmp += Force.Matrix[MatrixNumber][ k ][l+m]
    388               * Force.Matrix[MatrixNumber][ k ][l+m];
    389       if ((fabs(tmp) > MYEPSILON) && (tmp < stored)) {  // current force is greater than stored
    390         for (int m=NDIM;m--;)
    391           Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l+m]  = Force.Matrix[MatrixNumber][ k ][l+m];
    392         stored = tmp;
    393       }
    394     }
    395   }
     371        for (int l=Force.ColumnCounter;l--;)
     372                Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
     373        for (int l=5;l<Force.ColumnCounter;l+=3) {
     374                double stored = 0;
     375                int k=0;
     376                do {
     377                        for (int m=NDIM;m--;) {
     378                                stored += Force.Matrix[MatrixNumber][ k ][l+m]
     379                                                        * Force.Matrix[MatrixNumber][ k ][l+m];
     380                                Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l+m]       = Force.Matrix[MatrixNumber][ k ][l+m];
     381                        }
     382                        k++;
     383                } while ((fabs(stored) < MYEPSILON) && (k<Force.RowCounter[MatrixNumber]));
     384                for (;k<Force.RowCounter[MatrixNumber];k++) {
     385                        double tmp = 0;
     386                        for (int m=NDIM;m--;)
     387                                tmp += Force.Matrix[MatrixNumber][ k ][l+m]
     388                                                        * Force.Matrix[MatrixNumber][ k ][l+m];
     389                        if ((fabs(tmp) > MYEPSILON) && (tmp < stored)) {        // current force is greater than stored
     390                                for (int m=NDIM;m--;)
     391                                        Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l+m]       = Force.Matrix[MatrixNumber][ k ][l+m];
     392                                stored = tmp;
     393                        }
     394                }
     395        }
    396396};
    397397
     
    399399 * Results are stored in the matrix ForceMatrix::MatrixCounter of \a Force.
    400400 * \param Force ForceMatrix class containing matrix values
    401   * \param MatrixNumber the index for the ForceMatrix::matrix array
     401        * \param MatrixNumber the index for the ForceMatrix::matrix array
    402402 */
    403403void CreateMeanForce(class MatrixContainer &Force, int MatrixNumber)
    404404{
    405   int divisor = 0;
    406   for (int l=Force.ColumnCounter;l--;)
    407     Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
    408   for (int l=5;l<Force.ColumnCounter;l+=3) {
    409     double tmp = 0;
    410     for (int k=Force.RowCounter[MatrixNumber];k--;) {
    411       double norm = 0.;
    412       for (int m=NDIM;m--;)
    413         norm += Force.Matrix[MatrixNumber][ k ][l+m]
    414               * Force.Matrix[MatrixNumber][ k ][l+m];
    415       tmp += sqrt(norm);
    416       if (fabs(norm) > MYEPSILON) divisor++;
    417     }
    418     tmp /= (double)divisor;
    419     Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = tmp;
    420   }
     405        int divisor = 0;
     406        for (int l=Force.ColumnCounter;l--;)
     407                Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
     408        for (int l=5;l<Force.ColumnCounter;l+=3) {
     409                double tmp = 0;
     410                for (int k=Force.RowCounter[MatrixNumber];k--;) {
     411                        double norm = 0.;
     412                        for (int m=NDIM;m--;)
     413                                norm += Force.Matrix[MatrixNumber][ k ][l+m]
     414                                                        * Force.Matrix[MatrixNumber][ k ][l+m];
     415                        tmp += sqrt(norm);
     416                        if (fabs(norm) > MYEPSILON) divisor++;
     417                }
     418                tmp /= (double)divisor;
     419                Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = tmp;
     420        }
    421421};
    422422
     
    428428void CreateMaximumForce(class MatrixContainer &Force, int MatrixNumber)
    429429{
    430   for (int l=5;l<Force.ColumnCounter;l+=3) {
    431     double stored = 0;
    432     for (int k=Force.RowCounter[MatrixNumber];k--;) {
    433       double tmp = 0;
    434       for (int m=NDIM;m--;)
    435         tmp += Force.Matrix[MatrixNumber][ k ][l+m]
    436               * Force.Matrix[MatrixNumber][ k ][l+m];
    437       if (tmp > stored) {  // current force is greater than stored
    438         for (int m=NDIM;m--;)
    439           Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l+m]  = Force.Matrix[MatrixNumber][ k ][l+m];
    440         stored = tmp;
    441       }
    442     }
    443   }
     430        for (int l=5;l<Force.ColumnCounter;l+=3) {
     431                double stored = 0;
     432                for (int k=Force.RowCounter[MatrixNumber];k--;) {
     433                        double tmp = 0;
     434                        for (int m=NDIM;m--;)
     435                                tmp += Force.Matrix[MatrixNumber][ k ][l+m]
     436                                                        * Force.Matrix[MatrixNumber][ k ][l+m];
     437                        if (tmp > stored) {     // current force is greater than stored
     438                                for (int m=NDIM;m--;)
     439                                        Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l+m]       = Force.Matrix[MatrixNumber][ k ][l+m];
     440                                stored = tmp;
     441                        }
     442                }
     443        }
    444444};
    445445
     
    450450void CreateSameForce(class MatrixContainer &Force, int MatrixNumber)
    451451{
    452   // does nothing
     452        // does nothing
    453453};
    454454
     
    460460void CreateVectorSumForce(class MatrixContainer &Force, int MatrixNumber)
    461461{
    462   for (int l=Force.ColumnCounter;l--;)
    463     Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
    464   for (int l=0;l<Force.ColumnCounter;l++) {
    465     for (int k=Force.RowCounter[MatrixNumber];k--;)
    466       Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] += Force.Matrix[MatrixNumber][k][l];
    467   }
     462        for (int l=Force.ColumnCounter;l--;)
     463                Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] = 0.;
     464        for (int l=0;l<Force.ColumnCounter;l++) {
     465                for (int k=Force.RowCounter[MatrixNumber];k--;)
     466                        Force.Matrix[MatrixNumber][ Force.RowCounter[MatrixNumber] ][l] += Force.Matrix[MatrixNumber][k][l];
     467        }
    468468};
    469469
     
    472472 * \param *key position of key
    473473 * \param *logscale axis for logscale
    474  * \param *extraline extra set lines if desired 
     474 * \param *extraline extra set lines if desired
    475475 * \param mxtics small tics at ...
    476476 * \param xtics large tics at ...
    477  * \param *xlabel label for x axis 
     477 * \param *xlabel label for x axis
    478478 * \param *ylabel label for y axis
    479  */ 
     479 */
    480480void CreatePlotHeader(ofstream &output, const char *prefix, const int keycolumns, const char *key, const char *logscale, const char *extraline, const int mxtics, const int xtics, const char *xlabel, const char *ylabel)
    481481{
    482   //output << "#!/home/heber/build/pyxplot/pyxplot" << endl << endl;
    483   output << "reset" << endl;
    484   output << "set keycolumns "<< keycolumns << endl;
    485   output << "set key " << key << endl;
    486   output << "set mxtics "<< mxtics << endl;
    487   output << "set xtics "<< xtics << endl;
    488   if (logscale != NULL)
    489     output << "set logscale " << logscale << endl;
    490   if (extraline != NULL)
    491     output << extraline << endl;
    492   output << "set xlabel '" << xlabel << "'" << endl;
    493   output << "set ylabel '" << ylabel << "'" << endl;
    494   output << "set terminal eps color" << endl;
    495   output << "set output '"<< prefix << ".eps'" << endl;
     482        //output << "#!/home/heber/build/pyxplot/pyxplot" << endl << endl;
     483        output << "reset" << endl;
     484        output << "set keycolumns "<< keycolumns << endl;
     485        output << "set key " << key << endl;
     486        output << "set mxtics "<< mxtics << endl;
     487        output << "set xtics "<< xtics << endl;
     488        if (logscale != NULL)
     489                output << "set logscale " << logscale << endl;
     490        if (extraline != NULL)
     491                output << extraline << endl;
     492        output << "set xlabel '" << xlabel << "'" << endl;
     493        output << "set ylabel '" << ylabel << "'" << endl;
     494        output << "set terminal eps color" << endl;
     495        output << "set output '"<< prefix << ".eps'" << endl;
    496496};
    497497
     
    506506 * \param mxtics small tics at ...
    507507 * \param xtics large tics at ...
    508  * \param xlabel label for x axis 
     508 * \param xlabel label for x axis
    509509 * \param xlabel label for x axis
    510510 * \param *xrange xrange
     
    514514 * \param (*CreatePlotLines) function reference that writes a single plot line
    515515 * \return true if file was written successfully
    516  */   
     516 */
    517517bool CreatePlotOrder(class MatrixContainer &Matrix, const class KeySetsContainer &KeySet, const char *dir, const char *prefix, const int keycolumns, const char *key, const char *logscale, const char *extraline, const int mxtics, const int xtics, const char *xlabel, const char *ylabel, const char *xrange, const char *yrange, const char *xargument, const char *uses, void (*CreatePlotLines)(ofstream &, class MatrixContainer &, const char *, const char *, const char *))
    518518{
    519   stringstream filename;
    520   ofstream output;
    521 
    522   filename << prefix << ".pyx";
    523   if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
    524   CreatePlotHeader(output, prefix, keycolumns, key, logscale, extraline, mxtics, xtics, xlabel, ylabel);
    525   output << "plot " << xrange << " " << yrange << " \\" << endl;
    526   CreatePlotLines(output, Matrix, prefix, xargument, uses);
    527   output.close(); 
    528   return true;
     519        stringstream filename;
     520        ofstream output;
     521
     522        filename << prefix << ".pyx";
     523        if (!OpenOutputFile(output, dir, filename.str().c_str())) return false;
     524        CreatePlotHeader(output, prefix, keycolumns, key, logscale, extraline, mxtics, xtics, xlabel, ylabel);
     525        output << "plot " << xrange << " " << yrange << " \\" << endl;
     526        CreatePlotLines(output, Matrix, prefix, xargument, uses);
     527        output.close();
     528        return true;
    529529};
    530530
     
    538538void AbsEnergyPlotLine(ofstream &output, class MatrixContainer &Energy, const char *prefix, const char *xargument, const char *uses)
    539539{
    540   stringstream line(Energy.Header);
    541   string token;
    542 
    543   getline(line, token, '\t');
    544   for (int i=2; i<= Energy.ColumnCounter;i++) {
    545     getline(line, token, '\t');
    546     while (token[0] == ' ') // remove leading white spaces
    547       token.erase(0,1);
    548     output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(abs($" << i+2 << ")) " << uses;
    549     if (i != (Energy.ColumnCounter))
    550       output << ", \\";
    551     output << endl;
    552   }
     540        stringstream line(Energy.Header);
     541        string token;
     542
     543        getline(line, token, '\t');
     544        for (int i=2; i<= Energy.ColumnCounter;i++) {
     545                getline(line, token, '\t');
     546                while (token[0] == ' ') // remove leading white spaces
     547                        token.erase(0,1);
     548                output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(abs($" << i+2 << ")) " << uses;
     549                if (i != (Energy.ColumnCounter))
     550                        output << ", \\";
     551                output << endl;
     552        }
    553553};
    554554
     
    562562void EnergyPlotLine(ofstream &output, class MatrixContainer &Energy, const char *prefix, const char *xargument, const char *uses)
    563563{
    564   stringstream line(Energy.Header);
    565   string token;
    566 
    567   getline(line, token, '\t');
    568   for (int i=1; i<= Energy.ColumnCounter;i++) {
    569     getline(line, token, '\t');
    570     while (token[0] == ' ') // remove leading white spaces
    571       token.erase(0,1);
    572     output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":" << i+2 << " " << uses;
    573     if (i != (Energy.ColumnCounter))
    574       output << ", \\";
    575     output << endl;
    576   }
     564        stringstream line(Energy.Header);
     565        string token;
     566
     567        getline(line, token, '\t');
     568        for (int i=1; i<= Energy.ColumnCounter;i++) {
     569                getline(line, token, '\t');
     570                while (token[0] == ' ') // remove leading white spaces
     571                        token.erase(0,1);
     572                output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":" << i+2 << " " << uses;
     573                if (i != (Energy.ColumnCounter))
     574                        output << ", \\";
     575                output << endl;
     576        }
    577577};
    578578
     
    586586void ForceMagnitudePlotLine(ofstream &output, class MatrixContainer &Force, const char *prefix, const char *xargument, const char *uses)
    587587{
    588   stringstream line(Force.Header);
    589   string token;
    590 
    591   getline(line, token, '\t');
    592   getline(line, token, '\t');
    593   getline(line, token, '\t');
    594   getline(line, token, '\t');
    595   getline(line, token, '\t');
    596   for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
    597     getline(line, token, '\t');
    598     while (token[0] == ' ') // remove leading white spaces
    599       token.erase(0,1);
    600     token.erase(token.length(), 1);  // kill residual index char (the '0')
    601     output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(sqrt($" << i+1 << "*$" << i+1 << "+$" << i+2 << "*$" << i+2 << "+$" << i+3 << "*$" << i+3 << ")) " << uses;
    602     if (i != (Force.ColumnCounter-1))
    603       output << ", \\";
    604     output << endl;
    605     getline(line, token, '\t');
    606     getline(line, token, '\t');
    607   }
     588        stringstream line(Force.Header);
     589        string token;
     590
     591        getline(line, token, '\t');
     592        getline(line, token, '\t');
     593        getline(line, token, '\t');
     594        getline(line, token, '\t');
     595        getline(line, token, '\t');
     596        for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
     597                getline(line, token, '\t');
     598                while (token[0] == ' ') // remove leading white spaces
     599                        token.erase(0,1);
     600                token.erase(token.length(), 1); // kill residual index char (the '0')
     601                output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(sqrt($" << i+1 << "*$" << i+1 << "+$" << i+2 << "*$" << i+2 << "+$" << i+3 << "*$" << i+3 << ")) " << uses;
     602                if (i != (Force.ColumnCounter-1))
     603                        output << ", \\";
     604                output << endl;
     605                getline(line, token, '\t');
     606                getline(line, token, '\t');
     607        }
    608608};
    609609
     
    617617void AbsFirstForceValuePlotLine(ofstream &output, class MatrixContainer &Force, const char *prefix, const char *xargument, const char *uses)
    618618{
    619   stringstream line(Force.Header);
    620   string token;
    621 
    622   getline(line, token, '\t');
    623   getline(line, token, '\t');
    624   getline(line, token, '\t');
    625   getline(line, token, '\t');
    626   getline(line, token, '\t');
    627   for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
    628     getline(line, token, '\t');
    629     while (token[0] == ' ') // remove leading white spaces
    630       token.erase(0,1);
    631     token.erase(token.length(), 1);  // kill residual index char (the '0')
    632     output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(abs($" << i+1 << ")) " << uses;
    633     if (i != (Force.ColumnCounter-1))
    634       output << ", \\";
    635     output << endl;
    636     getline(line, token, '\t');
    637     getline(line, token, '\t');
    638   }
     619        stringstream line(Force.Header);
     620        string token;
     621
     622        getline(line, token, '\t');
     623        getline(line, token, '\t');
     624        getline(line, token, '\t');
     625        getline(line, token, '\t');
     626        getline(line, token, '\t');
     627        for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
     628                getline(line, token, '\t');
     629                while (token[0] == ' ') // remove leading white spaces
     630                        token.erase(0,1);
     631                token.erase(token.length(), 1); // kill residual index char (the '0')
     632                output << "'" << prefix << ".dat' title '" << token << "' using " << xargument << ":(abs($" << i+1 << ")) " << uses;
     633                if (i != (Force.ColumnCounter-1))
     634                        output << ", \\";
     635                output << endl;
     636                getline(line, token, '\t');
     637                getline(line, token, '\t');
     638        }
    639639};
    640640
     
    648648void BoxesForcePlotLine(ofstream &output, class MatrixContainer &Force, const char *prefix, const char *xargument, const char *uses)
    649649{
    650   stringstream line(Force.Header);
    651   char *fillcolor[5] = {"black", "red", "blue", "green", "cyan"};
    652   string token;
    653 
    654   getline(line, token, '\t');
    655   getline(line, token, '\t');
    656   getline(line, token, '\t');
    657   getline(line, token, '\t');
    658   getline(line, token, '\t');
    659   for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
    660     getline(line, token, '\t');
    661     while (token[0] == ' ') // remove leading white spaces
    662       token.erase(0,1);
    663     token.erase(token.length(), 1);  // kill residual index char (the '0')
    664     output << "'" << prefix << ".dat' title '" << token << "' using ($" << xargument << "+" << fixed << setprecision(1) << (double)((i-7)/3)*0.2 << "):(sqrt($" << i+1 << "*$" << i+1 << "+$" << i+2 << "*$" << i+2 << "+$" << i+3 << "*$" << i+3 << ")) " << uses << " " << fillcolor[(i-7)/3];
    665     if (i != (Force.ColumnCounter-1))
    666       output << ", \\";
    667     output << endl;
    668     getline(line, token, '\t');
    669     getline(line, token, '\t');
    670   }
     650        stringstream line(Force.Header);
     651        char *fillcolor[5] = {"black", "red", "blue", "green", "cyan"};
     652        string token;
     653
     654        getline(line, token, '\t');
     655        getline(line, token, '\t');
     656        getline(line, token, '\t');
     657        getline(line, token, '\t');
     658        getline(line, token, '\t');
     659        for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
     660                getline(line, token, '\t');
     661                while (token[0] == ' ') // remove leading white spaces
     662                        token.erase(0,1);
     663                token.erase(token.length(), 1); // kill residual index char (the '0')
     664                output << "'" << prefix << ".dat' title '" << token << "' using ($" << xargument << "+" << fixed << setprecision(1) << (double)((i-7)/3)*0.2 << "):(sqrt($" << i+1 << "*$" << i+1 << "+$" << i+2 << "*$" << i+2 << "+$" << i+3 << "*$" << i+3 << ")) " << uses << " " << fillcolor[(i-7)/3];
     665                if (i != (Force.ColumnCounter-1))
     666                        output << ", \\";
     667                output << endl;
     668                getline(line, token, '\t');
     669                getline(line, token, '\t');
     670        }
    671671};
    672672
     
    680680void BoxesFirstForceValuePlotLine(ofstream &output, class MatrixContainer &Force, const char *prefix, const char *xargument, const char *uses)
    681681{
    682   stringstream line(Force.Header);
    683   char *fillcolor[5] = {"black", "red", "blue", "green", "cyan"};
    684   string token;
    685 
    686   getline(line, token, '\t');
    687   getline(line, token, '\t');
    688   getline(line, token, '\t');
    689   getline(line, token, '\t');
    690   getline(line, token, '\t');
    691   for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
    692     getline(line, token, '\t');
    693     while (token[0] == ' ') // remove leading white spaces
    694       token.erase(0,1);
    695     token.erase(token.length(), 1);  // kill residual index char (the '0')
    696     output << "'" << prefix << ".dat' title '" << token << "' using ($" << xargument << "+" << fixed << setprecision(1) << (double)((i-7)/3)*0.2 << "):" << i+1 << " " << uses << " " << fillcolor[(i-7)/3];
    697     if (i != (Force.ColumnCounter-1))
    698       output << ", \\";
    699     output << endl;
    700     getline(line, token, '\t');
    701     getline(line, token, '\t');
    702   }
    703 };
     682        stringstream line(Force.Header);
     683        char *fillcolor[5] = {"black", "red", "blue", "green", "cyan"};
     684        string token;
     685
     686        getline(line, token, '\t');
     687        getline(line, token, '\t');
     688        getline(line, token, '\t');
     689        getline(line, token, '\t');
     690        getline(line, token, '\t');
     691        for (int i=7; i< Force.ColumnCounter;i+=NDIM) {
     692                getline(line, token, '\t');
     693                while (token[0] == ' ') // remove leading white spaces
     694                        token.erase(0,1);
     695                token.erase(token.length(), 1); // kill residual index char (the '0')
     696                output << "'" << prefix << ".dat' title '" << token << "' using ($" << xargument << "+" << fixed << setprecision(1) << (double)((i-7)/3)*0.2 << "):" << i+1 << " " << uses << " " << fillcolor[(i-7)/3];
     697                if (i != (Force.ColumnCounter-1))
     698                        output << ", \\";
     699                output << endl;
     700                getline(line, token, '\t');
     701                getline(line, token, '\t');
     702        }
     703};
Note: See TracChangeset for help on using the changeset viewer.