#include <FastHFShowerLibrary.h>

Classes
struct	Hit

Public Member Functions
	FastHFShowerLibrary (edm::ParameterSet const &p)

std::vector< Hit >	getHits (const G4ThreeVector &p, const G4ThreeVector &v, int parCode, double parEnergy, bool &ok, double weight, bool onlyLong=false)

const std::map< CaloHitID, float > &	getHitsMap ()

void const	initHFShowerLibrary (const edm::EventSetup &)

void	initRun ()

bool	isItinFidVolume (G4ThreeVector &)

void	recoHFShowerLibrary (const FSimTrack &myTrack)

	~FastHFShowerLibrary ()

Protected Member Functions
void	extrapolate (int, double)

std::vector< double >	getDDDArray (const std::string &, const DDsvalues_type &, int &)

void	getRecord (int, int)

void	interpolate (int, double)

void	loadEventInfo (TBranch *)

bool	rInside (double r)

void	storePhoton (int j)

Private Attributes
int	anuePDG

int	anumuPDG

int	anutauPDG

bool	applyFidCut

double	backProb

double	cFibre

double	dphi

TBranch *	emBranch

int	emPDG

int	epPDG

int	etaPDG

int	evtPerBin

FastHFFibre *	fibre

int	gammaPDG

int	geantinoPDG

std::vector< double >	gpar

TBranch *	hadBranch

TFile *	hf

std::map< CaloHitID, float >	hitMap

float	libVers

float	listVersion

bool	newForm

int	nMomBin

int	npe

int	nuePDG

HcalNumberingFromDDD *	numberingFromDDD

HcalNumberingScheme *	numberingScheme

int	numuPDG

int	nutauPDG

HFShowerPhotonCollection	pe

HFShowerPhotonCollection *	photo

HFShowerPhotonCollection	photon

int	pi0PDG

std::vector< double >	pmom

double	probMax

double	rMax

double	rMin

int	totEvents

bool	verbose

Detailed Description

Definition at line 38 of file FastHFShowerLibrary.h.

Constructor & Destructor Documentation

FastHFShowerLibrary::FastHFShowerLibrary ( edm::ParameterSet const & p )

Definition at line 28 of file FastHFShowerLibrary.cc.

References applyFidCut, backProb, cFibre, emBranch, event(), evtPerBin, edm::hlt::Exception, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), GeV, hadBranch, hf, i, info(), libVers, listVersion, loadEventInfo(), newForm, nMomBin, pmom, probMax, AlCaHLTBitMon_QueryRunRegistry::string, and totEvents.

                                                                   : fibre(0),
                                                                       hf(0),
                                 emBranch(0),
                                 hadBranch(0),
                                                           numberingScheme(0), 
                                                          numberingFromDDD(0), 
                                       npe(0),
                                                                     photo(0) {
   
   edm::ParameterSet m_HS   = p.getParameter<edm::ParameterSet>("HFShowerLibrary");
   edm::FileInPath fp       = m_HS.getParameter<edm::FileInPath>("FileName");
   std::string pTreeName    = fp.fullPath();
   std::string emName       = m_HS.getParameter<std::string>("TreeEMID");
   std::string hadName      = m_HS.getParameter<std::string>("TreeHadID");
   std::string branchEvInfo = m_HS.getUntrackedParameter<std::string>("BranchEvt","HFShowerLibraryEventInfos_hfshowerlib_HFShowerLibraryEventInfo");
   std::string branchPre    = m_HS.getUntrackedParameter<std::string>("BranchPre","HFShowerPhotons_hfshowerlib_");
   std::string branchPost   = m_HS.getUntrackedParameter<std::string>("BranchPost","_R.obj");
 
   probMax                  = m_HS.getParameter<double>("ProbMax");
   backProb                 = m_HS.getParameter<double>("BackProbability");  
   verbose                  = m_HS.getUntrackedParameter<bool>("Verbosity",false);
   applyFidCut              = m_HS.getParameter<bool>("ApplyFiducialCut");
   cFibre                   = c_light*(m_HS.getParameter<double>("CFibre"));
 
   if (pTreeName.find(".") == 0) pTreeName.erase(0,2);
   const char* nTree = pTreeName.c_str();
   hf                = TFile::Open(nTree);
 
   if (!hf->IsOpen()) { 
     edm::LogError("FastCalorimetry") << "HFShowerLibrary: opening " << nTree 
                          << " failed";
     throw cms::Exception("Unknown", "HFShowerLibrary") 
       << "Opening of " << pTreeName << " fails\n";
   } else {
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: opening " << nTree 
                         << " successfully"; 
   }
 
   newForm = (branchEvInfo == "");
   TTree* event(0);
   if (newForm) event = (TTree *) hf ->Get("HFSimHits");
   else         event = (TTree *) hf ->Get("Events");
   if (event) {
     TBranch *evtInfo(0);
     if (!newForm) {
       std::string info = branchEvInfo + branchPost;
       evtInfo          = event->GetBranch(info.c_str());
     }
     if (evtInfo || newForm) {
       loadEventInfo(evtInfo);
     } else {
       edm::LogError("FastCalorimetry") << "HFShowerLibrary: HFShowerLibrayEventInfo"
                        << " Branch does not exist in Event";
       throw cms::Exception("Unknown", "HFShowerLibrary")
     << "Event information absent\n";
     }
   } else {
     edm::LogError("FastCalorimetry") << "HFShowerLibrary: Events Tree does not "
                          << "exist";
     throw cms::Exception("Unknown", "HFShowerLibrary")
       << "Events tree absent\n";
   }
   
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Library " << libVers 
                       << " ListVersion "    << listVersion 
                       << " Events Total " << totEvents << " and "
                       << evtPerBin << " per bin";
 
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Energies (GeV) with " 
                       << nMomBin    << " bins";
 
   for (int i=0; i<nMomBin; i++) 
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: pmom[" << i << "] = "
                         << pmom[i]/GeV << " GeV";
 
   std::string nameBr = branchPre + emName + branchPost;
   emBranch         = event->GetBranch(nameBr.c_str());
   if (verbose) emBranch->Print();
   nameBr           = branchPre + hadName + branchPost;
   hadBranch        = event->GetBranch(nameBr.c_str());
   if (verbose) hadBranch->Print();
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary:Branch " << emName 
                       << " has " << emBranch->GetEntries() 
                       << " entries and Branch " << hadName 
                       << " has " << hadBranch->GetEntries() 
                       << " entries";
 
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary::No packing information -"
                       << " Assume x, y, z are not in packed form";
 
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Maximum probability cut off " 
                       << probMax << "  Back propagation of light prob. "
                                   << backProb ;
 }

FastHFShowerLibrary::~FastHFShowerLibrary ( )

Definition at line 123 of file FastHFShowerLibrary.cc.

References fibre, hf, numberingFromDDD, numberingScheme, and photo.

                                           {
   if (hf)               hf->Close();
   if (fibre)            delete fibre;
   if (photo)            delete photo;
   if (numberingFromDDD) delete numberingFromDDD;
   if (numberingScheme)  delete numberingScheme; 
 }

Member Function Documentation

void FastHFShowerLibrary::extrapolate	(	int	type,
		double	pin
	)

protected

Definition at line 621 of file FastHFShowerLibrary.cc.

References evtPerBin, getRecord(), j, LogDebug, newForm, nMomBin, npe, pe, photo, photon, pmom, alignCSCRings::r, storePhoton(), totEvents, and w.

Referenced by getHits().

                                                           {
 
   int nrec   = int(pin/pmom[nMomBin-1]);
   double w   = (pin - pmom[nMomBin-1]*nrec)/pmom[nMomBin-1];
   nrec++;
 #ifdef DebugLog
   LogDebug("FastCalorimetry") << "HFShowerLibrary:: Extrapolate for Energy " << pin 
                << " GeV with " << nMomBin << " momentum bins and " 
                << evtPerBin << " entries/bin -- total " << totEvents 
                << " using " << nrec << " records";
 #endif
 
   std::vector<int> irc(nrec);
 
   for (int ir=0; ir<nrec; ir++) {
     double r = G4UniformRand();
     irc[ir] = int(evtPerBin*0.5*r) +(nMomBin-1)*evtPerBin + 1;
     if (irc[ir]<1) {
       edm::LogWarning("FastCalorimetry") << "HFShowerLibrary:: Illegal irc[" << ir 
                   << "] = " << irc[ir] << " now set to 1";
       irc[ir] = 1;
     } else if (irc[ir] > totEvents) {
       edm::LogWarning("FastCalorimetry") << "HFShowerLibrary:: Illegal irc[" << ir 
                   << "] = " << irc[ir] << " now set to "
                   << totEvents;
       irc[ir] = totEvents;
 #ifdef DebugLog
     } else {
       LogDebug("FastCalorimetry") << "HFShowerLibrary::Extrapolation use irc[" 
                       << ir  << "] = " << irc[ir];
 #endif
     }
   }
 
   pe.clear(); 
   npe       = 0;
   int npold = 0;
   for (int ir=0; ir<nrec; ir++) {
     if (irc[ir]>0) {
       getRecord (type, irc[ir]);
       int nPhoton = (newForm) ? photo->size() : photon.size();
       npold      += nPhoton;
       for (int j=0; j<nPhoton; j++) {
     double r = G4UniformRand();
     if (ir != nrec-1 || r < w) {
       storePhoton (j);
     }
       }
 #ifdef DebugLog
       LogDebug("FastCalorimetry") << "HFShowerLibrary: Record [" << ir << "] = " 
                       << irc[ir] << " npold = " << npold;
 #endif
     }
   }
 #ifdef DebugLog
   LogDebug("FastCalorimetry") << "HFShowerLibrary:: uses " << npold << " photons";
 #endif
 
   if (npe > npold || npold == 0)
     edm::LogWarning("FastCalorimetry") << "HFShowerLibrary: Extrapolation Warning == "
                 << nrec << " records " << irc[0] << ", " 
                 << irc[1] << ", ... gives a buffer of " <<npold
                 << " photons and fills " << npe 
                 << " *****";
 #ifdef DebugLog
   else
     LogDebug("FastCalorimetry") << "HFShowerLibrary: Extrapolation == " << nrec
              << " records " << irc[0] << ", " << irc[1] 
              << ", ... gives a buffer of " << npold 
              << " photons and fills " << npe << " PE";
   for (int j=0; j<npe; j++)
     LogDebug("FastCalorimetry") << "Photon " << j << " " << pe[j];
 #endif
 }

std::vector< double > FastHFShowerLibrary::getDDDArray	(	const std::string &	str,
		const DDsvalues_type &	sv,
		int &	nmin
	)

protected

Definition at line 708 of file FastHFShowerLibrary.cc.

References DDfetch(), DDValue::doubles(), edm::hlt::Exception, LogDebug, and relativeConstraints::value.

Referenced by initHFShowerLibrary().

                                      {
 
 #ifdef DebugLog
   LogDebug("FastCalorimetry") << "HFShowerLibrary:getDDDArray called for " << str 
                       << " with nMin " << nmin;
 #endif
 
   DDValue value(str);
   if (DDfetch(&sv,value)) {
 #ifdef DebugLog
     LogDebug("FastCalorimetry") << "HFShowerLibrary:getDDDArray value " << value;
 #endif
 
     const std::vector<double> & fvec = value.doubles();
     int nval = fvec.size();
     if (nmin > 0) {
       if (nval < nmin) {
     edm::LogError("FastCalorimetry") << "HFShowerLibrary : # of " << str 
                   << " bins " << nval << " < " << nmin 
                   << " ==> illegal";
     throw cms::Exception("Unknown", "HFShowerLibrary")
       << "nval < nmin for array " << str << "\n";
       }
     } else {
       if (nval < 2) {
     edm::LogError("FastCalorimetry") << "HFShowerLibrary : # of " << str 
                   << " bins " << nval << " < 2 ==> illegal"
                   << " (nmin=" << nmin << ")";
     throw cms::Exception("Unknown", "HFShowerLibrary")
       << "nval < 2 for array " << str << "\n";
       }
     }
     nmin = nval;
 
     return fvec;
   } else {
     edm::LogError("FastCalorimetry") << "HFShowerLibrary : cannot get array " << str;
 
     throw cms::Exception("Unknown", "HFShowerLibrary") 
       << "cannot get array " << str << "\n";
   }
 }

std::vector< FastHFShowerLibrary::Hit > FastHFShowerLibrary::getHits	(	const G4ThreeVector &	p,
		const G4ThreeVector &	v,
		int	parCode,
		double	parEnergy,
		bool &	ok,
		double	weight,
		bool	onlyLong = `false`
	)

Definition at line 279 of file FastHFShowerLibrary.cc.

References funct::abs(), anuePDG, anumuPDG, anutauPDG, applyFidCut, FastHFFibre::attLength(), backProb, funct::cos(), FastHFShowerLibrary::Hit::depth, HLT_25ns14e33_v1_cff::depth, dphi, emPDG, epPDG, etaPDG, create_public_lumi_plots::exp, extrapolate(), fibre, gammaPDG, geantinoPDG, GeV, gpar, i, interpolate(), LogDebug, nMomBin, npe, nuePDG, numuPDG, nutauPDG, AlCaHLTBitMon_ParallelJobs::p, pe, pi0PDG, pmom, FastHFShowerLibrary::Hit::position, probMax, alignCSCRings::r, diffTwoXMLs::r1, diffTwoXMLs::r2, rInside(), funct::sin(), FastHFShowerLibrary::Hit::time, FastHFFibre::tShift(), detailsBasic3DVector::z, hit::z, and FastHFFibre::zShift().

Referenced by recoHFShowerLibrary().

                                                                            {
 
   std::vector<FastHFShowerLibrary::Hit> hit;
   ok = false;
   if (parCode == pi0PDG || parCode == etaPDG || parCode == nuePDG ||
       parCode == numuPDG || parCode == nutauPDG || parCode == anuePDG ||
       parCode == anumuPDG || parCode == anutauPDG || parCode == geantinoPDG) 
     return hit;
   ok = true;
 
   double tSlice = 0.1*hitPoint.mag()/29.98;
   double pz     = momDir.z(); 
   double zint   = hitPoint.z(); 
 
 // if particle moves from interaction point or "backwards (halo)
   int backward = 0;
   if (pz * zint < 0.) backward = 1;
   
   double sphi   = sin(momDir.phi());
   double cphi   = cos(momDir.phi());
   double ctheta = cos(momDir.theta());
   double stheta = sin(momDir.theta());
 
 #ifdef DebugLog
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: getHits " << parCode
                       << " of energy " << pin/GeV << " GeV"
                       << "  dir.orts " << momDir.x() << ", " <<momDir.y() 
                       << ", " << momDir.z() << "  Pos x,y,z = " 
                       << hitPoint.x() << "," << hitPoint.y() << "," 
                       << hitPoint.z() << "," 
                       << " sphi,cphi,stheta,ctheta  = " << sphi 
                       << ","  << cphi << ", " << stheta << "," << ctheta; 
 #endif    
 
   if (parCode == emPDG || parCode == epPDG || parCode == gammaPDG ) {
     if (pin<pmom[nMomBin-1]) {
       interpolate(0, pin);
     } else {
       extrapolate(0, pin);
     }
   } else {
     if (pin<pmom[nMomBin-1]) {
       interpolate(1, pin);
     } else {
       extrapolate(1, pin);
     }
   }
     
   int nHit = 0;
   FastHFShowerLibrary::Hit oneHit;
   for (int i = 0; i < npe; i++) {
     double zv = std::abs(pe[i].z()); // abs local z  
 
 #ifdef DebugLog
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Hit " << i << " " << pe[i] << " zv " << zv;
 #endif
 
     if (zv <= gpar[1] && pe[i].lambda() > 0 && 
     (pe[i].z() >= 0 || (zv > gpar[0] && (!onlyLong)))) {
       int depth = 1;
       if (onlyLong) {
       } else if (backward == 0) {    // fully valid only for "front" particles
     if (pe[i].z() < 0) depth = 2;// with "front"-simulated shower lib.
       } else {                       // for "backward" particles - almost equal
     double r = G4UniformRand();  // share between L and S fibers
         if (r > 0.5) depth = 2;
       } 
       
 
       // Updated coordinate transformation from local
       //  back to global using two Euler angles: phi and theta
       double pex = pe[i].x();
       double pey = pe[i].y();
 
       double xx = pex*ctheta*cphi - pey*sphi + zv*stheta*cphi; 
       double yy = pex*ctheta*sphi + pey*cphi + zv*stheta*sphi;
       double zz = -pex*stheta + zv*ctheta;
 
       G4ThreeVector pos  = hitPoint + G4ThreeVector(xx,yy,zz);
       zv = std::abs(pos.z()) - gpar[4] - 0.5*gpar[1];
       G4ThreeVector lpos = G4ThreeVector(pos.x(),pos.y(),zv);
 
       zv = fibre->zShift(lpos,depth,0);     // distance to PMT !
 
       double r  = pos.perp();
       double p  = fibre->attLength(pe[i].lambda());
       double fi = pos.phi();
       if (fi < 0) fi += twopi;
       int    isect = int(fi/dphi) + 1;
       isect        = (isect + 1) / 2;
       double dfi   = ((isect*2-1)*dphi - fi);
       if (dfi < 0) dfi = -dfi;
       double dfir  = r * sin(dfi);
 
 #ifdef DebugLog
       edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Position shift " << xx 
                           << ", " << yy << ", "  << zz << ": " << pos 
                           << " R " << r << " Phi " << fi << " Section " 
                           << isect << " R*Dfi " << dfir << " Dist " << zv;
 #endif
 
       zz           = std::abs(pos.z());
       double r1    = G4UniformRand();
       double r2    = G4UniformRand();
       double r3    = -9999.;
       if (backward)     r3    = G4UniformRand();
       if (!applyFidCut) dfir += gpar[5];
 
 #ifdef DebugLog
       edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: rLimits " << rInside(r)
                    << " attenuation " << r1 <<":" << exp(-p*zv) 
                    << " r2 " << r2 << " r3 " << r3 << " rDfi "  
                    << gpar[5] << " zz " 
                    << zz << " zLim " << gpar[4] << ":" 
                    << gpar[4]+gpar[1] << "\n"
                    << "  rInside(r) :" << rInside(r) 
                    << "  r1 <= exp(-p*zv) :" <<  (r1 <= exp(-p*zv))
                    << "  r2 <= probMax :"    <<  (r2 <= probMax*weight)
                    << "  r3 <= backProb :"   <<  (r3 <= backProb) 
                    << "  dfir > gpar[5] :"   <<  (dfir > gpar[5])
                    << "  zz >= gpar[4] :"    <<  (zz >= gpar[4])
                    << "  zz <= gpar[4]+gpar[1] :" 
                    << (zz <= gpar[4]+gpar[1]);   
 #endif
 
       if (rInside(r) && r1 <= exp(-p*zv) && r2 <= probMax*weight && 
       dfir > gpar[5] && zz >= gpar[4] && zz <= gpar[4]+gpar[1] && 
       r3 <= backProb && (depth != 2 || zz >= gpar[4]+gpar[0])) {
 
     oneHit.position = pos;
     oneHit.depth    = depth;
     oneHit.time     = (tSlice+(pe[i].t())+(fibre->tShift(lpos,depth,1)));
     hit.push_back(oneHit);
 #ifdef DebugLog
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Final Hit " << nHit 
                         <<" position " << (hit[nHit].position) 
                         << " Depth " << (hit[nHit].depth) <<" Time " 
                         << tSlice << ":" << pe[i].t() << ":" 
                         << fibre->tShift(lpos,depth,1) << ":" 
                         << (hit[nHit].time);
 #endif
     nHit++;
       }
 #ifdef DebugLog
       else  LogDebug("FastCalorimetry") << "HFShowerLibrary: REJECTED !!!";
 #endif
       if (onlyLong && zz >= gpar[4]+gpar[0] && zz <= gpar[4]+gpar[1]) {
     r1    = G4UniformRand();
     r2    = G4UniformRand();
     if (rInside(r) && r1 <= exp(-p*zv) && r2 <= probMax && dfir > gpar[5]){
       oneHit.position = pos;
       oneHit.depth    = 2;
       oneHit.time     = (tSlice+(pe[i].t())+(fibre->tShift(lpos,2,1)));
       hit.push_back(oneHit);
 #ifdef DebugLog
       edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Final Hit " << nHit 
                    << " position " << (hit[nHit].position) 
                    << " Depth " << (hit[nHit].depth) <<" Time "
                    << tSlice << ":" << pe[i].t() << ":" 
                    << fibre->tShift(lpos,2,1) << ":" 
                    << (hit[nHit].time);
 #endif
       nHit++;
     }
       }
     }
   }
 
 #ifdef DebugLog
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Total Hits " << nHit
                       << " out of " << npe << " PE";
 #endif
 
   if (nHit > npe && !onlyLong)
     edm::LogWarning("FastCalorimetry") << "HFShowerLibrary: Hit buffer " << npe 
                        << " smaller than " << nHit << " Hits";
   return hit;
 
 }

const std::map<CaloHitID,float>& FastHFShowerLibrary::getHitsMap ( )

inline

Definition at line 61 of file FastHFShowerLibrary.h.

References hitMap.

Referenced by CalorimetryManager::HDShowerSimulation(), and CalorimetryManager::reconstruct().

61 { return hitMap; };

FastHFShowerLibrary::hitMap

std::map< CaloHitID, float > hitMap

Definition: FastHFShowerLibrary.h:83

void FastHFShowerLibrary::getRecord	(	int	type,
		int	record
	)

protected

Definition at line 467 of file FastHFShowerLibrary.cc.

References emBranch, hadBranch, j, LogDebug, newForm, photo, photon, and totEvents.

Referenced by extrapolate(), and interpolate().

                                                         {
 
   int nrc     = record-1;
   photon.clear();
   photo->clear();
   if (type > 0) {
     if (newForm) {
       hadBranch->SetAddress(&photo);
       hadBranch->GetEntry(nrc+totEvents);
     } else {
       hadBranch->SetAddress(&photon);
       hadBranch->GetEntry(nrc);
     }
   } else {
     if (newForm) {
       emBranch->SetAddress(&photo);
     } else {
       emBranch->SetAddress(&photon);
     }
     emBranch->GetEntry(nrc);
   }
 
 #ifdef DebugLog
   int nPhoton = (newForm) ? photo->size() : photon.size();
   LogDebug("FastCalorimetry") << "HFShowerLibrary::getRecord: Record " << record
                        << " of type " << type << " with " << nPhoton 
                        << " photons";
   for (int j = 0; j < nPhoton; j++) 
     if (newForm) LogDebug("FastCalorimetry") << "Photon " << j << " " << photo->at(j);
     else         LogDebug("FastCalorimetry") << "Photon " << j << " " << photon[j];
 #endif
 }

void const FastHFShowerLibrary::initHFShowerLibrary ( const edm::EventSetup & iSetup )

Definition at line 131 of file FastHFShowerLibrary.cc.

References DDFilteredView::addFilter(), cFibre, dphi, DDSpecificsFilter::equals, edm::hlt::Exception, fibre, alcazmumu_cfi::filter, DDFilteredView::firstChild(), edm::EventSetup::get(), getDDDArray(), gpar, initRun(), DDFilteredView::mergedSpecifics(), mergeVDriftHistosByStation::name, HLT_25ns14e33_v1_cff::nEta, HLT_25ns14e33_v1_cff::nPhi, numberingFromDDD, numberingScheme, photo, rMax, rMin, DDSpecificsFilter::setCriteria(), AlCaHLTBitMon_QueryRunRegistry::string, and relativeConstraints::value.

Referenced by FamosManager::setupGeometryAndField().

                                                                                {
 
   edm::LogInfo("FastCalorimetry") << "initHFShowerLibrary::initialization"; 
 
   edm::ESTransientHandle<DDCompactView> cpv;
   iSetup.get<IdealGeometryRecord>().get(cpv);
 
   std::string name = "HcalHits";
   std::string attribute = "ReadOutName";
   std::string value     = name;
   DDSpecificsFilter filter;
   DDValue           ddv(attribute,value,0);
   filter.setCriteria(ddv,DDSpecificsFilter::equals);
   DDFilteredView fv(*cpv);
   fv.addFilter(filter);
   bool dodet = fv.firstChild();
   if (dodet) {
     DDsvalues_type sv(fv.mergedSpecifics());
     // Radius (minimum and maximum)
     int nR     = -1;
     std::vector<double> rTable = getDDDArray("rTable",sv,nR);
     rMin = rTable[0];
     rMax = rTable[nR-1];
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: rMIN " << rMin/cm 
                         << " cm and rMax " << rMax/cm;
     // Delta phi
     int nEta   = -1;
     std::vector<double> etaTable = getDDDArray("etaTable",sv,nEta);
     int nPhi   = nEta + nR - 2;
     std::vector<double> phibin   = getDDDArray("phibin",sv,nPhi);
     dphi       = phibin[nEta-1];
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: (Half) Phi Width of wedge " 
                         << dphi/deg;
 
     // Special Geometry parameters
     int ngpar = 7;
     gpar      = getDDDArray("gparHF",sv,ngpar);
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: " << ngpar << " gpar (cm)";
 
     for (int ig=0; ig<ngpar; ig++) 
       edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: gpar[" << ig << "] = "
                           << gpar[ig]/cm << " cm";
   } else {
     edm::LogError("FastCalorimetry") << "HFShowerLibrary: cannot get filtered "
                          << " view for " << attribute << " matching "
                          << name;
     throw cms::Exception("Unknown", "HFShowerLibrary")
       << "cannot match " << attribute << " to " << name <<"\n";
   }
 
   initRun();  
   fibre = new FastHFFibre(name, *cpv, cFibre);
   photo = new HFShowerPhotonCollection;
   numberingFromDDD = new HcalNumberingFromDDD(name, *cpv);
   numberingScheme  = new HcalNumberingScheme();
 }

void FastHFShowerLibrary::initRun ( )

Definition at line 188 of file FastHFShowerLibrary.cc.

References anuePDG, anumuPDG, anutauPDG, emPDG, epPDG, etaPDG, gammaPDG, geantinoPDG, nuePDG, numuPDG, nutauPDG, and pi0PDG.

Referenced by initHFShowerLibrary().

                                   {
 
   geantinoPDG = 0; gammaPDG = 22;
   emPDG   = 11; epPDG    = -11; nuePDG   = 12; anuePDG   = -12;
   numuPDG = 14; anumuPDG = -14; nutauPDG = 16; anutauPDG = -16;
   pi0PDG = 111; etaPDG   = 221;
 
 #ifdef DebugLog
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: Particle codes for e- = " 
                << emPDG << ", e+ = " << epPDG << ", gamma = " 
                << gammaPDG << ", pi0 = " << pi0PDG << ", eta = " 
                << etaPDG << ", geantino = " << geantinoPDG 
                << "\n        nu_e = " << nuePDG << ", nu_mu = " 
                << numuPDG << ", nu_tau = " << nutauPDG 
                << ", anti_nu_e = " << anuePDG << ", anti_nu_mu = " 
                << anumuPDG << ", anti_nu_tau = " << anutauPDG;
 #endif
 }

void FastHFShowerLibrary::interpolate	(	int	type,
		double	pin
	)

protected

Definition at line 529 of file FastHFShowerLibrary.cc.

References evtPerBin, getRecord(), GeV, j, LogDebug, newForm, nMomBin, npe, pe, photo, photon, pmom, alignCSCRings::r, storePhoton(), totEvents, and w.

Referenced by getHits().

                                                           {
 
 #ifdef DebugLog
   LogDebug("FastCalorimetry") << "HFShowerLibrary:: Interpolate for Energy " <<pin/GeV
                       << " GeV with " << nMomBin << " momentum bins and " 
                       << evtPerBin << " entries/bin -- total " << totEvents;
 #endif
 
   int irc[2];
   double w = 0.;
   double r = G4UniformRand();
 
   if (pin<pmom[0]) {
     w = pin/pmom[0];
     irc[1] = int(evtPerBin*r) + 1;
     irc[0] = 0;
   } else {
     for (int j=0; j<nMomBin-1; j++) {
       if (pin >= pmom[j] && pin < pmom[j+1]) {
     w = (pin-pmom[j])/(pmom[j+1]-pmom[j]);
     if (j == nMomBin-2) { 
       irc[1] = int(evtPerBin*0.5*r);
     } else {
       irc[1] = int(evtPerBin*r);
     }
     irc[1] += (j+1)*evtPerBin + 1;
     r = G4UniformRand();
     irc[0] = int(evtPerBin*r) + 1 + j*evtPerBin;
     if (irc[0]<0) {
       edm::LogWarning("FastCalorimetry") << "HFShowerLibrary:: Illegal irc[0] = "
                              << irc[0] << " now set to 0";
       irc[0] = 0;
     } else if (irc[0] > totEvents) {
       edm::LogWarning("FastCalorimetry") << "HFShowerLibrary:: Illegal irc[0] = "
                              << irc[0] << " now set to "<< totEvents;
       irc[0] = totEvents;
     }
       }
     }
   }
   if (irc[1]<1) {
     edm::LogWarning("FastCalorimetry") << "HFShowerLibrary:: Illegal irc[1] = " 
                        << irc[1] << " now set to 1";
 
     irc[1] = 1;
   } else if (irc[1] > totEvents) {
     edm::LogWarning("FastCalorimetry") << "HFShowerLibrary:: Illegal irc[1] = " 
                        << irc[1] << " now set to "<< totEvents;
 
     irc[1] = totEvents;
   }
 
 #ifdef DebugLog
   LogDebug("FastCalorimetry") << "HFShowerLibrary:: Select records " << irc[0] 
                       << " and " << irc[1] << " with weights " << 1-w 
                       << " and " << w;
 #endif
 
   pe.clear(); 
   npe       = 0;
   int npold = 0;
   for (int ir=0; ir < 2; ir++) {
     if (irc[ir]>0) {
       getRecord (type, irc[ir]);
       int nPhoton = (newForm) ? photo->size() : photon.size();
       npold      += nPhoton;
       for (int j=0; j<nPhoton; j++) {
     r = G4UniformRand();
     if ((ir==0 && r > w) || (ir > 0 && r < w)) {
       storePhoton (j);
     }
       }
     }
   }
 
   if ((npe > npold || (npold == 0 && irc[0] > 0)) && !(npe == 0 && npold == 0)) 
     edm::LogWarning("FastCalorimetry") << "HFShowerLibrary: Interpolation Warning =="
                 << " records " << irc[0] << " and " << irc[1]
                 << " gives a buffer of " << npold 
                 << " photons and fills " << npe << " *****";
 #ifdef DebugLog
   else
     LogDebug("FastCalorimetry") << "HFShowerLibrary: Interpolation == records " 
              << irc[0] << " and " << irc[1] << " gives a "
              << "buffer of " << npold << " photons and fills "
              << npe << " PE";
   for (int j=0; j<npe; j++)
     LogDebug("FastCalorimetry") << "Photon " << j << " " << pe[j];
 #endif
 
 }

bool FastHFShowerLibrary::isItinFidVolume ( G4ThreeVector & hitPoint )

Definition at line 262 of file FastHFShowerLibrary.cc.

References applyFidCut, archive::flag, and HFFibreFiducial::PMTNumber().

Referenced by recoHFShowerLibrary().

                                                                   {
   bool flag = true;
   if (applyFidCut) {
     int npmt = HFFibreFiducial:: PMTNumber(hitPoint);
 #ifdef DebugLog
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary::isItinFidVolume:#PMT= " 
                                     << npmt << " for hit point " << hitPoint;
 #endif
     if (npmt <= 0) flag = false;
   }
 #ifdef DebugLog
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary::isItinFidVolume: point " 
                                     << hitPoint << " return flag " << flag;
 #endif
   return flag;
 }

void FastHFShowerLibrary::loadEventInfo ( TBranch * branch )

protected

Definition at line 500 of file FastHFShowerLibrary.cc.

References evtPerBin, GeV, i, libVers, listVersion, nMomBin, pmom, and totEvents.

Referenced by FastHFShowerLibrary().

                                                        {
 
   if (branch) {
     std::vector<HFShowerLibraryEventInfo> eventInfoCollection;
     branch->SetAddress(&eventInfoCollection);
     branch->GetEntry(0);
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary::loadEventInfo loads "
                              << " EventInfo Collection of size "
                              << eventInfoCollection.size() << " records";
     totEvents   = eventInfoCollection[0].totalEvents();
     nMomBin     = eventInfoCollection[0].numberOfBins();
     evtPerBin   = eventInfoCollection[0].eventsPerBin();
     libVers     = eventInfoCollection[0].showerLibraryVersion();
     listVersion = eventInfoCollection[0].physListVersion();
     pmom        = eventInfoCollection[0].energyBins();
   } else {
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary::loadEventInfo loads "
                              << " EventInfo from hardwired numbers";
     nMomBin     = 16;
     evtPerBin   = 5000;
     totEvents   = nMomBin*evtPerBin;
     libVers     = 1.1;
     listVersion = 3.6;
     pmom        = {2,3,5,7,10,15,20,30,50,75,100,150,250,350,500,1000};
   }
   for (int i=0; i<nMomBin; i++) 
     pmom[i] *= GeV;
 }

void FastHFShowerLibrary::recoHFShowerLibrary ( const FSimTrack & myTrack )

Definition at line 207 of file FastHFShowerLibrary.cc.

References HLT_25ns14e33_v1_cff::depth, getHits(), HcalNumberingScheme::getUnitID(), hitMap, i, FSimTrack::id(), isItinFidVolume(), numberingFromDDD, numberingScheme, convertSQLiteXML::ok, FSimTrack::onVFcal(), cond::rpcobgas::time, tmp, CoreSimTrack::type(), HcalNumberingFromDDD::unitID(), RawParticle::vertex(), FSimTrack::vfcalEntrance(), histoStyle::weight, RawParticle::X(), RawParticle::Y(), and RawParticle::Z().

Referenced by CalorimetryManager::HDShowerSimulation(), and CalorimetryManager::reconstruct().

                                                                       {
 
 #ifdef DebugLog
   edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: recoHFShowerLibrary ";
 #endif 
 
   if(!myTrack.onVFcal()) {
 #ifdef DebugLog
     edm::LogInfo("FastCalorimetry") << "HFShowerLibrary: we should not be here ";
 #endif
   }
 
   hitMap.clear();
   double eGen  = 1000.*myTrack.vfcalEntrance().e();                // energy in [MeV]
   double delZv = (myTrack.vfcalEntrance().vertex().Z()>0.0) ? 50.0 : -50.0;
   G4ThreeVector vertex( 10.*myTrack.vfcalEntrance().vertex().X(),
                         10.*myTrack.vfcalEntrance().vertex().Y(),
                         10.*myTrack.vfcalEntrance().vertex().Z()+delZv); // in [mm]
 
   G4ThreeVector direction(myTrack.vfcalEntrance().Vect().X(),
                           myTrack.vfcalEntrance().Vect().Y(),
                           myTrack.vfcalEntrance().Vect().Z());
 
   bool ok;
   double weight = 1.0;                     // rad. damage 
   int parCode   = myTrack.type();
   std::vector<FastHFShowerLibrary::Hit> hits =
               getHits(vertex, direction, parCode, eGen, ok, weight, false);
 
   for (unsigned int i=0; i<hits.size(); ++i) {
     G4ThreeVector pos = hits[i].position;
     int depth         = hits[i].depth;
     double time       = hits[i].time;
 
     if (isItinFidVolume (pos)) {     
       int det = 5;
       int lay = 1;
       uint32_t id = 0;
       HcalNumberingFromDDD::HcalID tmp = numberingFromDDD->unitID(det, pos, depth, lay);
       id = numberingScheme->getUnitID(tmp);
 
       CaloHitID current_id(id,time,myTrack.id());
       std::map<CaloHitID,float>::iterator cellitr;
       cellitr = hitMap.find(current_id);
       if(cellitr==hitMap.end()) {
          hitMap.insert(std::pair<CaloHitID,float>(current_id,1.0));
       } else {
          cellitr->second += 1.0;
       }
     }  // end of isItinFidVolume check 
 
   } // end loop over hits
 
 }

bool FastHFShowerLibrary::rInside ( double r )

protected

Definition at line 461 of file FastHFShowerLibrary.cc.

References rMax, and rMin.

Referenced by getHits().

                                           {
 
   if (r >= rMin && r <= rMax) return true;
   else                        return false;
 }

void FastHFShowerLibrary::storePhoton ( int j )

protected

Definition at line 696 of file FastHFShowerLibrary.cc.

References LogDebug, newForm, npe, pe, photo, and photon.

Referenced by extrapolate(), and interpolate().

                                            {
 
   if (newForm) pe.push_back(photo->at(j));
   else         pe.push_back(photon[j]);
 #ifdef DebugLog
   LogDebug("FastCalorimetry") << "HFShowerLibrary: storePhoton " << j << " npe " 
                       << npe << " " << pe[npe];
 #endif
 
   npe++;
 }