Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

PhotonCalibratedPeComputer Class Reference

#include <PhotonCalibratedPeComputer.h>

Inheritance diagram for PhotonCalibratedPeComputer:

PhotonTransportModule List of all members.

Public Member Functions

 PhotonCalibratedPeComputer ()
virtual ~PhotonCalibratedPeComputer ()
virtual void Configure (const Registry &r)
virtual void Reset (const VldContext &cx)
virtual Bool_t ComputePhotons (const DigiScintHit *inHit, Double_t inBluePrescale, Double_t inWlsPrescale, Double_t inGreenPrescale, PhotonCount &outCount)

Private Member Functions

virtual Double_t ClipFraction (double x) const
 ClassDef (PhotonCalibratedPeComputer, 0)

Private Attributes

Double_t fOverallLightOutput
Double_t fBirksConstant
Double_t fGeVPerMip
Int_t fScintillatorClipping
Double_t fMaxPePerHit

Constructor & Destructor Documentation

PhotonCalibratedPeComputer::PhotonCalibratedPeComputer (   ) 
 

Definition at line 16 of file PhotonCalibratedPeComputer.cxx.

References PhotonConfiguration().

00017 {
00018   Configure(PhotonConfiguration());
00019 }

PhotonCalibratedPeComputer::~PhotonCalibratedPeComputer (   )  [virtual]
 

Definition at line 23 of file PhotonCalibratedPeComputer.cxx.

00024 {
00025  
00026 }

Member Function Documentation

PhotonCalibratedPeComputer::ClassDef (  PhotonCalibratedPeComputer  ,
0 
)  [private]
 

Double_t PhotonCalibratedPeComputer::ClipFraction (  double  x  )  const [private, virtual]
 

Definition at line 280 of file PhotonCalibratedPeComputer.cxx.

Referenced by ComputePhotons().

00281 {
00282   // Returns the integral of the double-lorentzian blue photon profile
00283   // from -inf to x.
00284   
00285   // Same parameters as in PhotonFastBlueModel.
00286 
00287   const double kLorentzWidth1 = 0.017853;
00288   const double kArea1 = 8.178052;
00289   const double kLorentzWidth2 = 0.001309;
00290   const double kArea2 = 0.864215;
00291 
00292   const double kArea = kArea1+kArea2;
00293 
00294   double frac = kArea1 * atan(x/kLorentzWidth1)
00295     + kArea2 * atan(x/kLorentzWidth2);
00296   frac = frac / (kArea * TMath::Pi());
00297   frac += 0.5;
00298 
00299   // Bonus extra fudge-factor:
00300   // We don't want this effect to be seen outside ~20 cm 
00301   // (where the effect is already down around 1%)
00302   // So, supress by an additional amount:
00303   if (x>0.2) frac = 1-(1-frac)*exp(-(x-0.2)*20.);
00304 
00305   return frac;
00306 }

Bool_t PhotonCalibratedPeComputer::ComputePhotons (  const DigiScintHit *  inHit,
Double_t  inBluePrescale,
Double_t  inWlsPrescale,
Double_t  inGreenPrescale,
PhotonCount &  outCount
)  [virtual]
 

Reimplemented from PhotonTransportModule.

Definition at line 59 of file PhotonCalibratedPeComputer.cxx.

References PlexStripEndId::AsString(), ClipFraction(), DigiScintHit::DE(), Calibrator::DecalAttenCorrected(), Calibrator::DecalDrift(), Calibrator::DecalMIP(), Calibrator::DecalStripToStrip(), DigiScintHit::DS(), fBirksConstant, fMaxPePerHit, fOverallLightOutput, UgliStripHandle::GetHalfLength(), PhotonCount::GetPe_Neg(), PhotonCount::GetPe_Pos(), Calibrator::GetPhotoElectrons(), UgliGeomHandle::GetStripHandle(), Calibrator::Instance(), UgliStripHandle::IsMirrored(), MSG, UgliStripHandle::PartialLength(), PlexStripEndId::SetEnd(), PhotonCount::SetPeNegPos(), DigiScintHit::StripEndId(), UgliStripHandle::WlsBypass(), DigiScintHit::X1(), and DigiScintHit::X2().

00064 {
00065   // Number of photons made:
00066   if (!inHit) return false;
00067   double dE = inHit->DE();
00068   double dx = inHit->DS();
00069   double dEdx = dE/dx;
00070 
00071   // Get position.
00072   double x = 0.5*(inHit->X1() + inHit->X2());
00073   
00074   double clippedfrac = 1.0;
00075 
00076   UgliGeomHandle ugli(fContext);
00077   UgliStripHandle ustrip=ugli.GetStripHandle(inHit->StripEndId());
00078 
00079   if (fScintillatorClipping) {
00080     // Clip the photon response by shaving off photons lost out the
00081     // ends of the scintillator either at the planed edges or the bore hole.
00082 
00083     // To do this, we cheat a bit: we assume that ALL strip ends play a role
00084     // in diminishing the return. 
00085    
00086     double distToEnd;
00087     
00088     distToEnd = ustrip.GetHalfLength()-x; // +ve end.
00089     clippedfrac *= ClipFraction( distToEnd ); // End1
00090 
00091     distToEnd = x+ustrip.GetHalfLength(); // -ve end.
00092     clippedfrac *= ClipFraction( distToEnd ); // End2
00093 
00094     if (ustrip.WlsBypass() > 0.) {
00095       // There is a coil bypass, so look at coil stripends.
00096       // Where is the coil center?
00097       double coilpos = 0.5*(ustrip.PartialLength(StripEnd::kNegative) - ustrip.PartialLength(StripEnd::kPositive));
00098       if (x<coilpos) {
00099         distToEnd = (-x) - ( ustrip.GetHalfLength() - ustrip.PartialLength(StripEnd::kNegative));
00100         clippedfrac *= ClipFraction( distToEnd );
00101 
00102       } else {
00103         distToEnd = (x) - ( ustrip.GetHalfLength() - ustrip.PartialLength(StripEnd::kPositive));
00104         clippedfrac *= ClipFraction( distToEnd );
00105       }
00106     }
00107   }
00108 
00109   if (clippedfrac < 0.49) {
00110     MSG("Photon",Msg::kError) << "Error: saw clipfrac less than 0.5: " << clippedfrac<< endl;
00111     clippedfrac = 0.5;
00112   }
00113  
00114   PlexStripEndId end[2];
00115   end[0] = inHit->StripEndId(); end[0].SetEnd(StripEnd::kNegative);
00116   end[1] = inHit->StripEndId(); end[1].SetEnd(StripEnd::kPositive);
00117 
00118   // Energy->MIP
00119   // Apply normalization and birks' constant.
00120   double emips = 
00121     clippedfrac                       // Light that didnt' leak out ends of strip                   
00122     * fOverallLightOutput             // Tuning parameter.
00123     * dE/(1. + fBirksConstant*dEdx)   // Birk's law
00124     / fGeVPerMip;                     // Convert from energy to mips
00125 
00126   
00127   // Useful for deugging:
00128   //MSG("Photon",Msg::kInfo) << "dE: " << dE << "  fGeVPerMip: " << fGeVPerMip << "  Birk suppress: " << 1.0/(1. + fBirksConstant*dEdx) << "  clippfrac: " << clippedfrac
00129   //     << " gives emips = " << emips << endl;
00130 
00131   const Calibrator& cal = Calibrator::Instance();
00132   // MIP->SigMap.  These two values should be the same (in principle)
00133   // but we may as well be pedantic.
00134   double sigmap[2] = {0,0};
00135   double sigcorr[2]= {0,0};
00136   double siglin[2] = {0,0};
00137   double adc[2]    = {0,0};
00138   double meanpe[2] = {0,0};
00139   double pe[2]     = {0,0};
00140   int npe[2]       = {0,0};
00141 
00142   if (ustrip.IsMirrored(StripEnd::kNegative) || ustrip.IsMirrored(StripEnd::kPositive) ) {
00143     // All of the energy goes to one end.
00144   } else {
00145     // The energy is split so that light goes both ways.
00146     emips *= 0.5;
00147   }
00148 
00149   for (int i=0;i<2;i++) {
00150     if ( ustrip.IsMirrored(end[i].GetEnd()) ) {
00151       // This end doesn't read out. Set the output to zero.
00152       sigmap[i] = sigcorr[i] = siglin[i] = pe[i] = 0;
00153       npe[i] = 0;
00154     } else {
00155       // Go through the decalibrtation chain:
00156 
00157       // Mips -> sigmaps.
00158       sigmap[i] =  cal.DecalMIP( emips , end[i]);
00159 
00160       // SigMap->SigCorr
00161       sigcorr[i] = cal.DecalAttenCorrected(sigmap[i],x,end[i]);
00162 
00163       // SigCorr->SigLin (siglin = ADCs for no drift, no linearity)
00164       siglin[i] = cal.DecalStripToStrip(sigcorr[i],end[i]);
00165 
00166       // SigLin->ADCs (Correct for drift. Linearity is applied during PMT simulation.)
00167       adc[i] = cal.DecalDrift(siglin[i],end[i]);
00168 
00169       // SigLin->Pe
00170       meanpe[i] = cal.GetPhotoElectrons(adc[i],end[i]);
00171 
00172       if (isnan(meanpe[i])) {
00173         MSG("Photon",Msg::kError) 
00174           << "Got an expectation value of NaN photons in CalibratedPeComputer.\n"
00175           << "\tHit on stripend " << end[i].AsString() << endl
00176           << "\tPosition: " << x << " m (in strip coords) " << endl
00177           << "\t" << dE/Munits::MeV  << " MeV,  "<< endl
00178           << "\t" << emips      << " MIP,   "<< endl
00179           << "\t" << sigmap[0]  << " / " << sigmap[1]  << " SigMap, "<< endl
00180           << "\t" << sigcorr[0] << " / " << sigcorr[1] << " SigCorr, "<< endl
00181           << "\t" << siglin[0]  << " / " << siglin[1]  << " SigLin, "<< endl
00182           << "\t" << adc[0]  << " / " << adc[1]  << " ADC, "<< endl
00183           << "\t" << meanpe[0]  << " / " << meanpe[1]  << " PE" << endl;
00184         outCount.SetPeNegPos(0,0); // Explicitly do nothing.  
00185         return false;              // and give up on this hit.
00186       }
00187 
00188       // Set to a poisson number about this mean.
00189       pe[i] = fRandom->PoissonD(meanpe[i]);
00190 
00191       // Handle some error conditions:
00192       if (isnan(pe[i])) {
00193         MSG("Photon",Msg::kError)
00194           << "Got a random poisson value of NaN photons in CalibratedPeComputer.\n"
00195           << "\tHit on stripend " << end[i].AsString() << endl
00196           << "\tPosition: " << x << " m (in strip coords) " << endl
00197           << "\t" << dE/Munits::MeV  << " MeV,  "<< endl
00198           << "\t" << emips      << " MIP,   "<< endl
00199           << "\t" << sigmap[0]  << " / " << sigmap[1]  << " SigMap, "<< endl
00200           << "\t" << sigcorr[0] << " / " << sigcorr[1] << " SigCorr, "<< endl
00201           << "\t" << siglin[0]  << " / " << siglin[1]  << " SigLin, "<< endl
00202           << "\t" << adc[0]  << " / " << adc[1]  << " ADC, "<< endl
00203           << "\t" << meanpe[0]  << " / " << meanpe[1]  << " PE" << endl;
00204         outCount.SetPeNegPos(0,0); // Explicitly do nothing.  
00205         return false;              // and give up on this hit.
00206       }
00207         
00208       if (pe[i] < 0) {
00209         MSG("Photon",Msg::kError) 
00210           << "Got -ve PEs on end " << i << endl
00211           << "\tHit on stripend " << end[i].AsString() << endl
00212           << "\tPosition: " << x << " m (in strip coords) " << endl
00213           << "\t" << dE/Munits::MeV  << " MeV,  "<< endl
00214           << "\t" << emips      << " MIP,   "<< endl
00215           << "\t" << sigmap[0]  << " / " << sigmap[1]  << " SigMap, "<< endl
00216           << "\t" << sigcorr[0] << " / " << sigcorr[1] << " SigCorr, "<< endl
00217           << "\t" << siglin[0]  << " / " << siglin[1]  << " SigLin, "<< endl
00218           << "\t" << adc[0]  << " / " << adc[1]  << " ADC, "<< endl
00219           << "\t" << meanpe[0]  << " / " << meanpe[1]  << " PE" << endl;
00220         pe[i] = 0;
00221       }
00222   
00223       if (pe[i] > fMaxPePerHit) {
00224         MSG("Photon",Msg::kError) 
00225           << "Got huge pulse on end " << i << endl
00226           << "\tHit on stripend " << end[i].AsString() << endl
00227           << "\tPosition: " << x << " m (strip coords) " << endl
00228           << "\t" << dE/Munits::MeV  << " MeV,  "<< endl
00229           << "\t" << emips      << " MIP,   "<< endl
00230           << "\t" << sigmap[0]  << " / " << sigmap[1]  << " SigMap, "<< endl
00231           << "\t" << sigcorr[0] << " / " << sigcorr[1] << " SigCorr, "<< endl
00232           << "\t" << siglin[0]  << " / " << siglin[1]  << " SigLin, "<< endl
00233           << "\t" << adc[0]  << " / " << adc[1]  << " ADC, "<< endl
00234           << "\t" << meanpe[0]  << " / " << meanpe[1]  << " PE" << endl
00235           << "\t ... Truncating to " << fMaxPePerHit << " pes." << endl;
00236         pe[i] = fMaxPePerHit;
00237       }
00238       npe[i] = (int) pe[i];
00239     }
00240   }
00241 
00242 
00243   // Set the output.
00244   outCount.SetPeNegPos( npe[0], npe[1] );
00245 
00246   /*
00247    // Consistency check. Disabled in production.
00248   float recomip[2];
00249   for (int i=0;i<2;i++) {
00250     float gain, width;
00251     cal.DecalGainAndWidth(gain, width, end[i]);
00252     float q = (float)(meanpe[i])*gain;
00253     q = cal.GetStripToStripCorrected(q,end[i]);
00254     q = cal.GetAttenCorrected(q,x,end[i]);
00255     q = cal.GetMIP(q,end[i]);
00256     recomip[i] = q;
00257   }
00258   MSG("Photon",Msg::kInfo) 
00259     << "Input mip: " << emips << "  Output: " << recomip[0] << "  " << recomip[1] << endl;
00260   */
00261 
00262   MSG("Photon",Msg::kDebug) << "CalibratedPeComputer: " 
00263                            << dE/Munits::MeV  << " MeV,  "
00264                            << emips << " MIP,   "
00265                            << sigmap[0] + sigmap[1] << " SigMap, "
00266                            << sigcorr[0] + sigcorr[1] << " SigCorr, "
00267                            << siglin[0] + siglin[1] << " SigLin, "
00268                            << meanpe[0] + meanpe[1] << " PE" 
00269                            << endl;
00270   MSG("Photon",Msg::kDebug) << "                      " 
00271                            << outCount.GetPe_Neg() 
00272                            << " (" << pe[0] << ") pe Neg   " 
00273                            << outCount.GetPe_Pos() 
00274                            << " (" << pe[1] << ") pe Pos" 
00275                            << endl;
00276 
00277   return true;
00278 }

void PhotonCalibratedPeComputer::Configure (  const Registry &  r  )  [virtual]
 

Reimplemented from PhotonTransportModule.

Definition at line 31 of file PhotonCalibratedPeComputer.cxx.

References fBirksConstant, fGeVPerMip, fMaxPePerHit, fOverallLightOutput, fScintillatorClipping, Registry::Get(), and MSG.

00032 {
00033   MSG("Photon",Msg::kDebug) << "Configuring CalibratedPeComputer." << std::endl;
00034 
00035   bool ok = true;
00036   ok = ok && r.Get("OverallLightOutput",fOverallLightOutput);
00037   ok = ok && r.Get("BirksConstant",fBirksConstant);
00038   ok = ok && r.Get("GeVPerMip",fGeVPerMip);
00039   ok = ok && r.Get("ScintillatorClipping",fScintillatorClipping);
00040   ok = ok && r.Get("MaxPePerHit",fMaxPePerHit);
00041   if (!ok) MSG("Photon",Msg::kError)
00042              << "Trouble configuring PhotonCalibratedPeComputer." << endl;
00043 }

void PhotonCalibratedPeComputer::Reset (  const VldContext &  cx  )  [virtual]
 

Reimplemented from PhotonTransportModule.

Definition at line 48 of file PhotonCalibratedPeComputer.cxx.

00049 {
00050 }

Member Data Documentation

Double_t PhotonCalibratedPeComputer::fBirksConstant [private]
 

Definition at line 32 of file PhotonCalibratedPeComputer.h.

Referenced by ComputePhotons(), and Configure().

Double_t PhotonCalibratedPeComputer::fGeVPerMip [private]
 

Definition at line 33 of file PhotonCalibratedPeComputer.h.

Referenced by Configure().

Double_t PhotonCalibratedPeComputer::fMaxPePerHit [private]
 

Definition at line 35 of file PhotonCalibratedPeComputer.h.

Referenced by ComputePhotons(), and Configure().

Double_t PhotonCalibratedPeComputer::fOverallLightOutput [private]
 

Definition at line 31 of file PhotonCalibratedPeComputer.h.

Referenced by ComputePhotons(), and Configure().

Int_t PhotonCalibratedPeComputer::fScintillatorClipping [private]
 

Definition at line 34 of file PhotonCalibratedPeComputer.h.

Referenced by Configure().

The documentation for this class was generated from the following files:
Generated on Mon Feb 15 11:10:01 2010 for loon by  doxygen 1.3.9.1