--- branches/development/src/applications/staticProps/SpatialStatistics.cpp	2013/05/15 15:09:35	1874
+++ trunk/src/applications/staticProps/SpatialStatistics.cpp	2014/10/22 12:23:59	2026
@@ -43,6 +43,10 @@
 #include "applications/staticProps/SpatialStatistics.hpp"
 #include "io/DumpReader.hpp"
 #include "primitives/Molecule.hpp"
+#ifdef _MSC_VER
+#define isnan(x) _isnan((x))
+#define isinf(x) (!_finite(x) && !_isnan(x))
+#endif
 
 namespace OpenMD {
   
@@ -137,7 +141,7 @@ namespace OpenMD {
       processStuntDouble( sd, bin );
 
       dynamic_cast<Accumulator *>(counts_->accumulator[bin])->add(1);
-    }      
+    }
   }
   
 
@@ -183,7 +187,7 @@ namespace OpenMD {
       }
         
       outStream << "#######################################################\n";
-      outStream << "# Standard Deviations in those quantities follow:\n";
+      outStream << "# 95% confidence intervals in those quantities follow:\n";
       outStream << "#######################################################\n";
       
       for (int j = 0; j < nBins_; j++) {
@@ -196,7 +200,7 @@ namespace OpenMD {
             
             int n = outputData->accumulator[j]->count();
             if (n != 0) {
-              writeStdDev( outStream, outputData, j );
+              writeErrorBars( outStream, outputData, j );
             }
           }
           outStream << std::endl;
@@ -221,15 +225,13 @@ namespace OpenMD {
     int n = dat->accumulator[bin]->count();
     if (n == 0) return;
 
-    RealType r;
-    Vector3d v;
-
     if( dat->dataType == odtReal ) {
+      RealType r;
       dynamic_cast<Accumulator*>(dat->accumulator[bin])->getAverage(r);      
       if (isinf(r) || isnan(r) ) {      
         sprintf( painCave.errMsg,
                  "SpatialStatistics detected a numerical error writing:\n"
-                 "\t%s for bin %d",
+                 "\t%s for bin %u",
                  dat->title.c_str(), bin);
         painCave.isFatal = 1;
         simError();
@@ -238,13 +240,14 @@ namespace OpenMD {
       os << "\t" << r;      
 
     } else if ( dat->dataType == odtVector3 ) {
+      Vector3d v;
       dynamic_cast<VectorAccumulator*>(dat->accumulator[bin])->getAverage(v);
       if (isinf(v[0]) || isnan(v[0]) || 
           isinf(v[1]) || isnan(v[1]) || 
           isinf(v[2]) || isnan(v[2]) ) {      
         sprintf( painCave.errMsg,
                  "SpatialStatistics detected a numerical error writing:\n"
-                 "\t%s for bin %d",
+                 "\t%s for bin %u",
                  dat->title.c_str(), bin);
         painCave.isFatal = 1;
         simError();
@@ -254,21 +257,19 @@ namespace OpenMD {
     }
   }
 
-  void SpatialStatistics::writeStdDev(ostream& os, OutputData* dat, 
+  void SpatialStatistics::writeErrorBars(ostream& os, OutputData* dat, 
                                     unsigned int bin) {
     assert(int(bin) < nBins_);
     int n = dat->accumulator[bin]->count();
     if (n == 0) return;
 
-    RealType r;
-    Vector3d v;
-
     if( dat->dataType == odtReal ) {
-      dynamic_cast<Accumulator*>(dat->accumulator[bin])->getStdDev(r);      
+      RealType r;
+      dynamic_cast<Accumulator*>(dat->accumulator[bin])->get95percentConfidenceInterval(r);      
       if (isinf(r) || isnan(r) ) {      
         sprintf( painCave.errMsg,
                  "SpatialStatistics detected a numerical error writing:\n"
-                 "\tstandard deviation of %s for bin %d",
+                 "\tstandard deviation of %s for bin %u",
                  dat->title.c_str(), bin);
         painCave.isFatal = 1;
         simError();
@@ -277,13 +278,14 @@ namespace OpenMD {
       os << "\t" << r;      
 
     } else if ( dat->dataType == odtVector3 ) {
-      dynamic_cast<VectorAccumulator*>(dat->accumulator[bin])->getStdDev(v);
+      Vector3d v;
+      dynamic_cast<VectorAccumulator*>(dat->accumulator[bin])->get95percentConfidenceInterval(v);
       if (isinf(v[0]) || isnan(v[0]) || 
           isinf(v[1]) || isnan(v[1]) || 
           isinf(v[2]) || isnan(v[2]) ) {      
         sprintf( painCave.errMsg,
                  "SpatialStatistics detected a numerical error writing:\n"
-                 "\tstandard deviation of %s for bin %d",
+                 "\tstandard deviation of %s for bin %u",
                  dat->title.c_str(), bin);
         painCave.isFatal = 1;
         simError();
@@ -321,12 +323,12 @@ namespace OpenMD {
   }
 
   SlabStatistics::~SlabStatistics() {
-    delete z_;
   }
 
 
   void SlabStatistics::processFrame(int istep) {
     RealType z;
+
     hmat_ = currentSnapshot_->getHmat();
     for (int i = 0; i < nBins_; i++) {
       z = (((RealType)i + 0.5) / (RealType)nBins_) * hmat_(2,2);
@@ -349,11 +351,31 @@ namespace OpenMD {
 
   ShellStatistics::ShellStatistics(SimInfo* info, const string& filename, 
                                    const string& sele, int nbins) : 
-    SpatialStatistics(info, filename, sele, nbins){
+    SpatialStatistics(info, filename, sele, nbins), coordinateOrigin_(V3Zero) {
     
-    coordinateOrigin_ = V3Zero;
     binWidth_ = 1.0;
 
+    Globals* simParams = info->getSimParams();
+    RNEMDParameters* rnemdParams = simParams->getRNEMDParameters();
+    bool hasCoordinateOrigin = rnemdParams->haveCoordinateOrigin();
+    
+    if (hasCoordinateOrigin) {
+        std::vector<RealType> co = rnemdParams->getCoordinateOrigin();
+        if (co.size() != 3) {
+          sprintf(painCave.errMsg,
+                  "RNEMD: Incorrect number of parameters specified for coordinateOrigin.\n"
+                  "\tthere should be 3 parameters, but %lu were specified.\n", 
+                  co.size());
+          painCave.isFatal = 1;
+          simError();      
+        }
+        coordinateOrigin_.x() = co[0];
+        coordinateOrigin_.y() = co[1];
+        coordinateOrigin_.z() = co[2];
+    } else {
+      coordinateOrigin_ = V3Zero;
+    }
+    
     r_ = new OutputData;
     r_->units =  "Angstroms";
     r_->title =  "R";
@@ -371,10 +393,9 @@ namespace OpenMD {
   }
 
   ShellStatistics::~ShellStatistics() {
-    delete r_;
   }
 
-  int ShellStatistics::getBin(Vector3d pos) {    
+  int ShellStatistics::getBin(Vector3d pos) { 
     Vector3d rPos = pos - coordinateOrigin_;
     return int(rPos.length() / binWidth_);
   }