Refactor sim.py into a reusable Simulation class that is called by the

main function().  Also cleanup more event data structure names and add
an nphoton value that is preserved even if you prune off all the
actual photon vertices.

2 files changed, 30 insertions, 28 deletions

diff --git a/fileio/root.C b/fileio/root.C
index 8f7d0a9..7052a71 100644
--- a/fileio/root.C
+++ b/fileio/root.C
@@ -4,10 +4,10 @@
 #include <string>
 
 struct Photon {
-  double t;
-  TVector3 pos;
-  TVector3 dir;
-  TVector3 pol;
+  double time;
+  TVector3 position;
+  TVector3 direction;
+  TVector3 polarization;
   double wavelength; // nm
   unsigned int history;
   int last_hit_triangle;
@@ -15,19 +15,20 @@ struct Photon {
 
 struct Track {
   std::string particle;
-  double t;
-  TVector3 pos;
-  TVector3 dir;
+  double time;
+  TVector3 position;
+  TVector3 direction;
   double start_time;
   double total_energy;
 };
 
 struct MC {
   std::string particle;
-  TVector3 gen_pos;
-  TVector3 gen_dir;
-  double gen_total_e;
+  TVector3 gen_position;
+  TVector3 gen_direction;
+  double gen_total_energy;
 
+  int nphoton;
   
   std::vector<Track> subtrack;
   std::vector<Photon> photon_start;
@@ -36,10 +37,10 @@ struct MC {
 };
 
 struct Channel {
-  Channel() : channel_id(-1), t(-9999.0), q(-9999.0) { };
+  Channel() : channel_id(-1), time(-9999.0), charge(-9999.0) { };
   int channel_id;
-  double t;
-  double q;
+  double time;
+  double charge;
   unsigned int mc_history;
 };
   
@@ -65,8 +66,8 @@ struct Event {
 
       if (channel_id < nentries) {
 	hit[channel_id] = 1;
-	time[channel_id] = channel[i].t;
-	charge[channel_id] = channel[i].q;
+	time[channel_id] = channel[i].time;
+	charge[channel_id] = channel[i].charge;
       }
     }
   }
@@ -83,10 +84,10 @@ void fill_photons(Event *ev, bool start,
   
   for (unsigned int i=0; i < nphotons; i++) {
     Photon &photon = photons[i];
-    photon.t = t0[i];
-    photon.pos.SetXYZ(pos[3*i], pos[3*i + 1], pos[3*i + 2]);
-    photon.dir.SetXYZ(dir[3*i], dir[3*i + 1], dir[3*i + 2]);
-    photon.pol.SetXYZ(pol[3*i], pol[3*i + 1], pol[3*i + 2]);
+    photon.time = t0[i];
+    photon.position.SetXYZ(pos[3*i], pos[3*i + 1], pos[3*i + 2]);
+    photon.direction.SetXYZ(dir[3*i], dir[3*i + 1], dir[3*i + 2]);
+    photon.polarization.SetXYZ(pol[3*i], pol[3*i + 1], pol[3*i + 2]);
     photon.wavelength = wavelength[i];
     if (histories)
       photon.history = histories[i];
@@ -101,18 +102,18 @@ void fill_photons(Event *ev, bool start,
 }
 		  
 
-void fill_hits(Event *ev, unsigned int nchannels, float *t,
-	       float *q, unsigned int *history)
+void fill_hits(Event *ev, unsigned int nchannels, float *time,
+	       float *charge, unsigned int *history)
 {
   ev->channel.resize(0);
   ev->nhit = 0;
   Channel ch;
   for (unsigned int i=0; i < nchannels; i++) {
-    if (t[i] < 1e8) {
+    if (time[i] < 1e8) {
       ev->nhit++;
       ch.channel_id = i;
-      ch.t = t[i];
-      ch.q = q[i];
+      ch.time = time[i];
+      ch.charge = charge[i];
       ch.mc_history = history[i];
       ev->channel.push_back(ch);
     }
diff --git a/fileio/root.py b/fileio/root.py
index 5bebda0..77c23ff 100644
--- a/fileio/root.py
+++ b/fileio/root.py
@@ -20,9 +20,10 @@ class RootWriter(object):
         self.ev.event_id = pyev.event_id
 
         self.ev.mc.particle = pyev.particle_name
-        self.ev.mc.gen_pos.SetXYZ(*pyev.gen_position)
-        self.ev.mc.gen_dir.SetXYZ(*pyev.gen_direction)
+        self.ev.mc.gen_position.SetXYZ(*pyev.gen_position)
+        self.ev.mc.gen_direction.SetXYZ(*pyev.gen_direction)
         self.ev.mc.gen_total_energy = pyev.gen_total_energy
+        self.ev.mc.nphoton = pyev.nphoton
 
         if pyev.photon_start is not None:
             photons = pyev.photon_start
@@ -46,8 +47,8 @@ class RootWriter(object):
             self.ev.mc.subtrack.resize(len(pyev.subtracks))
             for i, subtrack in enumerate(pyev.subtracks):
                 self.ev.mc.subtrack[i].name = subtrack.particle_name
-                self.ev.mc.subtrack[i].pos.SetXYZ(*subtrack.position)
-                self.ev.mc.subtrack[i].dir.SetXYZ(*subtrack.direction)
+                self.ev.mc.subtrack[i].position.SetXYZ(*subtrack.position)
+                self.ev.mc.subtrack[i].direction.SetXYZ(*subtrack.direction)
                 self.ev.mc.subtrack[i].start_time = subtrack.start_time
                 self.ev.mc.subtrack[i].total_energy = subtrack.total_energy