LightInjector Class Reference

#include <LightInjector.h>

Inheritance diagram for LightInjector:
JobCInputModule JobCModule

List of all members.

Public Member Functions

 LightInjector ()
 ~LightInjector ()
void BeginJob ()
const RegistryDefaultConfig () const
void Config (const Registry &r)
Int_t GetCurrentRun () const
Int_t GetCurrentSnarl () const
JobCResult Get ()
JobCResult Next (int n=1)
JobCResult Prev (int n=1)
JobCResult GoTo (int run, int snarl, int searchDir=0)
JobCResult Ana (const MomNavigator *mom)
JobCResult GoTo (const VldContext &)

Private Attributes

Int_t fDetector
Double_t fRate
Int_t fDate
Int_t fTime
Int_t fPlane
Int_t fStrip
Int_t fEnd
Double_t fPulseHeight
Double_t fPulseWidth
VldTimeStamp fPulseTime
VldContext fContext
PlexStripEndId fPulseSeid
int fCurrentRun
int fCurrentSnarl

Static Private Attributes

static PlexPixelSpotId fPsid
static TFile * fFile = 0
static TH1 * fHisto [65]

Detailed Description

Definition at line 17 of file LightInjector.h.


Constructor & Destructor Documentation

LightInjector::LightInjector (  ) 

Definition at line 51 of file LightInjector.cxx.

00051                              : 
00052   fCurrentRun(1000),
00053   fCurrentSnarl(0)
00054 { }

LightInjector::~LightInjector (  ) 

Definition at line 58 of file LightInjector.cxx.

00059 { 
00060 }


Member Function Documentation

JobCResult LightInjector::Ana ( const MomNavigator mom  )  [virtual]

Implement this for read only access to the MomNavigator

Reimplemented from JobCModule.

Definition at line 244 of file LightInjector.cxx.

References VldContext::AsString(), fFile, fHisto, Form(), fPsid, MomNavigator::FragmentIter(), RawDigit::GetADC(), RawDigit::GetChannel(), RawDigitDataBlock::GetDatumIter(), RawChannelId::GetElecType(), PlexMuxBoxId::GetElecType(), RawRecord::GetRawBlockIter(), RawRecord::GetRawHeader(), RecMinosHdr::GetVldContext(), RawChannelId::IsSameChannel(), Msg::kDebug, JobCResult::kPassed, ElecType::kQIE, Msg::kVerbose, MSG, pixels, and PlexPixelSpotId::SetPixel().

00245 {
00246   // Makes histograms of all the pixels
00247   // in the pulsed PMT.
00248   MSG("DetSim",Msg::kVerbose) << "LI Ana:" << endl;
00249   
00250   if(fFile==0) {
00251     TFile* save = gFile;
00252     fFile = new TFile("mcli.root","RECREATE");
00253     fFile->cd();
00254     for(int c=0;c<65;c++) {
00255       fHisto[c] = new TH1F(Form("h%02d",c),
00256                            Form("Pixel %02d",c),
00257                            20500,-500,20000 );
00258     }
00259     if(save) save->cd();
00260   }
00261  
00262 
00263   int ndigits_tot = 0;
00264   int ndigits = 0;
00265 
00266   // Get the RawDigitDataBlock.
00267   TObject* tobj;
00268   TIter    fragiter = mom->FragmentIter();
00269   while( ( tobj = fragiter.Next() ) ) {
00270     RawRecord* rawRec = dynamic_cast<RawRecord*>(tobj);
00271     if(rawRec) {
00272       TIter itr = rawRec->GetRawBlockIter();
00273       RawDataBlock* rawBlk;
00274       while ( ( rawBlk = dynamic_cast<RawDataBlock*>(itr.Next()) ) ) {
00275         RawDigitDataBlock* rddb = dynamic_cast<RawDigitDataBlock*>(rawBlk);
00276         if(rddb) {
00277           
00278           cout << "Ana context: " << rawRec->GetRawHeader()->GetVldContext().AsString() << endl;
00279           PlexHandle plex(rawRec->GetRawHeader()->GetVldContext());
00280 
00281           int npixels = 16;
00282           if(fPsid.GetElecType()==ElecType::kQIE) npixels = 64;
00283           
00284           PlexPixelSpotId pixels[65];
00285           RawChannelId    channels[65];
00286           Double_t        values[65];
00287           for(int i=0;i<npixels;i++){
00288             pixels[i] = fPsid;
00289             pixels[i].SetPixel(i);
00290             channels[i] = plex.GetRawChannelId(pixels[i]);
00291             values[i]=0;
00292           }
00293 
00294           TIter iter = rddb->GetDatumIter();
00295           while( (tobj=iter.Next() ) ) {
00296             
00297             // Process a digit.     
00298             RawDigit*        rawDigit = dynamic_cast<RawDigit *>(tobj);
00299             
00300             if(rawDigit) {
00301               // Try to add it.
00302               // We do it this way to 
00303               // make sure we get all buckets.
00304               ndigits_tot++;
00305               for(int i=0;i<npixels;i++) {
00306                 if(rawDigit->GetChannel().IsSameChannel(channels[i])) {
00307                   ndigits++;
00308                   values[i] += rawDigit->GetADC();
00309                   if(rawDigit->GetChannel().GetElecType()==ElecType::kQIE)
00310                     values[i] -= 50; // Remove digital offset for each bucket.
00311                 }
00312               }    
00313             }
00314 
00315           }
00316 
00317 
00318           MSG("DetSim",Msg::kDebug) << "Got " << ndigits << " digits; total in event = " << ndigits_tot << endl;
00319           // Fill histograms.
00320           for(int i=0;i<npixels;i++) {
00321             fHisto[i]->Fill(values[i]);
00322           }
00323           
00324 
00325 
00326         }       
00327       }
00328     }
00329   } // It's so ugly!
00330   
00331 
00332   return JobCResult::kPassed; // kNoDecision, kFailed, etc.
00333 }

void LightInjector::BeginJob ( void   )  [virtual]

Implement for notification of begin of job

Reimplemented from JobCModule.

Definition at line 64 of file LightInjector.cxx.

00065 {
00066 }

void LightInjector::Config ( const Registry r  )  [virtual]

Return the actual configuration. If your module directly pulls its configuration from the fConfig Registry, you don't need to override this. Override if you have local config variables.

Reimplemented from JobCModule.

Definition at line 101 of file LightInjector.cxx.

References Detector::CharToEnum(), StripEnd::CharToEnum(), fContext, fDate, fDetector, fEnd, fPlane, fPulseHeight, fPulseSeid, fPulseTime, fPulseWidth, fRate, fStrip, fTime, Registry::Get(), and SimFlag::kReroot.

00102 {
00103 //======================================================================
00104 // Configure the module given the Registry r
00105 //======================================================================
00106   int    tmpi;
00107   double tmpd;
00108 
00109   if (r.Get("detectorType",tmpi)) { fDetector = tmpi; }
00110   if (r.Get("rate",tmpd)) { fRate = tmpd; }
00111   if (r.Get("date",tmpi)) { fDate = tmpi; }
00112   if (r.Get("time",tmpi)) { fTime = tmpi; }
00113   if (r.Get("plane",tmpi)) { fPlane = tmpi; }
00114   if (r.Get("strip",tmpi)) { fStrip = tmpi; }
00115   if (r.Get("end",tmpi)) { fEnd = tmpi; }
00116   if (r.Get("pulseHeight",tmpd)) { fPulseHeight = tmpd; }
00117   if (r.Get("pulseWidth",tmpd)) { fPulseWidth = tmpd; }
00118 
00119   fPulseTime = VldTimeStamp(fDate,fTime,0),
00120   fContext =   VldContext(Detector::CharToEnum((char)fDetector),
00121                           SimFlag::kReroot,
00122                           fPulseTime),
00123   fPulseSeid = PlexStripEndId(Detector::CharToEnum((char)fDetector),
00124                           fPlane,
00125                           fStrip,
00126                           StripEnd::CharToEnum((char)fEnd));
00127 }

const Registry & LightInjector::DefaultConfig ( void   )  const [virtual]

Get the default configuration registry. This should normally be overridden. One useful idiom is to implement it like:

const Registry& MyModule::DefaultConfig() const { static Registry cfg; // never is destroyed if (cfg.Size()) return cfg; // already filled it // set defaults: cfg.Set("TheAnswer",42); cfg.Set("Units","unknown"); return cfg; }

Reimplemented from JobCModule.

Definition at line 71 of file LightInjector.cxx.

References JobCModule::GetName(), Registry::LockValues(), Registry::Set(), and Registry::UnLockValues().

00072 {
00073 //======================================================================
00074 // Get the default configuration for this module
00075 //======================================================================
00076   static Registry r; // Default configuration for module
00077 
00078   // Set name of config
00079   std::string name = this->GetName();
00080   name += ".config.default";
00081   r.SetName(name.c_str());
00082 
00083   // Set values in configuration
00084   r.UnLockValues();
00085   r.Set("detectorType", 2);
00086   r.Set("rate",         1.0);
00087   r.Set("date",         20030101);
00088   r.Set("time",         000000);
00089   r.Set("plane",        2);
00090   r.Set("strip",        10);
00091   r.Set("end",          2);
00092   r.Set("pulseHeight",  10.0);
00093   r.Set("pulseWidth",   0.0);
00094   r.LockValues();
00095 
00096   return r;
00097 }

JobCResult LightInjector::Get (  )  [virtual]

Reimplemented from JobCInputModule.

Definition at line 131 of file LightInjector.cxx.

References VldTimeStamp::Add(), RecDataRecord< T >::AdoptTemporary(), VldContext::AsString(), PlexPixelSpotId::AsString(), PlexStripEndId::AsString(), Detector::CharToEnum(), fContext, fCurrentRun, fCurrentSnarl, fDate, fDetector, fPsid, fPulseHeight, fPulseSeid, fPulseTime, fPulseWidth, fRate, fTime, JobCInputModule::GetMom(), PlexHandle::GetPixelSpotId(), RecRecordImp< T >::GetTempTags(), gSystem(), hostname, PlexPixelSpotId::IsValid(), JobCResult::kAOK, Msg::kError, JobCResult::kError, DigiSignal::kGenuine, Msg::kInfo, SimFlag::kMC, MSG, Munits::ns, Registry::Set(), and SimSnarlHeader::TrimCodename().

Referenced by GoTo(), Next(), and Prev().

00132 {
00133   MomNavigator* mom = this->GetMom();
00134   assert(mom);
00135   mom -> Clear(); // Moving on so clear contents of Mom
00136 
00137 
00138   // Make a SimSnarl.
00139   Int_t   trigbits = 0;
00140   Short_t subrun = 0;
00141   Short_t runtype = 0;
00142   Int_t   errcode = 0;
00143   Int_t   timeframe = -1;
00144 
00145   // Increment the time by the pulse rate.
00146   int secs;
00147   double frac = frexp(fRate*fCurrentSnarl,&secs);
00148   int nsecs = int(frac*1e9);
00149   VldTimeStamp offset(secs,nsecs);
00150   fPulseTime = VldTimeStamp(fDate,fTime,0);
00151   fPulseTime.Add(offset);
00152 
00153   // Reset the context.
00154   fContext = VldContext(Detector::CharToEnum((char)fDetector),
00155                         SimFlag::kMC,
00156                         fPulseTime);
00157 
00158   // Get a PlexHandle to find the psid.
00159   PlexHandle plex(fContext);
00160   fPsid = plex.GetPixelSpotId(fPulseSeid);
00161   if(!fPsid.IsValid()) {
00162     MSG("DetSim",Msg::kError) << "Invalid psid " << fPsid.AsString()
00163                               << " for strip " << fPulseSeid.AsString()
00164                               << " at context " << fContext.AsString()
00165                               << endl;
00166     return JobCResult::kError;
00167   }
00168 
00169   VldTimeStamp now;
00170   std::string  codename = SimSnarlHeader::TrimCodename("$Name:  $");
00171   std::string  hostname(gSystem->HostName());
00172 
00173   SimSnarlHeader simheader(fContext,fCurrentRun,subrun,runtype,
00174                            errcode,fCurrentSnarl,trigbits,timeframe,-1,
00175                            now,codename,hostname);
00176                                                 
00177   SimSnarlRecord* simsnarl = new SimSnarlRecord(simheader);
00178 
00179   simsnarl->GetTempTags().Set("stream","SimSnarl");  // no idea....
00180     
00181   // Give it to mom.
00182   mom -> AdoptFragment(simsnarl);
00183 
00184   // Add the DigiPE list to the simsnarl.
00185   TObjArray* PeList = new TObjArray(0);
00186   PeList->SetName("DigiListPe");
00187   PeList->SetOwner(true);
00188 
00189   static double tot_pulses = 0;
00190   static double tot_npe = 0;
00191   
00192   // Find out how many of them to make.
00193   gRandom->SetSeed(fCurrentSnarl);
00194   Int_t npe = gRandom->Poisson(fPulseHeight);
00195 
00196   // Make the photoelectron objects.
00197   for(int i=0;i<npe;i++) {
00198     // Time is 100ns +- pulse width.
00199     double time = gRandom->Gaus(100.0*Munits::ns, fPulseWidth);
00200     
00201     DigiPE*  pe = new DigiPE(time, fPsid, DigiSignal::kGenuine);
00202     PeList->Add(pe);
00203   }
00204 
00205   tot_pulses +=1;
00206   tot_npe += npe;
00207   MSG("DetSim",Msg::kInfo) << "SnarlNo: " << fCurrentSnarl << "  Pulse: " << npe << "avg: " << tot_npe/tot_pulses << " pe into " << fPsid.AsString()  << endl;
00208 
00209   simsnarl->AdoptTemporary(PeList);
00210 
00211   return JobCResult::kAOK;
00212 }

Int_t LightInjector::GetCurrentRun (  )  const [inline, virtual]

Return the current run #. This should not be overridden. Use this in BeginRun().

Reimplemented from JobCInputModule.

Definition at line 33 of file LightInjector.h.

References fCurrentRun.

00033 { return fCurrentRun; };

Int_t LightInjector::GetCurrentSnarl (  )  const [inline, virtual]

Reimplemented from JobCInputModule.

Definition at line 34 of file LightInjector.h.

References fCurrentSnarl.

00034 { return fCurrentSnarl; };

JobCResult LightInjector::GoTo ( const VldContext  )  [inline, virtual]

Reimplemented from JobCInputModule.

Definition at line 43 of file LightInjector.h.

References Get().

00043 { return this->Get();};

JobCResult LightInjector::GoTo ( int  run,
int  snarl,
int  searchDir = 0 
) [virtual]

Reimplemented from JobCInputModule.

Definition at line 236 of file LightInjector.cxx.

References fCurrentSnarl, and Get().

00237 {
00238   fCurrentSnarl = snarl;
00239   return this->Get();
00240 }

JobCResult LightInjector::Next ( int  n = 1  )  [virtual]

Reimplemented from JobCInputModule.

Definition at line 216 of file LightInjector.cxx.

References fCurrentSnarl, and Get().

00217 {
00218 //======================================================================
00219 // Advance the position in the stream n record sets. Load the records
00220 // at the last position
00221 //======================================================================
00222   fCurrentSnarl+=n;
00223   return this->Get();
00224 }

JobCResult LightInjector::Prev ( int  n = 1  )  [virtual]

Reimplemented from JobCInputModule.

Definition at line 228 of file LightInjector.cxx.

References fCurrentSnarl, and Get().

00229 {
00230   if(fCurrentSnarl>0) fCurrentSnarl-=n;
00231   return this->Get();
00232 }


Member Data Documentation

Definition at line 71 of file LightInjector.h.

Referenced by Config(), and Get().

Definition at line 79 of file LightInjector.h.

Referenced by Get(), and GetCurrentRun().

Definition at line 80 of file LightInjector.h.

Referenced by Get(), GetCurrentSnarl(), GoTo(), Next(), and Prev().

Int_t LightInjector::fDate [private]

Definition at line 61 of file LightInjector.h.

Referenced by Config(), and Get().

Int_t LightInjector::fDetector [private]

Definition at line 43 of file LightInjector.h.

Referenced by Config(), and Get().

Int_t LightInjector::fEnd [private]

Definition at line 65 of file LightInjector.h.

Referenced by Config().

TFile * LightInjector::fFile = 0 [static, private]

Definition at line 76 of file LightInjector.h.

Referenced by Ana().

TH1 * LightInjector::fHisto [static, private]

Definition at line 77 of file LightInjector.h.

Referenced by Ana().

Int_t LightInjector::fPlane [private]

Definition at line 63 of file LightInjector.h.

Referenced by Config().

Definition at line 75 of file LightInjector.h.

Referenced by Ana(), and Get().

Double_t LightInjector::fPulseHeight [private]

Definition at line 66 of file LightInjector.h.

Referenced by Config(), and Get().

Definition at line 72 of file LightInjector.h.

Referenced by Config(), and Get().

Definition at line 70 of file LightInjector.h.

Referenced by Config(), and Get().

Double_t LightInjector::fPulseWidth [private]

Definition at line 67 of file LightInjector.h.

Referenced by Config(), and Get().

Double_t LightInjector::fRate [private]

Definition at line 60 of file LightInjector.h.

Referenced by Config(), and Get().

Int_t LightInjector::fStrip [private]

Definition at line 64 of file LightInjector.h.

Referenced by Config().

Int_t LightInjector::fTime [private]

Definition at line 62 of file LightInjector.h.

Referenced by Config(), and Get().


The documentation for this class was generated from the following files:

Generated on 2 Nov 2017 for loon by  doxygen 1.6.1