PhotonDefaultModel Class Reference

#include <PhotonDefaultModel.h>

Inheritance diagram for PhotonDefaultModel:
PhotonTransportModule

List of all members.

Public Member Functions

 PhotonDefaultModel ()
virtual void Configure (const Registry &r)
virtual Bool_t ComputePhotons (const DigiScintHit *inHit, Double_t inBluePrescale, Double_t inWlsPrescale, Double_t inGreenPrescale, PhotonCount &outCount)
virtual Bool_t ScintPhotonToFibreHit (const DigiScintHit *hit, const UgliStripHandle &inStrip, DigiPhoton *&outPhoton)
virtual Bool_t FibreHitToGreenPhoton (const UgliStripHandle &inStrip, StripEnd::StripEnd_t inDirection, DigiPhoton *inBluePhoton, DigiPhoton *&outGreenPhoton)
virtual Bool_t GreenPhotonToPe (const UgliStripHandle &inStrip, DigiPhoton *inGreenPhoton, DigiPE *&outPe, StripEnd::StripEnd_t &outEnd)
virtual void MakeNoise (std::vector< DigiPE * > &peList, Double_t t_start, Double_t t_end)
virtual void MakeNoiseOld (std::vector< DigiPE * > &peList, Double_t t_start, Double_t t_end)
virtual void MakeAfterpulses (const std::vector< DigiPE * > &inPeList, std::vector< DigiPE * > &outPeList)
virtual void Print (Option_t *option="") const

Protected Member Functions

virtual void Reset (const VldContext &cx)

Private Member Functions

Float_t GetStripEndRate (PlexStripEndId seid, PlexHandle ph)
 ClassDef (PhotonDefaultModel, 1)

Private Attributes

Double_t fOverallLightOutput
Double_t fBirksConstant
Double_t fScintDecayTime
Double_t fClearFibreIndex
Double_t fFibreRadius
Double_t fFibreFlurorDecayTime
Double_t fOuterCladdingIndex
Double_t fInnerCladdingIndex
Double_t fFibreIndex
Double_t fFibreVelocityFudge
Double_t fFibreCoreRadius
Double_t fReflectorFibreReflec
Double_t fQuantumEfficiency
Double_t fFibreAttenN1
Double_t fFibreAttenN2
Double_t fFibreAttenLength1
Double_t fFibreAttenLength2
Double_t fClearFibreAttenN1
Double_t fClearFibreAttenN2
Double_t fClearFibreAttenLength1
Double_t fClearFibreAttenLength2
Double_t fDarkNoiseRate
Double_t fGreenNoiseRate
Double_t fNoiseWindow
map< UShort_t, Float_t > chanhitsmap
Double_t fPoissonNoiseMean
Double_t fExpTailFraction
Int_t fDetectorNoiseRate
Int_t fUseSimpleNoiseModel

Detailed Description

Definition at line 23 of file PhotonDefaultModel.h.


Constructor & Destructor Documentation

PhotonDefaultModel::PhotonDefaultModel (  ) 

Definition at line 19 of file PhotonDefaultModel.cxx.

References PhotonConfiguration().

00020 {
00021   Configure(PhotonConfiguration()); // Set up default values.
00022 }


Member Function Documentation

PhotonDefaultModel::ClassDef ( PhotonDefaultModel  ,
 
) [private]

Reimplemented from PhotonTransportModule.

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

Reimplemented from PhotonTransportModule.

Definition at line 106 of file PhotonDefaultModel.cxx.

References DigiScintHit::DE(), DigiScintHit::DS(), fBirksConstant, fOverallLightOutput, PhotonTransportModule::fRandom, and PhotonCount::SetBlueTotal().

00111 {
00112   // Number of photons made:
00113   if(!inHit) return false;
00114   double dE = inHit->DE();
00115   double dx = inHit->DS();
00116 
00117   double dEdx = dE/dx;
00118   
00119   //get the number of photons
00120   //from the PDG 2002 Section 27 on Particle Detectors you get
00121   //about 1 photon/100eV of energy deposited.
00122   // .. so, 10photons/keV = 10^7 photons/GeV
00123 
00124   // From rough comparison with REROOT results, the total 
00125   // constant for this algorithm should be ~27500
00126 
00127   // Here's what I can account for:
00128   const double kLightYield = ( 1/100.  // photons/eV
00129                                * 1e9 // GeV/ev
00130                                * inBluePrescale
00131                                * inWlsPrescale
00132                                * inGreenPrescale
00133                                );
00134 
00135   // This fudge factor comes from a comparison of using GMINOS FLS digits
00136   // and this method.  It has no physical justification.
00137   const double kFudgeFactor = 2.0;
00138                         
00139   // Apply normalization and birks' constant.
00140   double meanBlue = ( kLightYield * kFudgeFactor * fOverallLightOutput 
00141                       * dE/(1. + fBirksConstant*dEdx)   );
00142 
00143   // Set to a poisson number about this mean.
00144   int nblue =  fRandom->Poisson(meanBlue);
00145     
00146   //MSG("Photon",Msg::kDebug) << "DefaultComputer: " << dE/Munits::MeV 
00147   //                        << " Mev -> " << meanBlue << " expected blue -> " 
00148   //                        << nblue << " blue photons.  " 
00149   //                        << inBluePrescale << "  " 
00150   //                        << inWlsPrescale << "  "
00151   //                        << inGreenPrescale << "  "
00152   //                        << endl;
00153   
00154 
00155   outCount.SetBlueTotal( nblue );
00156   return true;
00157 }

void PhotonDefaultModel::Configure ( const Registry r  )  [virtual]

Reimplemented from PhotonTransportModule.

Definition at line 27 of file PhotonDefaultModel.cxx.

References fBirksConstant, fClearFibreAttenLength1, fClearFibreAttenLength2, fClearFibreAttenN1, fClearFibreAttenN2, fClearFibreIndex, fDarkNoiseRate, fDetectorNoiseRate, fExpTailFraction, fFibreAttenLength1, fFibreAttenLength2, fFibreAttenN1, fFibreAttenN2, fFibreCoreRadius, fFibreFlurorDecayTime, fFibreIndex, fFibreRadius, fFibreVelocityFudge, fGreenNoiseRate, fInnerCladdingIndex, fNoiseWindow, fOuterCladdingIndex, fOverallLightOutput, fPoissonNoiseMean, fQuantumEfficiency, fReflectorFibreReflec, fScintDecayTime, fUseSimpleNoiseModel, Registry::Get(), Msg::kDebug, MSG, PhotonTransportModule::SetBluePrescaleFactor(), PhotonTransportModule::SetGreenPrescaleFactor(), and PhotonTransportModule::SetWlsPrescaleFactor().

00028 {
00029   // The other classses in the model need configuring:
00030   ScintPhoton::Configure(r);
00031 
00032   MSG("Photon",Msg::kDebug) << "Configuring DefaultModel." << std::endl;
00033 
00034   double tmpd;
00035   int tmpi;
00036   // FIXME!
00037   if (r.Get("OverallLightOutput",tmpd))   fOverallLightOutput  = tmpd;
00038   if (r.Get("BirksConstant",tmpd))        fBirksConstant  = tmpd;
00039   if (r.Get("ScintDecayTime",tmpd))       fScintDecayTime = tmpd;
00040   if (r.Get("ClearFibreIndex",tmpd))      fClearFibreIndex  = tmpd;
00041   if (r.Get("FibreRadius",tmpd))          fFibreRadius  = tmpd;
00042   if (r.Get("FibreFlurorDecayTime",tmpd)) fFibreFlurorDecayTime  = tmpd;
00043 
00044   if (r.Get("OuterCladdingIndex",tmpd))   fOuterCladdingIndex  = tmpd;
00045   if (r.Get("InnerCladdingIndex",tmpd))   fInnerCladdingIndex  = tmpd;
00046   if (r.Get("FibreIndex",tmpd))           fFibreIndex  = tmpd;
00047   if (r.Get("FibreVelocityFudge",tmpd))   fFibreVelocityFudge  = tmpd;
00048   if (r.Get("FibreCoreRadius",tmpd))      fFibreCoreRadius  = tmpd;
00049   if (r.Get("ReflectorFibreReflec",tmpd)) fReflectorFibreReflec  = tmpd;
00050   if (r.Get("QuantumEfficiency",tmpd))    fQuantumEfficiency = tmpd;
00051 
00052   if (r.Get("FibreAttenN1",tmpd))        fFibreAttenN1 = tmpd;
00053   if (r.Get("FibreAttenN2",tmpd))        fFibreAttenN2 = tmpd;
00054   if (r.Get("FibreAttenLength1",tmpd))   fFibreAttenLength1 = tmpd;
00055   if (r.Get("FibreAttenLength2",tmpd))   fFibreAttenLength2 = tmpd;
00056   if (r.Get("ClearFibreAttenN1",tmpd))   fClearFibreAttenN1 = tmpd;
00057   if (r.Get("ClearFibreAttenN2",tmpd))   fClearFibreAttenN2 = tmpd;
00058   if (r.Get("ClearFibreAttenLength1",tmpd)) fClearFibreAttenLength1 = tmpd;
00059   if (r.Get("ClearFibreAttenLength2",tmpd)) fClearFibreAttenLength2 = tmpd ;
00060 
00061   if (r.Get("DarkNoiseRate",tmpd))        fDarkNoiseRate = tmpd; 
00062   if (r.Get("GreenNoiseRate",tmpd))       fGreenNoiseRate = tmpd;
00063   if (r.Get("NoiseWindow",tmpd))          fNoiseWindow = tmpd;
00064   if (r.Get("PoissonNoiseMean",tmpd))     fPoissonNoiseMean = tmpd;
00065   if (r.Get("ExpTailFraction",tmpd))      fExpTailFraction = tmpd;
00066   if (r.Get("FDNoiseRate",tmpi))          fDetectorNoiseRate = tmpi;
00067 
00068   if (r.Get("UseSimpleNoiseModel",tmpi))  fUseSimpleNoiseModel = tmpi;
00069 
00070   // Set the prescale factors.
00071   SetBluePrescaleFactor(1.0); // No prescaling.
00072   SetWlsPrescaleFactor(0.2);  // We always accept a contained ray, no matter what.
00073   SetGreenPrescaleFactor(fQuantumEfficiency); // Phototube photocathode acceptance.
00074 
00075 
00076 }

Bool_t PhotonDefaultModel::FibreHitToGreenPhoton ( const UgliStripHandle inStrip,
StripEnd::StripEnd_t  inDirection,
DigiPhoton inBluePhoton,
DigiPhoton *&  outGreenPhoton 
) [virtual]

Reimplemented from PhotonTransportModule.

Definition at line 245 of file PhotonDefaultModel.cxx.

References fFibreCoreRadius, fFibreFlurorDecayTime, fFibreIndex, fOuterCladdingIndex, PhotonTransportModule::fRandom, UgliStripHandle::IsMirrored(), StripEnd::kNegative, StripEnd::kPositive, n, DigiPhoton::ParentHit(), MuELoss::R, DigiPhoton::T(), and DigiPhoton::X().

00250 {
00251   // Simple model. No wavelength dependencies.
00252   // Return a contained green photon.
00253   // This code currently does not do counting correctly; it just keeps
00254   // simulating until it gets a green photon that is contained in the fibre.
00255   // This implicitly assumes there is no correllation between the blue photon
00256   // and the green photon in any way: one blue photon is as good as another.
00257   // This is both computationally efficient and pretty true.
00258   // The blue photon is assumed to uniformly in the green fibre.
00259 
00260 
00261   outGreenPhoton = NULL;
00262   if(!inBluePhoton) return false;
00263 
00264   // Loop until we have a valid solution.
00265   TVector3 dir;
00266   TVector3 x0;
00267   while(true) {
00268 
00269     // Find the new direction of the green photon.
00270     // Go in the forward or backward direction, depending on the 
00271     // given parameter.
00272     
00273     // Without loss of generality, rotate the system coordinates about the X-axis
00274     // such that the photon start position it on the Z axis between 0 and R.
00275     
00276     // Chose a position.
00277     // FIXME: really, the position should be chosen by 
00278    // computing the depth the blue photon reaches in the fibre.
00279     // But frankly, who cares?
00280 
00281     // Chose radial position as volume-weighted.
00282     double R = fFibreCoreRadius;
00283     double r = sqrt(fRandom->Rndm())*R;
00284     x0.SetXYZ(inBluePhoton->X(),0,r);
00285 
00286     // Choose a direction.
00287     double cosTheta  = 2.*fRandom->Rndm()-1.0;
00288     double sinTheta = sqrt(1.0-cosTheta*cosTheta); //TMath::Sin(TMath::ACos(cosTheta));
00289     double phi = 2.*TMath::Pi()*fRandom->Rndm();
00290     dir.SetXYZ( cosTheta, 
00291                 sinTheta*TMath::Cos(phi),
00292                 sinTheta*TMath::Sin(phi) );
00293     
00294     // Find the position of first reflection:
00295     double xv = r*dir.z();         // xv = x0.y()*dir.y() + x0.z()*dir.z();
00296     double v2 = sinTheta*sinTheta; //  double v2 = dir.y()*dir.y() + dir.z()*dir.z();
00297     if(v2==0) return false;
00298     
00299     double lambda = (-xv +  sqrt(xv*xv + v2*(R*R-r*r)) )/v2;
00300     
00301     TVector3 xhit = x0 + lambda*dir;
00302     
00303     // The surface normal at this place.
00304     TVector3 n(0, -xhit.y()/R, -xhit.z()/R);
00305 
00306     // Find angle of incidence.
00307     double cosAngle =  n.Dot(dir);
00308     double sinAngle = sqrt(1.-cosAngle*cosAngle); //sin(acos(cosAngle));
00309     if(sinAngle > (fOuterCladdingIndex/fFibreIndex)) {
00310       // Solution is good!
00311       break;
00312     } else {
00313       // Return if we want to do counting properly, or if
00314       // we have some sort of correlation between green and blue.
00315       // Otherwise, why waste a perfectly good blue photon? Just
00316       // keep trying till it bloody well DOES go down the fibre.
00317       // return false;
00318     }
00319   }
00320 
00321   // Choose photon direction:
00322   if(inStrip.IsMirrored(StripEnd::kPositive)||inStrip.IsMirrored(StripEnd::kNegative)) {
00323     // If the strip is mirrored, we _don't_ want to make an a priori decision
00324     // since either direction can make a hit on the correct side.
00325 
00326   } else {
00327 
00328     // Let's ensure it's going in the demanded direction:
00329     switch(inDirection) {
00330     case StripEnd::kPositive:
00331       dir.SetX(fabs(dir.X()));
00332       break;
00333       
00334     case StripEnd::kNegative:
00335       dir.SetX(-fabs(dir.X()));
00336       break;
00337       
00338     default:
00339       // Do nothing.
00340       break;    
00341     }
00342   }
00343   
00344   // There's another (unlikely) condition: the photon bounces around exactly sideways
00345   // and so never gets to the end of the tube.
00346   if(dir.X() ==0) return false;
00347 
00348   // Otherwise.. we have a winner.
00349 
00350   // Set up the time correctly now.
00351   Double_t decayTime = fRandom->Exp(fFibreFlurorDecayTime);
00352 
00353   outGreenPhoton = new DigiPhoton( inBluePhoton->ParentHit(),
00354                                    0, // Wavelength is generic.
00355                                    inBluePhoton->T() + decayTime,
00356                                    x0,
00357                                    dir
00358                                    );
00359                                    
00360 
00361   return true;
00362 }

Float_t PhotonDefaultModel::GetStripEndRate ( PlexStripEndId  seid,
PlexHandle  ph 
) [private]

Definition at line 591 of file PhotonDefaultModel.cxx.

References chanhitsmap, RawChannelId::GetChAdd(), RawChannelId::GetCrate(), PlexHandle::GetRawChannelId(), Msg::kDebug, and MSG.

Referenced by MakeNoise().

00593 {
00594   // We could perhaps speed this up by storing the rate by *strip end*
00595   // in the config file, instead of the rate by channel (saves asking
00596   // the plex)
00597   RawChannelId rcid = ph.GetRawChannelId(seid);
00598   UShort_t chadd = rcid.GetChAdd();
00599   Int_t crate = rcid.GetCrate();
00600   MSG("Photon",Msg::kDebug) << "chanhitsmap[" << chadd+5329*crate 
00601                             << "] = " << chanhitsmap[chadd+5329*crate] 
00602                             << endl;
00603   return chanhitsmap[chadd+5329*crate];
00604 }

Bool_t PhotonDefaultModel::GreenPhotonToPe ( const UgliStripHandle inStrip,
DigiPhoton inGreenPhoton,
DigiPE *&  outPe,
StripEnd::StripEnd_t outEnd 
) [virtual]

Reimplemented from PhotonTransportModule.

Definition at line 371 of file PhotonDefaultModel.cxx.

References Munits::c_light, UgliStripHandle::ClearFiber(), fClearFibreAttenLength1, fClearFibreAttenLength2, fClearFibreAttenN1, fClearFibreAttenN2, fClearFibreIndex, PhotonTransportModule::fContext, fFibreAttenLength1, fFibreAttenLength2, fFibreAttenN1, fFibreAttenN2, fFibreIndex, fFibreVelocityFudge, PhotonTransportModule::fRandom, fReflectorFibreReflec, DigiPhoton::GetCosX(), UgliStripHandle::GetHalfLength(), PlexHandle::GetPixelSpotId(), UgliStripHandle::GetSEId(), UgliStripHandle::IsMirrored(), StripEnd::kNegative, StripEnd::kPositive, DigiPhoton::ParentHit(), UgliStripHandle::PartialLength(), PlexStripEndId::SetEnd(), DigiPhoton::T(), UgliStripHandle::WlsBypass(), UgliStripHandle::WlsPigtail(), and DigiPhoton::X().

00375 {
00376   outPe = NULL;
00377   if(!inGreenPhoton) return false;
00378 
00379   // First, figure out which way we're going.
00380   StripEnd::StripEnd_t dir = 
00381     (inGreenPhoton->GetCosX()>0) ? StripEnd::kPositive : StripEnd::kNegative;
00382 
00383   StripEnd::StripEnd_t readout = dir; // Assume we're reading out the end it's going.
00384 
00385   // Is the end we've gotten to mirrored? 
00386   if(inStrip.IsMirrored(dir)) {
00387     // Ok, then we're going to need to flip it.
00388     // Does this photon survive the flip?
00389     if(fRandom->Rndm() > fReflectorFibreReflec) {
00390       return false;
00391     }
00392 
00393     // We survived the reflection.
00394     readout = (dir==StripEnd::kPositive) ? StripEnd::kNegative : StripEnd::kPositive;
00395   }
00396 
00397   // The photon is reflecting it's way down the strip.
00398   // Figure out how far it has to go in each medium.
00399 
00400   double distToCenter = (inGreenPhoton->GetCosX()>0) ? (-inGreenPhoton->X())
00401     : (inGreenPhoton->X());
00402 
00403   Double_t greenDist =
00404     distToCenter // Distance to center of strip.
00405     +inStrip.GetHalfLength()         // Distance to end of strip;
00406     +inStrip.WlsPigtail(readout); // Pigtail on readout end.
00407 
00408   //MSG("Photon",Msg::kVerbose) << "Distance to center of strip: " << distToCenter 
00409   //                          << "  end= " << readout 
00410   //                          << "  halflen=" << inStrip.GetHalfLength()
00411   //                          << "  x="<<inGreenPhoton->X() 
00412   //                          << "  dx=" << inGreenPhoton->GetCosX() << endl;
00413  
00414   // Deal with WLS bypass, if any
00415   double bypassExtra = 0;
00416   if(inStrip.WlsBypass() >0) {    
00417     bypassExtra = ( inStrip.PartialLength(StripEnd::kNegative) + 
00418                     inStrip.PartialLength(StripEnd::kPositive) +
00419                     inStrip.WlsBypass() - 
00420                     inStrip.GetHalfLength()*2.0 );
00421 
00422     // place in the middle of the break;
00423     float x_break =0.5* (inStrip.PartialLength(StripEnd::kNegative) 
00424                          - inStrip.PartialLength(StripEnd::kPositive) );
00425 
00426     if(inGreenPhoton->X() < x_break) {
00427       // If the photon is at -ve end and going +ve, it goes through the bypass
00428       if(dir==StripEnd::kPositive)  greenDist+= bypassExtra;
00429     } else {
00430       // If the photon is at +ve end and going -ve, it goes through the bypass
00431       if(dir==StripEnd::kNegative)  greenDist+= bypassExtra;
00432     }
00433   }
00434 
00435   if(readout!=dir) {
00436     // i.e. we're going towards the mirror end but being read out on the non-mirror end
00437     greenDist+= inStrip.GetHalfLength()*2.0; // Gotta double back...
00438     greenDist+= bypassExtra;                 // ... through the ENTIRE green fibre
00439     greenDist+= inStrip.WlsPigtail(dir)*2.0; // And we have to go through the pigtail, if any
00440                                              // (N.B. This is for caldet, which has pigtails before the reflector)
00441   }
00442 
00443 
00444   // Find the distance down the clear strip.
00445   Double_t clearDist = inStrip.ClearFiber(readout);
00446 
00447   // But... we don't travel in a straight line!
00448   // We bounce around the fibre a lot. This makes the path length longer.
00449   greenDist /= fabs(inGreenPhoton->GetCosX());
00450   clearDist /= fabs(inGreenPhoton->GetCosX());
00451   
00452   //MSG("Photon",Msg::kVerbose) << "After pathlength correction: greenDist = " << greenDist
00453   //                          << "  clearDist = " << clearDist << endl;
00454 
00455   // Compute the probablility we just lose it.
00456   // FIXME: could use full treatment, and at the very least
00457   // shouldn't be hard-coded.
00458   // This treatment comes straight out of GMINOS: init/scintfiber_postffr.F
00459   double norm_g = 1.0/(fFibreAttenN1+fFibreAttenN2);
00460 
00461   // 1 is no attenuation, 0 is complete loss
00462   // So, a small number is a small attenuation prob.
00463   double greenAtten = norm_g * ( fFibreAttenN1 * exp(-greenDist/fFibreAttenLength1) +
00464                                  fFibreAttenN2 * exp(-greenDist/fFibreAttenLength2) );
00465 
00466   //MSG("Photon",Msg::kVerbose) << "Green atten prob = " << greenAtten << endl;
00467   if(fRandom->Rndm() > greenAtten) return false;
00468 
00469   // Clear fibres:
00470   // Again, these numbers are from the GMINOS defaults.
00471   double norm_c =  1.0/(fClearFibreAttenN1+fClearFibreAttenN2);
00472 
00473   // 1 is no attenuation, 0 is complete loss
00474   double clearAtten = norm_c * ( fClearFibreAttenN1 * exp(-clearDist/fClearFibreAttenLength1) +
00475                                  fClearFibreAttenN2 * exp(-clearDist/fClearFibreAttenLength2) );
00476   
00477   //MSG("Photon",Msg::kVerbose) << "Clear atten prob = " << clearAtten << endl;
00478   if(fRandom->Rndm() > clearAtten) return false;
00479 
00480   // ASSUME: no loss at connector interfaces
00481 
00482   // PMT: ASSUME:
00483   // FIXME: Assume 100% quantum efficiency, or efficiency dealt with
00484   // in photon computation.
00485 
00486   // Build the resultant digipe.
00487   // Find the time of arrival.
00488   double time = inGreenPhoton->T()
00489     + greenDist * fFibreIndex * fFibreVelocityFudge / Munits::c_light
00490     + clearDist * fClearFibreIndex * fFibreVelocityFudge / Munits::c_light;
00491 
00492   PlexStripEndId seid =  inStrip.GetSEId();
00493   seid.SetEnd(readout);
00494 
00495   PlexHandle plex(fContext);
00496   PlexPixelSpotId psid = plex.GetPixelSpotId(seid);
00497 
00498   outPe = new DigiPE( time, psid, inGreenPhoton->ParentHit() );
00499   outEnd = readout;
00500 
00501   return true;  
00502 }

void PhotonDefaultModel::MakeAfterpulses ( const std::vector< DigiPE * > &  inPeList,
std::vector< DigiPE * > &  outPeList 
) [virtual]

The "Afterpulser" must impliment this: Create PE for each afterpulse.

Default model: no afterpulsing.

Reimplemented from PhotonTransportModule.

Definition at line 686 of file PhotonDefaultModel.cxx.

00688 {
00695 
00696   return;
00697 }

void PhotonDefaultModel::MakeNoise ( std::vector< DigiPE * > &  peList,
Double_t  t_start,
Double_t  t_end 
) [virtual]

Apply green fibre light to get a new list of PEs.

Reimplemented from PhotonTransportModule.

Definition at line 505 of file PhotonDefaultModel.cxx.

References PhotonTransportModule::fContext, fDetectorNoiseRate, fExpTailFraction, fNoiseWindow, fPoissonNoiseMean, PhotonTransportModule::fRandom, fUseSimpleNoiseModel, VldContext::GetDetector(), GetStripEndRate(), PlexPixelSpotId::IsValid(), DigiSignal::kDarkNoise, Msg::kDebug, Detector::kFar, DigiSignal::kFibreLight, MakeNoiseOld(), MSG, and Munits::ns.

00508 {
00513 
00514 
00515   // So far, this only works for the far detector, so if in ND (or the
00516   // use specifically requested it), use old method:
00517   if ( (fContext.GetDetector() != Detector::kFar) 
00518        || fUseSimpleNoiseModel ) {
00519     MakeNoiseOld(peList, t_start, t_end);
00520     return;
00521   }
00522   
00523   // Work out our window.
00524   double t1 = t_start - fNoiseWindow*0.5;
00525   double t2 = t_end   + fNoiseWindow*0.5;
00526   double dt = t2-t1;
00527 
00528   MSG("Photon",Msg::kDebug) << "Computing noise for window at" << t1 
00529                            << " s for " << (t2-t1)/Munits::ns << " ns."
00530                            << std::endl;
00531 
00532   PlexHandle plex(fContext);
00533   
00534   const std::vector<PlexStripEndId>& stripEnds = plex.GetAllStripEnds();
00535   
00536   // How many fibre hits would we expect?
00537   Double_t fibreHits = dt * fDetectorNoiseRate;
00538   // Now find the number that DID happen...
00539   Int_t nGreenHits = fRandom->Poisson(fibreHits);
00540   MSG("Photon",Msg::kDebug) 
00541     << "Creating " << nGreenHits << " green WLS noise hits.\n";
00542   
00543   // Create the hits.
00544   for(int i=0; i<nGreenHits; i++) {
00545     // Pick a strip end at random.
00546     Int_t whichSeid = fRandom->Integer((int)stripEnds.size());
00547     // Decide whether or not to put this hit on this channel based on
00548     // its normalized noise rate
00549     while ( fRandom->Uniform(0,1) > 
00550             GetStripEndRate(stripEnds[whichSeid] , plex) )
00551       whichSeid = fRandom->Integer((int)stripEnds.size());
00552 
00553     PlexPixelSpotId psid = plex.GetPixelSpotId(stripEnds[whichSeid]);
00554     if(psid.IsValid()) {
00555       Int_t npe = 0;
00556       // Discard the 0 bin
00557       while (npe==0) npe=fRandom->Poisson(fPoissonNoiseMean);
00558       Double_t hitTime = fRandom->Uniform(t1,t2);
00559       for (int j=0; j<npe; ++j) {
00560         DigiPE* pe = new DigiPE( hitTime,
00561                                  psid,
00562                                  DigiSignal::kDarkNoise ); // No Scint Hit associated with the PE.
00563         peList.push_back(pe);
00564       }
00565     }
00566   }
00567   // Do the same for the exponential tail
00568   Int_t nExpHits = fRandom->Poisson(fibreHits*fExpTailFraction);
00569   for (Int_t i=0; i<nExpHits; ++i) {
00570     // Pick a strip end at random.
00571     Int_t whichSeid = fRandom->Integer((int)stripEnds.size());
00572     // Decide whether or not to put this hit on this channel based on
00573     // its normalized noise rate
00574     while ( fRandom->Uniform(0,1) > GetStripEndRate(stripEnds[whichSeid] , plex))
00575       whichSeid = fRandom->Integer((int)stripEnds.size());
00576     
00577     PlexPixelSpotId psid = plex.GetPixelSpotId(stripEnds[whichSeid]);
00578 
00579     Int_t npe = 4+TMath::Nint(fRandom->Exp(0.7));
00580     Double_t hitTime = fRandom->Uniform(t1,t2);
00581     for (int j=0; j<npe; ++j) {
00582       DigiPE* pe = new DigiPE( hitTime,
00583                                psid,
00584                                DigiSignal::kFibreLight ); // No Scint Hit associated with the PE.
00585       peList.push_back(pe);
00586     }
00587   }
00588 }

void PhotonDefaultModel::MakeNoiseOld ( std::vector< DigiPE * > &  peList,
Double_t  t_start,
Double_t  t_end 
) [virtual]

Apply dark noise and green fibre light to get a new list of PEs.

Definition at line 607 of file PhotonDefaultModel.cxx.

References PhotonTransportModule::fContext, fDarkNoiseRate, fGreenNoiseRate, fNoiseWindow, PhotonTransportModule::fRandom, PlexMuxBoxId::GetElecType(), PlexPixelSpotId::IsValid(), DigiSignal::kDarkNoise, Msg::kDebug, DigiSignal::kFibreLight, ElecType::kQIE, ElecType::kVA, MSG, Munits::ns, PlexPixelSpotId::SetPixel(), and PlexPixelSpotId::SetSpot().

Referenced by MakeNoise().

00610 {
00615 
00616   // Work out our window.
00617   double t1 = t_start - fNoiseWindow*0.5;
00618   double t2 = t_end   + fNoiseWindow*0.5;
00619   double dt = t2-t1;
00620 
00621   MSG("Photon",Msg::kDebug) << "Computing noise using simple noise"
00622                             << " model for window at" << t1 
00623                             << " s for " << (t2-t1)/Munits::ns << " ns."
00624                             << std::endl;
00625 
00626   PlexHandle plex(fContext);
00627   // Get a list of PMTs from the Plex.
00628   const std::vector<PlexPixelSpotId>& tubes = plex.GetAllTubes();
00629   
00630   // How many PMT hits would we expect?
00631   Double_t pmtHits = (double)tubes.size() * dt * fDarkNoiseRate;
00632   // Now find the number that DID happen...
00633   Int_t nPmtHits = fRandom->Poisson(pmtHits);
00634   MSG("Photon",Msg::kDebug) << "Creating " << nPmtHits 
00635                             << " dark noise hits.\n";
00636   
00637   // Create the hits.
00638   for(int i=0; i<nPmtHits; i++) {
00639     // Pick a tube at random.
00640     Int_t whichPmt = fRandom->Integer((int)tubes.size());
00641     // Pick a pixelspot.
00642     PlexPixelSpotId psid = tubes[whichPmt];
00643     if(psid.GetElecType()==ElecType::kQIE) {
00644       psid.SetPixel(fRandom->Integer(64));
00645       psid.SetSpot(1);      
00646     }
00647     if(psid.GetElecType()==ElecType::kVA) {
00648       psid.SetPixel(fRandom->Integer(16));
00649       psid.SetSpot(fRandom->Integer(8)+1);      
00650     }
00651     if(psid.IsValid()) {
00652       DigiPE* pe = new DigiPE( fRandom->Uniform(t1,t2),
00653                                psid,
00654                                DigiSignal::kDarkNoise ); // No Scint Hit associated with the PE.
00655       peList.push_back(pe);
00656     }
00657   }
00658 
00659   // Get list of strip ends from the Plex.  
00660   const std::vector<PlexStripEndId>& stripEnds = plex.GetAllStripEnds();
00661 
00662   // How many fibre hits would we expect?
00663   Double_t fibreHits = (double)stripEnds.size() * dt * fGreenNoiseRate;
00664   // Now find the number that DID happen...
00665   Int_t nGreenHits = fRandom->Poisson(fibreHits);
00666   MSG("Photon",Msg::kDebug) << "Creating " << nGreenHits 
00667                             << " green WLS noise hits.\n";
00668   
00669   // Create the hits.
00670   for(int i=0; i<nGreenHits; i++) {
00671     // Pick a strip end at random.
00672     Int_t whichSeid = fRandom->Integer((int)stripEnds.size());
00673     
00674     PlexPixelSpotId psid = plex.GetPixelSpotId(stripEnds[whichSeid]);
00675 
00676     if(psid.IsValid()) {  
00677       DigiPE* pe =  new DigiPE( fRandom->Uniform(t1,t2),
00678                                 psid,
00679                                 DigiSignal::kFibreLight ); // No Scint Hit associated with the PE.
00680       peList.push_back(pe);
00681     }
00682   }
00683 }

void PhotonDefaultModel::Print ( Option_t *  option = ""  )  const [virtual]

Reimplemented from PhotonTransportModule.

Definition at line 702 of file PhotonDefaultModel.cxx.

References fBirksConstant, fOverallLightOutput, fQuantumEfficiency, fReflectorFibreReflec, Msg::kInfo, and MSG.

00703 {
00704   if(option[0]=='c') {
00705     MSG("Photon",Msg::kInfo) << "-<PhotonComputer>     PhotonDefaultModel:-----" << endl;
00706     MSG("Photon",Msg::kInfo) << "  OverallLightOutput:   " << fOverallLightOutput << endl;
00707     MSG("Photon",Msg::kInfo) << "  Birk's Constant:      " << fBirksConstant << endl;   
00708   }
00709   if(option[0]=='b') {
00710     MSG("Photon",Msg::kInfo) << "-<BluePhotonTracker>  PhotonDefaultModel:-----" << endl;
00711     //FIXME
00712   }
00713   if(option[0]=='w') {
00714     MSG("Photon",Msg::kInfo) << "-<WlsFibreModel>      PhotonDefaultModel:-----" << endl;  
00715     // FIXME
00716   }
00717   if(option[0]=='g') {
00718     MSG("Photon",Msg::kInfo) << "-<GreenPhotonTracker> PhotonDefaultModel:-----" << endl;  
00719     MSG("Photon",Msg::kInfo) << "  ReflectorFibreReflec: " << fReflectorFibreReflec << endl; 
00720     MSG("Photon",Msg::kInfo) << "  Quantum Efficiency:   " << fQuantumEfficiency    << endl; 
00721     
00722   }
00723   if(option[0]=='n') {
00724     MSG("Photon",Msg::kInfo) << "-<NoiseMaker>          PhotonDefaultModel:-----" << endl;  
00725     //FIXME
00726   }
00727 }

void PhotonDefaultModel::Reset ( const VldContext cx  )  [protected, virtual]

Reimplemented from PhotonTransportModule.

Definition at line 79 of file PhotonDefaultModel.cxx.

References chanhitsmap, PhotonTransportModule::fContext, fUseSimpleNoiseModel, ChannelNoiseRates::GetChAdd(), ChannelNoiseRates::GetCrate(), VldContext::GetDetector(), ChannelNoiseRates::GetNormRate(), DbiResultPtr< T >::GetNumRows(), DbiResultPtr< T >::GetRow(), Detector::kFar, Msg::kFatal, and MSG.

00080 {
00081   // Ensure tables are up-to-date.
00082   
00083   // VldContext should be set by now, so can set up noise rates
00084 
00085   // Set up the channel noise rates if in far detector
00086   if ( (fContext.GetDetector() == Detector::kFar) 
00087        && fUseSimpleNoiseModel==0 ) {
00088     DbiResultPtr<ChannelNoiseRates> resPtr(fContext);
00089     if (!resPtr.GetNumRows())
00090       MSG("Photon", Msg::kFatal) << "Got no DB rows with noise rates!" 
00091                                  << endl;
00092     
00093     for (UInt_t irow=0; irow < resPtr.GetNumRows(); ++irow) {
00094       const ChannelNoiseRates* cnr = resPtr.GetRow(irow);
00095       chanhitsmap[cnr->GetChAdd()+5329*cnr->GetCrate()]=cnr->GetNormRate();
00096     }
00097   }
00098 }

Bool_t PhotonDefaultModel::ScintPhotonToFibreHit ( const DigiScintHit hit,
const UgliStripHandle inStrip,
DigiPhoton *&  outPhoton 
) [virtual]

Reimplemented from PhotonTransportModule.

Definition at line 170 of file PhotonDefaultModel.cxx.

References PhotonTransportModule::fRandom, fScintDecayTime, ScintPhoton::GeomError(), UgliStripHandle::GetHalfLength(), ScintPhoton::InFibre(), StripEnd::kNegative, StripEnd::kPositive, UgliStripHandle::PartialLength(), ScintPhoton::Propagate(), ScintPhoton::Reset(), DigiScintHit::T1(), DigiScintHit::T2(), UgliStripHandle::WlsBypass(), DigiPhoton::X(), DigiScintHit::X1(), DigiScintHit::X2(), DigiScintHit::Y1(), DigiScintHit::Y2(), DigiScintHit::Z1(), and DigiScintHit::Z2().

00173 {
00174   outPhoton = NULL;
00175   // This default method uses the ScintPhoton class to do the dirty work.
00176 
00177   
00178   if(!hit) return false;
00179 
00180   // Make the photon.
00181   ScintPhoton*  photon = new ScintPhoton();  
00182 
00183   // Find the location/time for the blue photon.
00184 
00185   Double_t delta = fRandom->Rndm();
00186   TVector3 x1(hit->X1(),hit->Y1(),hit->Z1());
00187   TVector3 x2(hit->X2(),hit->Y2(),hit->Z2());
00188 
00189   TVector3 p = x1 + delta*(x2-x1);
00190   Double_t t = hit->T1() + delta * (hit->T2() - hit->T1())
00191     + fRandom->Exp(fScintDecayTime);
00192   
00193 
00194   // The point of this loop is to ensure that we actually manage to 
00195   // simulate a photon successfully.. it doesn't just disappear, but
00196   // is successfully tracked to SOMEWHERE.
00197   bool ok;
00198   do {
00199     TVector3 dir(0,0,0); // let the scint photon do it.
00200     photon->Reset(fRandom,
00201                   hit,
00202                   inStrip,
00203                   p, t,
00204                   dir);
00205     
00206     photon->Propagate();
00207     ok = ! photon->GeomError(); // Geometry errors are bad.
00208   } while(!ok);
00209 
00210   // Successful photon.. may or may not have hit green.
00211   outPhoton = (DigiPhoton*) photon;
00212 
00213   // Simple check for intra-strip light leakage. Not completely correct; this SHOULD
00214   // be part of the photon tracking code.. but simple and straightforward.
00215   if(inStrip.WlsBypass() > 0.) {
00216     float x = outPhoton->X();
00217     // Deal with the case where the hit was near the edge of a two-part strip.
00218     // See if the hit is contained or not.
00219     if( ( (x + inStrip.GetHalfLength()) > inStrip.PartialLength(StripEnd::kNegative) ) 
00220         && ( (inStrip.GetHalfLength() - x) > inStrip.PartialLength(StripEnd::kPositive) )
00221         ) {
00222 
00223       // The photon got leaked out into the coil hole bypass
00224       return false;
00225     }
00226   }
00227 
00228   return (photon->InFibre()); 
00229 }


Member Data Documentation

map<UShort_t, Float_t> PhotonDefaultModel::chanhitsmap [private]

Definition at line 97 of file PhotonDefaultModel.h.

Referenced by GetStripEndRate(), and Reset().

Definition at line 68 of file PhotonDefaultModel.h.

Referenced by ComputePhotons(), Configure(), and Print().

Definition at line 88 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 89 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 86 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 87 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 70 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 91 of file PhotonDefaultModel.h.

Referenced by Configure(), and MakeNoiseOld().

Definition at line 101 of file PhotonDefaultModel.h.

Referenced by Configure(), and MakeNoise().

Definition at line 100 of file PhotonDefaultModel.h.

Referenced by Configure(), and MakeNoise().

Definition at line 84 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 85 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 82 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 83 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 78 of file PhotonDefaultModel.h.

Referenced by Configure(), and FibreHitToGreenPhoton().

Definition at line 72 of file PhotonDefaultModel.h.

Referenced by Configure(), and FibreHitToGreenPhoton().

Double_t PhotonDefaultModel::fFibreIndex [private]

Definition at line 76 of file PhotonDefaultModel.h.

Referenced by Configure(), FibreHitToGreenPhoton(), and GreenPhotonToPe().

Definition at line 71 of file PhotonDefaultModel.h.

Referenced by Configure().

Definition at line 77 of file PhotonDefaultModel.h.

Referenced by Configure(), and GreenPhotonToPe().

Definition at line 92 of file PhotonDefaultModel.h.

Referenced by Configure(), and MakeNoiseOld().

Definition at line 75 of file PhotonDefaultModel.h.

Referenced by Configure().

Definition at line 93 of file PhotonDefaultModel.h.

Referenced by Configure(), MakeNoise(), and MakeNoiseOld().

Definition at line 74 of file PhotonDefaultModel.h.

Referenced by Configure(), and FibreHitToGreenPhoton().

Definition at line 67 of file PhotonDefaultModel.h.

Referenced by ComputePhotons(), Configure(), and Print().

Definition at line 99 of file PhotonDefaultModel.h.

Referenced by Configure(), and MakeNoise().

Definition at line 80 of file PhotonDefaultModel.h.

Referenced by Configure(), and Print().

Definition at line 79 of file PhotonDefaultModel.h.

Referenced by Configure(), GreenPhotonToPe(), and Print().

Definition at line 69 of file PhotonDefaultModel.h.

Referenced by Configure(), and ScintPhotonToFibreHit().

Definition at line 102 of file PhotonDefaultModel.h.

Referenced by Configure(), MakeNoise(), and Reset().


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

Generated on 22 Nov 2017 for loon by  doxygen 1.6.1