FillStrip.cxx

Go to the documentation of this file.

// $Id: FillStrip.cxx,v 1.17 2008/01/31 22:18:17 rustem Exp $

// ROOT
#include "TClonesArray.h"

// MINOS
#include "CandNtupleSR/NtpSREvent.h"
#include "CandNtupleSR/NtpSRShower.h"
#include "CandNtupleSR/NtpSRSlice.h"
#include "CandNtupleSR/NtpSRStrip.h"
#include "CandNtupleSR/NtpSRTrack.h"
#include "MessageService/MsgService.h"
#include "MinosObjectMap/MomNavigator.h"
#include "Plex/PlexStripEndId.h"
#include "StandardNtuple/NtpStRecord.h"
#include "Validity/VldContext.h"

// Local
#include "PhysicsNtuple/Default.h"
#include "PhysicsNtuple/Factory.h"

// Local
#include "FillStrip.h"

CVSID("$Id: FillStrip.cxx,v 1.17 2008/01/31 22:18:17 rustem Exp $");

REGISTER_ANP_OBJECT(AlgStore,FillStrip)

using namespace std;

//---------------------------------------------------------------------------------------------
Anp::FillStrip::FillStrip()
{
}

//---------------------------------------------------------------------------------------------
Anp::FillStrip::~FillStrip()
{
}

//---------------------------------------------------------------------------------------------
bool Anp::FillStrip::Run(Record &record, TObject *ptr)
{
   //
   // Fill information for all reconstructed strips
   //

   NtpStRecord *ntprec = dynamic_cast<NtpStRecord *>(ptr);
   if(!ntprec)
   {
      const MomNavigator *mom = dynamic_cast<const MomNavigator *> (ptr);
      if(mom)
      {
         ntprec = dynamic_cast<NtpStRecord *>(mom -> GetFragment("NtpStRecord")); 
      }
      else
      {
         MSG("FillAlg", Msg::kError) << "Failed to find MomNavigator pointer" << endl;
         return false;
      }      
   }

   if(!ntprec)
   {
      MSG("FillAlg", Msg::kError) << "Failed to get NtpStRecord pointer" << endl;
      return false;
   }

   const TClonesArray *strip_array = ntprec -> stp;
   if(!strip_array)
   {
      MSG("FillAlg", Msg::kWarning) << "NtpStRecord does not have valid strip array." << endl;
      return false;
   }
      
   const int nstrip = strip_array -> GetEntries();
   if(nstrip != strip_array -> GetEntriesFast())
   {
      MSG("FillAlg", Msg::kWarning) << "NtpSRStrip TClonesArray.has holes" << endl;
   }

   std::vector<Anp::Strip> &strip_vec = record.GetStrip();

   // fill vector with empty strips
   strip_vec.clear();
   strip_vec.insert(strip_vec.begin(), nstrip, Anp::Strip());

   const Detector::Detector_t detector =  ntprec -> GetHeader().GetVldContext().GetDetector();

   for(int i_strip = 0; i_strip != nstrip; ++i_strip)
   {            
      const NtpSRStrip *ntp_strip  = dynamic_cast<const NtpSRStrip *>(strip_array -> At(i_strip));
      if(!ntp_strip)
      {
         MSG("FillAlg", Msg::kWarning) << "Could not get NtpSRStrip" << endl;
         continue;
      }      
      if(i_strip != ntp_strip -> index)
      {  
         MSG("FillAlg", Msg::kWarning) << "Wrong strip and TClonesArray index" << endl;
         continue;
      }

      Anp::Strip &strip = strip_vec[i_strip];
      
      strip.index     = ntp_strip -> index;
      strip.strip     = ntp_strip -> strip;
      strip.plane     = ntp_strip -> plane;
      strip.ndigit    = ntp_strip -> ndigit;
      strip.planeview = ntp_strip -> planeview;
      strip.tpos      = ntp_strip -> tpos;
      strip.zpos      = ntp_strip -> z;      

      if(detector == Detector::kNear)
      {
         if(ntp_strip -> pmtindex0 > 0)
         {
            MSG("FillAlg", Msg::kWarning) << "Near detector has east strip readout " << endl;
            continue;
         }
         if(ntp_strip -> pmtindex1 < 0 || !(ntp_strip -> ph1.raw > 0.0))
         {
            MSG("FillAlg", Msg::kWarning) << "Near detector does not have west strip readout " << endl;
            continue;
         }

         strip.adc       = ntp_strip -> ph1.raw;
         strip.pe        = ntp_strip -> ph1.pe;
         strip.siglin    = ntp_strip -> ph1.siglin;
         strip.sigcor0   = 0.0;
         strip.sigcor1   = ntp_strip -> ph1.sigcor;
         strip.time0     = 0.0;
         strip.time1     = ntp_strip -> time1;
         strip.pmtindex0 = -1;
         strip.pmtindex1 = ntp_strip -> pmtindex1;
      }
      else if(detector == Detector::kFar)
      {
         if(!(ntp_strip -> pmtindex0 >= 0) && !(ntp_strip -> pmtindex1 >= 0))
         {
            MSG("FillAlg", Msg::kWarning) << "Far detector does not have valid strip readout " << endl
                                          << *ntp_strip;
            continue;
         }

         if(!(ntp_strip -> ph0.raw > 0.0) && !(ntp_strip -> ph1.raw > 0.0))
         {
            MSG("FillAlg", Msg::kWarning) << "Far detector does not have valid strip readout " << endl
                                          << *ntp_strip;
            continue;
         }

         strip.adc     = 0.0;
         strip.pe      = 0.0;
         strip.siglin  = 0.0;
         strip.sigcor0 = 0.0;
         strip.sigcor1 = 0.0;
         
         if(ntp_strip -> pmtindex0 > 0 || 
            (ntp_strip -> pmtindex0 == 0 && ntp_strip -> ph0.raw > 0.0))
         {
            strip.adc    += ntp_strip -> ph0.raw;
            strip.pe     += ntp_strip -> ph0.pe;
            strip.siglin += ntp_strip -> ph0.siglin;

            strip.sigcor0   = ntp_strip -> ph0.sigcor;
            strip.time0     = ntp_strip -> time0;
            strip.pmtindex0 = ntp_strip -> pmtindex0;
         }

         if(ntp_strip -> pmtindex1 > 0 ||
            (ntp_strip -> pmtindex1 == 0 && ntp_strip -> ph1.raw > 0.0))
         {
            strip.adc    += ntp_strip -> ph1.raw;
            strip.pe     += ntp_strip -> ph1.pe;
            strip.siglin += ntp_strip -> ph1.siglin;

            strip.sigcor1   = ntp_strip -> ph1.sigcor;
            strip.time1     = ntp_strip -> time1;
            strip.pmtindex1 = ntp_strip -> pmtindex1;
         }
      }
      else
      {
         MSG("FillAlg", Msg::kWarning) << "Unknown detector type" << endl;
      }

      const PlexStripEndId seid = GetSEID(strip, detector);

      if(seid.GetPlane() != strip.plane || seid.GetStrip() != strip.strip)
      {
         MSG("FillAlg", Msg::kError) << "Mismatch between Plex and Strip" << endl;
      }

      strip.seid = seid.GetEncoded();
   }

   GetEvtStrip(strip_vec, *ntprec);
   GetShwStrip(strip_vec, *ntprec);
   GetSlcStrip(strip_vec, *ntprec);
   GetTrkStrip(strip_vec, *ntprec);

   return true;
}

//-----------------------------------------------------------------------------
void Anp::FillStrip::Config(const Registry &)
{
   //
   // Configure self
   //
}

//---------------------------------------------------------------------------------------------
void Anp::FillStrip::GetEvtStrip(SVec &svec, const NtpStRecord &ntprec) const
{
   const TClonesArray *event_array = ntprec.evt;   
   if(!event_array)
   {
      MSG("FillAlg", Msg::kWarning) << "NtpStRecord does not have valid event array" << endl;
      return;
   }
   
   const short size = svec.size();
   const int nevent = event_array -> GetEntries();

   for(int i_event = 0; i_event != nevent; ++i_event)
   {      
      const NtpSREvent* ntpevt = dynamic_cast<const NtpSREvent *> (event_array -> At(i_event));
      if(!ntpevt)
      {
         MSG("FillAlg", Msg::kWarning) << "Failed to find NtpSREvent at "<< i_event << endl;
         continue; 
      }

      const short index = ntpevt -> index;
      if(index != i_event)
      {
         MSG("FillAlg", Msg::kError) << "NtpSREvent "<< index << " has mistmatched index" << endl;
         continue; 
      }

      if(!ntpevt -> stp)
      {
         MSG("FillAlg", Msg::kWarning) << "NtpSREvent "<< index << " doesn't have strip array" << endl;
         continue;
      }
      
      for(short i_strip = 0; i_strip < ntpevt -> nstrip; ++i_strip)
      {
         const short strip_index = ntpevt -> stp[i_strip];

         if(strip_index < 0 || strip_index >= size)
         {
            MSG("FillAlg", Msg::kWarning) << "Strip index is out of range: " << strip_index << endl;
            continue;
         }

         Strip &strip = svec[strip_index];       
         if(strip.index < 0 || strip.index != strip_index)
         {
            MSG("FillAlg", Msg::kWarning) << "Mismateched strip index" << endl;
         }
         else
         {
            strip.evt.push_back(index);
         }
      }
   }
}

//---------------------------------------------------------------------------------------------
void Anp::FillStrip::GetSlcStrip(SVec &svec, const NtpStRecord &ntprec) const
{
   const TClonesArray *slice_array = ntprec.slc;
   if(!slice_array)
   {
      MSG("FillAlg", Msg::kWarning) << "NtpStRecord does not have valid slice array." << endl;
      return;
   }

   const short size = svec.size();   
   const int nslice = slice_array -> GetEntries();

   for(int i_slice = 0; i_slice != nslice; ++i_slice)
   {      
      const NtpSRSlice* ntp_slice = dynamic_cast<const NtpSRSlice *> (slice_array -> At(i_slice));
      if(!ntp_slice)
      {
         MSG("FillAlg", Msg::kWarning) << "Failed to find NtpSRSlice at "<< i_slice << endl;
         continue; 
      }

      const short index = ntp_slice -> index;
      if(index != i_slice)
      {
         MSG("FillAlg", Msg::kError) << "NtpSRSlice "<< index << " has mistmatched index" << endl;
         continue; 
      }

      if(!ntp_slice -> stp)
      {
         MSG("FillAlg", Msg::kWarning) << "NtpSRSlice "<< index << " doesn't have strip array" << endl;
         continue;
      }
      
      for(short i_strip = 0; i_strip < ntp_slice -> nstrip; ++i_strip)
      {
         const short strip_index = ntp_slice -> stp[i_strip];

         if(strip_index < 0 || strip_index >= size)
         {
            MSG("FillAlg", Msg::kWarning) << "Strip index is out of range: " << strip_index << endl;
            continue;
         }

         Strip &strip = svec[strip_index];       
         if(strip.index < 0 || strip.index != strip_index)
         {
            MSG("FillAlg", Msg::kWarning) << "Mismateched strip index" << endl;
         }
         else
         {
            strip.slc.push_back(index);
         }
      }
   }
}

//---------------------------------------------------------------------------------------------
void Anp::FillStrip::GetShwStrip(SVec &svec, const NtpStRecord &ntprec) const
{
   const TClonesArray *shower_array = ntprec.shw;   
   if(!shower_array)
   {
      MSG("FillAlg", Msg::kWarning) << "NtpStRecord does not have valid shower array." << endl;
      return;
   }
   
   const short size = svec.size();
   const int nshower = shower_array -> GetEntries();

   for(int i_shower = 0; i_shower != nshower; ++i_shower)
   {      
      const NtpSRShower* ntpshw = dynamic_cast<const NtpSRShower *> (shower_array -> At(i_shower));
      if(!ntpshw)
      {
         MSG("FillAlg", Msg::kWarning) << "Failed to find NtpSRShower at "<< i_shower << endl;
         continue; 
      }

      const short index = ntpshw -> index;
      if(index != i_shower)
      {
         MSG("FillAlg", Msg::kError) << "NtpSRShower "<< index << " has mistmatched index" << endl;
         continue; 
      }

      if(!ntpshw -> stp)
      {
         MSG("FillAlg", Msg::kWarning) << "NtpSRShower "<< index << " doesn't have strip array" << endl;
         continue;
      }
      
      for(short i_strip = 0; i_strip < ntpshw -> nstrip; ++i_strip)
      {
         const short strip_index = ntpshw -> stp[i_strip];
         
         if(strip_index < 0 || strip_index >= size)
         {
            MSG("FillAlg", Msg::kWarning) << "Strip index is out of range: " << strip_index << endl;
            continue;
         }

         Strip &strip = svec[strip_index];       
         if(strip.index < 0 || strip.index != strip_index)
         {
            MSG("FillAlg", Msg::kWarning) << "Mismateched strip index" << endl;
         }
         else
         {
            strip.shw.push_back(index);
         }
      }
   }
}

//---------------------------------------------------------------------------------------------
void Anp::FillStrip::GetTrkStrip(SVec &svec, const NtpStRecord &ntprec) const
{
   const TClonesArray *track_array = ntprec.trk;   
   if(!track_array)
   {
      MSG("FillAlg", Msg::kWarning) << "NtpStRecord does not have valid track array." << endl;
      return;
   }
   
   const short size = svec.size();
   const int ntrack = track_array -> GetEntries();

   for(int i_track = 0; i_track != ntrack; ++i_track)
   {      
      const NtpSRTrack* ntptrk = dynamic_cast<const NtpSRTrack *> (track_array -> At(i_track));
      if(!ntptrk)
      {
         MSG("FillAlg", Msg::kWarning) << "Failed to find NtpSRTrack at "<< i_track << endl;
         continue; 
      }

      const short index = ntptrk -> index;
      if(index != i_track)
      {
         MSG("FillAlg", Msg::kError) << "NtpSRTrack "<< index << " has mistmatched index" << endl;
         continue; 
      }

      if(!ntptrk -> stp)
      {
         MSG("FillAlg", Msg::kWarning) << "NtpSRTrack "<< index << " doesn't have strip array" << endl;
         continue;
      }
      
      for(short i_strip = 0; i_strip < ntptrk -> nstrip; ++i_strip)
      {
         const short strip_index = ntptrk -> stp[i_strip];
         
         if(strip_index < 0 || strip_index >= size)
         {
            MSG("FillAlg", Msg::kWarning) << "Strip index is out of range: " << strip_index << endl;
            continue;
         }

         Strip &strip = svec[strip_index];
         if(strip.index < 0 || strip.index != strip_index)
         {
            MSG("FillAlg", Msg::kError) << "Mismateched strip index" << endl;
            continue;
         }

         strip.trk.push_back(index);

         Strip::TrackInfo info;

         info.vpos = ntptrk -> stpv[i_strip];
         info.upos = ntptrk -> stpu[i_strip];
            
         info.sigmap_east = ntptrk -> stpph0sigmap[i_strip];
         info.sigmap_west = ntptrk -> stpph1sigmap[i_strip];
         info.mip_east = ntptrk -> stpph0mip[i_strip];
         info.mip_west = ntptrk -> stpph1mip[i_strip];

         if(!strip.trk_info.insert(Strip::TrackInfoMap::value_type(index, info)).second)
         {
            MSG("FillAlg", Msg::kWarning) << "Failed to insert TrackInfo structure" << endl;
         }
      }
   }
}

//-----------------------------------------------------------------------------
const PlexStripEndId Anp::FillStrip::GetSEID(const Strip &strip, const Detector::Detector_t det) const
{
   StripEnd::StripEnd_t end =  StripEnd::kWest;
   if(det == Detector::kFar)
   {
      if((strip.pmtindex0 >= 0 && strip.sigcor0 > 0.0) && 
         (strip.pmtindex1 >= 0 && strip.sigcor1 > 0.0))
      {
         end = StripEnd::kWhole;
      }
      else if(strip.pmtindex0 >= 0 && strip.sigcor0 > 0.0)
      {
         end = StripEnd::kEast;
      }
      else if(strip.pmtindex1 >= 0 && strip.sigcor1 > 0.0)
      {
         end = StripEnd::kWest;
      }
      else
      {
         cerr << "FillStrip::GetSEID - failed to determine strip end" << endl;
      }
   }
   
   return PlexStripEndId(det, strip.plane, strip.strip, end);
}

---