--- trunk/src/applications/staticProps/GofR.cpp 2005/02/10 14:15:52 309 +++ trunk/src/applications/staticProps/GofR.cpp 2005/02/17 16:21:07 360 @@ -46,30 +46,33 @@ namespace oopse { namespace oopse { -GofR::GofR(SimInfo* info, const std::string& filename, const std::string& sele1, const std::string& sele2, double len) - : RadialDistrFunc(info, filename, sele1, sele2, len){ +GofR::GofR(SimInfo* info, const std::string& filename, const std::string& sele1, const std::string& sele2, double len, int nrbins) + : RadialDistrFunc(info, filename, sele1, sele2), len_(len), nRBins_(nrbins){ - histogram_.resize(nbins_); - avgGofr_.resize(nbins_); + deltaR_ = len_ /nRBins_; + + histogram_.resize(nRBins_); + avgGofr_.resize(nRBins_); + + setOutputName(getPrefix(filename) + ".gofr"); } void GofR::preProcess() { - avgGofr_.resize(nbins_); std::fill(avgGofr_.begin(), avgGofr_.end(), 0.0); } void GofR::initalizeHistogram() { - npairs_ = 0; std::fill(histogram_.begin(), histogram_.end(), 0); } void GofR::processHistogram() { + int nPairs = getNPairs(); double volume = info_->getSnapshotManager()->getCurrentSnapshot()->getVolume(); - double pairDensity = npairs_ /volume; - double pairConstant = ( 4.0 * PI * pairDensity ) / 3.0; + double pairDensity = nPairs /volume * 2.0; + double pairConstant = ( 4.0 * NumericConstant::PI * pairDensity ) / 3.0; for(int i = 0 ; i < histogram_.size(); ++i){ @@ -96,9 +99,10 @@ void GofR::collectHistogram(StuntDouble* sd1, StuntDou double distance = r12.length(); - int whichBin = distance / deltaR_; - histogram_[whichBin] ++; - npairs_++; + if (distance < len_) { + int whichBin = distance / deltaR_; + histogram_[whichBin] += 2; + } } @@ -116,7 +120,9 @@ void GofR::writeRdf() { } else { - + sprintf(painCave.errMsg, "GofR: unable to open %s\n", outputFilename_.c_str()); + painCave.isFatal = 1; + simError(); } rdfStream.close();