initial draft of simIDE interface

Author: Jacob McLaughlin

sbndcode/BlipRecoSBND/Alg/BlipRecoAlg.cc

@@ -357,18 +357,31 @@ namespace blip {
     if (evt.getByLabel(fGeantProducer,pHandle))
       art::fill_ptr_vector(plist, pHandle);
 
-    // -- SimEnergyDeposits
-    art::Handle<std::vector<sim::IDE> > sedHandle;
-    std::vector<art::Ptr<sim::IDE> > sedlist;
-    if (evt.getByLabel(fSimDepProducer,sedHandle)){
-      art::fill_ptr_vector(sedlist, sedHandle);
-    }
-
+    // -- SimEnergyDeposits (usually dropped in reco
+    //art::Handle<std::vector<sim::EnergyDeposit> > sedHandle;
+    std::vector<sim::IDE > sedlist;
+    //if (evt.getByLabel(fSimDepProducer,sedHandle)){
+    //  art::fill_ptr_vector(sedlist, sedHandle);
+    // }
+    std::cout << " \n\n\n\n\n\n\n\n\n Look here \n\n\n\n\n " << std::endl;
     // -- SimChannels (usually dropped in reco)
     art::Handle<std::vector<sim::SimChannel> > simchanHandle;
     std::vector<art::Ptr<sim::SimChannel> > simchanlist;
-    if (evt.getByLabel(fSimChanProducer,simchanHandle)) 
+    std::cout << fSimChanProducer << " Look at my sim channel handle here " <<std::endl;
+    if (evt.getByLabel(fSimChanProducer,simchanHandle))
+      { 
       art::fill_ptr_vector(simchanlist, simchanHandle);
+      //Loop over channels to get full sedlist
+      for(int chIndex=0; chIndex<int(simchanlist.size()); chIndex++)
+	{
+	  std::vector< sim::IDE > TempChIDE = (*simchanlist[chIndex]).TrackIDsAndEnergies(0, 999999999);
+	  for(int ideIndex=0; ideIndex<int(TempChIDE.size()); ideIndex++)
+	    {
+	      //art::fill_ptr_vector(sedlist, simchanHandle.TrackIDsAndEnergies(0, 99999999));
+	      sedlist.push_back( TempChIDE[ideIndex] ); //may need to add a &
+	    }
+	}
+      }
 
     // -- hits (from input module, usually track-masked subset of gaushit)
     art::Handle< std::vector<recob::Hit> > hitHandle;

sbndcode/BlipRecoSBND/Utils/BlipUtils.cc

@@ -61,12 +61,12 @@ namespace BlipUtils {
     pinfo.depEnergy     = 0;
     pinfo.depElectrons  = 0;
     for(auto& sed : sedvec ) {
-      if( sed->trackID == part.TrackId() ) {
-        pinfo.depEnergy     += sed->energy;
-        pinfo.depElectrons  += sed->numElectrons;
+      if( -1*sed.trackID == part.TrackId() || sed.trackID == part.TrackId() ) {
+        pinfo.depEnergy     += sed.energy;
+        pinfo.depElectrons  += sed.numElectrons;
       }
     }
-    
+    std::cout << pinfo.depEnergy << "  " << pinfo.trackId << std::endl;
     return;
 
   }
@@ -86,7 +86,7 @@ namespace BlipUtils {
       //std::cout<<"Making true blip for "<<part.TrackId()<<" (PDG "<<part.PdgCode()<<", which deposited "<<pinfo[i].depEnergy<<"\n";
 
       // If this is a photon or neutron, don't even bother!
-      if( part.PdgCode() == 2112 || part.PdgCode() == 22 ) continue;
+      //if( part.PdgCode() == 2112 || part.PdgCode() == 22 ) continue;
 
       // If this is an electron that came from another electron, it would 
       // have already been grouped as part of the contiguous "blip" previously.
@@ -95,7 +95,9 @@ namespace BlipUtils {
 
       // Create the new blip
       blip::TrueBlip tb;
+      std::cout << " making a new blip " << std::endl;
       GrowTrueBlip(pinfo[i],tb);
+      std::cout << tb.Energy << std::endl;
       if( !tb.Energy ) continue;  
 
       // We want to loop through any contiguous electrons (produced
@@ -112,6 +114,7 @@ namespace BlipUtils {
 
       // Final check -- ensure there was non-negligible number 
       // of deposted ionization electrons
+      std::cout << tb.DepElectrons << " total electrons " << std::endl;
       if( tb.DepElectrons < 20 ) continue;
 
       // Calculate TPC-specific quantities
@@ -126,6 +129,7 @@ namespace BlipUtils {
       tb.DriftTime = tick_calc*clockData.TPCClock().TickPeriod() + clockData.TriggerOffsetTPC();
 
       tb.ID = trueblips.size();
+      std::cout << "going to push this into holder " << std::endl;
       trueblips.push_back(tb);
 
     }
@@ -139,7 +143,7 @@ namespace BlipUtils {
     simb::MCParticle& part = pinfo.particle;
 
     // Skip neutrons, photons
-    if( part.PdgCode() == 2112 || part.PdgCode() == 22 ) return;
+    // if( part.PdgCode() == 2112 || part.PdgCode() == 22 ) return;
 
     // Check that path length isn't zero
     if( !pinfo.pathLength ) return;

sbndcode/BlipRecoSBND/Utils/BlipUtils.h

@@ -38,7 +38,7 @@
 #include "TH1D.h"
 
 
-typedef std::vector<art::Ptr<sim::IDE>> SEDVec_t;
+typedef std::vector<sim::IDE> SEDVec_t;
 
 geo::View_t kViews[3]={geo::kU, geo::kV, geo::kW};
 

sbndcode/BlipRecoSBND/blipreco_configs.fcl

@@ -9,7 +9,7 @@ sbnd_blipalg:
   TrkProducer:          "pandoraTrack"       #// input recob::Tracks to use for blip reconstruction
   GeantProducer:        "largeant"      #// input sim::MCParticles (getting true particle info)
   SimEDepProducer:      "ionandscint"   #// input sim::SimEnergyDeposits (getting energy/electrons deposited)
-  SimChanProducer:      "simdrift"      #// label for sim::SimChannels (getting drifted charge; optional)
+  SimChanProducer:      "simtpc2d:simpleSC:DetSim"      #// label for sim::SimChannels (getting drifted charge; optional)
   MaxHitTrkLength:      5.0             #// hits in trks > this length will be vetoed [cm]
   DoHitFiltering:       false           #// filter hits based on amp, width, RMS
   HitClustWireRange:    1               #// clustering wire width