Merge pull request #174 from JeffersonLab/all_threads_get_their_own_b…

…eam_generator_rtj

All threads get their own beam generator rtj

src/AdaptiveSampler.cc

@@ -134,6 +134,7 @@
 //    of AdaptiveSampler.
 
 #include <assert.h>
+#include <cmath>
 #include <fstream>
 #include <sstream>
 #include <iostream>
@@ -834,18 +835,20 @@ int AdaptiveSampler::Cell::serialize(std::ofstream &ofs, bool optimized)
    if (nhit != 0)
       ofs << "nhit=" << nhit << std::endl;
    if (sum_wI != 0)
-      ofs << "sum_wI=" << sum_wI << std::endl;
+      ofs << "sum_wI=" << (std::isfinite(sum_wI)? sum_wI : 0) << std::endl;
    if (sum_wI2 != 0)
-      ofs << "sum_wI2=" << sum_wI2 << std::endl;
+      ofs << "sum_wI2=" << (std::isfinite(sum_wI2)? sum_wI2 : 0) << std::endl;
    if (sum_wI4 != 0)
-      ofs << "sum_wI4=" << sum_wI4 << std::endl;
+      ofs << "sum_wI4=" << (std::isfinite(sum_wI4)? sum_wI4 : 0) << std::endl;
    for (int i=0; i < 3*ndim; ++i) {
       if (sum_wI2d[i] != 0)
-         ofs << "sum_wI2d[" << i << "]=" << sum_wI2d[i] << std::endl;
+         ofs << "sum_wI2d[" << i << "]=" 
+             << (std::isfinite(sum_wI2d[i])? sum_wI2d[i] : 0) << std::endl;
    }
    for (int i=0; i < 3*ndim; ++i) {
       if (sum_wI4d[i] != 0)
-         ofs << "sum_wI4d[" << i << "]=" << sum_wI4d[i] << std::endl;
+         ofs << "sum_wI4d[" << i << "]=" 
+             << (std::isfinite(sum_wI4d[i])? sum_wI4d[i] : 0) << std::endl;
    }
    if (optimized)
       ofs << "subset=" << std::setprecision(20) << opt_subset << std::endl;

src/GlueXPrimaryGeneratorAction.cc

@@ -79,6 +79,7 @@ GlueXPrimaryGeneratorAction::GlueXPrimaryGeneratorAction()
    fParticleGun = new GlueXParticleGun();
 
    if (fSourceType == SOURCE_TYPE_HDDM) {
+      clone_photon_beam_generator();
       fPrimaryGenerator = new GlueXPrimaryGenerator(fHDDMistream);
       return;
    }
@@ -493,7 +494,7 @@ GlueXPrimaryGeneratorAction &GlueXPrimaryGeneratorAction::operator=(const
    else if (fSourceType == SOURCE_TYPE_PARTICLE_GUN) {
       fParticleGun->SetParticleDefinition(fGunParticle.partDef);
    }
-   if (fPhotonBeamGenerator) {
+   if (src.fPhotonBeamGenerator) {
       clone_photon_beam_generator();
    }
    fBeamvertex = src.fBeamvertex;

src/GlueXUserEventInformation.cc

@@ -93,13 +93,18 @@ long int *GlueXUserEventInformation::fStartingSeeds = 0;
 G4Mutex GlueXUserEventInformation::fMutex = G4MUTEX_INITIALIZER;
 
 GlueXUserEventInformation::GlueXUserEventInformation(hddm_s::HDDM *hddmevent)
- : fKeepEvent(true),
+ : fOutputRecord(0),
+   fKeepEvent(true),
    fNprimaries(0),
    fNvertices(0)
 {
+   fEventSequenceNo = HddmOutput::incrementEventNo();
    if (hddmevent == 0) {
+      int runNo = GetRunNo();
       fOutputRecord = new hddm_s::HDDM();
-      fOutputRecord->addPhysicsEvents();
+      hddm_s::PhysicsEventList pev = fOutputRecord->addPhysicsEvents();
+      pev(0).setEventNo(fEventSequenceNo);
+      pev(0).setRunNo(runNo);
    }
    else {
       fOutputRecord = hddmevent;
@@ -120,12 +125,6 @@ GlueXUserEventInformation::~GlueXUserEventInformation()
       }
       if (fKeepEvent) {
          hddm_s::PhysicsEventList pev = fOutputRecord->getPhysicsEvents();
-         int runno = HddmOutput::getRunNo();
-         if (runno > 0)
-            pev(0).setRunNo(runno);
-         if (pev(0).getEventNo() == 0) {
-            pev(0).setEventNo(HddmOutput::incrementEventNo());
-         }
          if (fWriteNoHitEvents || pev(0).getHitViews().size() > 0) {
             HddmOutput::WriteOutputHDDM(*fOutputRecord);
          }
@@ -452,19 +451,16 @@ void GlueXUserEventInformation::SetRandomSeeds()
 
    hddm_s::PhysicsEventList pev = fOutputRecord->getPhysicsEvents();
    hddm_s::ReactionList rea = pev(0).getReactions();
-   long int eventNo = GetEventSequenceNo();
    if (rea.size() == 0) {
       rea = pev(0).addReactions();
    }
    hddm_s::RandomList rnd = rea(0).getRandoms();
    if (rnd.size() > 0) {
       fEventSeeds[0] = rnd(0).getSeed1();
       fEventSeeds[1] = rnd(0).getSeed2();
-      if (pev(0).getEventNo() > 0)
-         eventNo = pev(0).getEventNo();
       G4Random::setTheSeeds(fEventSeeds);
 #if VERBOSE_RANDOMS
-      G4cout << "New event " << eventNo << " with starting seeds " 
+      G4cout << "New event " << pev(0).getEventNo() << " with starting seeds " 
              << fEventSeeds[0] << ", " << fEventSeeds[1] << G4endl;
 #endif
    }
@@ -486,7 +482,7 @@ void GlueXUserEventInformation::SetRandomSeeds()
       rnd(0).setSeed3(709975946 + pev(0).getEventNo());
       rnd(0).setSeed4(912931182 + pev(0).getEventNo());
 #if VERBOSE_RANDOMS
-      G4cout << "New event " << eventNo << " with starting seeds " 
+      G4cout << "New event " << pev(0).getEventNo() << " with starting seeds " 
              << fEventSeeds[0] << ", " << fEventSeeds[1] << G4endl;
 #endif
    }
@@ -495,7 +491,8 @@ void GlueXUserEventInformation::SetRandomSeeds()
 int GlueXUserEventInformation::GetRunNo()
 {
    if (fOutputRecord && fOutputRecord->getPhysicsEvents().size() > 0) {
-      return fOutputRecord->getPhysicsEvent().getRunNo();
+      hddm_s::PhysicsEventList pev = fOutputRecord->getPhysicsEvents();
+      return pev(0).getRunNo();
    }
    G4RunManager *runmgr = G4RunManager::GetRunManager();
    if (runmgr != 0 && runmgr->GetCurrentRun() != 0) {
@@ -504,13 +501,13 @@ int GlueXUserEventInformation::GetRunNo()
    return 0;
 }
 
-long int GlueXUserEventInformation::GetEventSequenceNo()
+long int GlueXUserEventInformation::GetEventNo()
 {
-   G4RunManager *runmgr = G4RunManager::GetRunManager();
-   if (runmgr != 0 && runmgr->GetCurrentEvent() != 0) {
-      return runmgr->GetCurrentEvent()->GetEventID();
+   if (fOutputRecord && fOutputRecord->getPhysicsEvents().size() > 0) {
+      hddm_s::PhysicsEventList pev = fOutputRecord->getPhysicsEvents();
+      return pev(0).getEventNo();
    }
-   return 0;
+   return fEventSequenceNo;
 }
 
 double GlueXUserEventInformation::GetBeamPhotonEnergy()

src/GlueXUserEventInformation.hh

@@ -45,7 +45,7 @@ class GlueXUserEventInformation: public G4VUserEventInformation
    void AddMCtrajectoryPoint(const G4Step &step, int save_option);
 
    int GetRunNo();
-   long int GetEventSequenceNo();
+   long int GetEventNo();
    double GetBeamPhotonEnergy();
    int GetGlueXTrackID(int g4ID);
    int GetGlueXTrackID(const G4Track *track);
@@ -93,6 +93,7 @@ class GlueXUserEventInformation: public G4VUserEventInformation
    std::map<std::string, std::fstream*> fDlogfile;
    std::map<std::string, int> fDlogreading;
    long int fEventSeeds[2];
+   long int fEventSequenceNo;
 
    static G4Mutex fMutex;
 };

src/HddmOutput.cc

@@ -54,6 +54,7 @@ HddmOutput& HddmOutput::operator=(HddmOutput &src)
 
 void HddmOutput::setRunNo(int runno)
 {
+   G4AutoLock barrier(&fMutex);
    fRunNo = runno;
 }
 

src/HddmOutput.hh

@@ -57,6 +57,7 @@ inline int HddmOutput::getEventNo()
 
 inline void HddmOutput::setEventNo(int eventno)
 {
+   G4AutoLock barrier(&fMutex);
    fEventNo = eventno;
 }