GflashProtonShowerProfile Class Reference

#include <GflashProtonShowerProfile.h>

Inheritance diagram for GflashProtonShowerProfile:

Public Member Functions
	GflashProtonShowerProfile (edm::ParameterSet parSet)
void	loadParameters ()
	~GflashProtonShowerProfile ()
Public Member Functions inherited from GflashHadronShowerProfile
std::vector< GflashHit > &	getGflashHitList ()
GflashShowino *	getGflashShowino ()
	GflashHadronShowerProfile (edm::ParameterSet parSet)
void	hadronicParameterization ()
void	initialize (int showerType, double energy, double globalTime, double charge, Gflash3Vector &position, Gflash3Vector &momentum)
virtual	~GflashHadronShowerProfile ()

Additional Inherited Members
Protected Member Functions inherited from GflashHadronShowerProfile
double	depthScale (double ssp, double ssp0, double length)
void	doCholeskyReduction (double **cc, double **vv, const int ndim)
double	fLnE1 (double einc, const double *par)
double	fTanh (double einc, const double *par)
double	gammaProfile (double alpha, double beta, double depth, double lengthUnit)
void	getFluctuationVector (double *lowTriangle, double *correlationVector)
int	getNumberOfSpots (Gflash::CalorimeterNumber kCalor)
double	hoProfile (double pathLength, double refDepth)
Gflash3Vector	locateHitPosition (GflashTrajectoryPoint &point, double lateralArm)
double	longitudinalProfile ()
double	medianLateralArm (double depth, Gflash::CalorimeterNumber kCalor)
void	setEnergyScale (double einc, Gflash3Vector ssp)
double	twoGammaProfile (double *par, double depth, Gflash::CalorimeterNumber kIndex)
Protected Attributes inherited from GflashHadronShowerProfile
double	averageSpotEnergy [Gflash::kNumberCalorimeter]
double	energyScale [Gflash::kNumberCalorimeter]
double	lateralPar [Gflash::kNumberCalorimeter][Gflash::Nrpar]
double	longEcal [Gflash::NPar]
double	longHcal [Gflash::NPar]
double	theBField
bool	theGflashHcalOuter
std::vector< GflashHit >	theGflashHitList
GflashHistogram *	theHisto
edm::ParameterSet	theParSet
GflashShowino *	theShowino

Detailed Description

Definition at line 6 of file GflashProtonShowerProfile.h.

Constructor & Destructor Documentation

GflashProtonShowerProfile::GflashProtonShowerProfile ( edm::ParameterSet  parSet )

inline

Definition at line 12 of file GflashProtonShowerProfile.h.

                                                     : 
    GflashHadronShowerProfile (parSet) {}; 

GflashProtonShowerProfile::~GflashProtonShowerProfile ( )

inline

Definition at line 14 of file GflashProtonShowerProfile.h.

{};

Member Function Documentation

void GflashProtonShowerProfile::loadParameters ( )

virtual

Reimplemented from GflashHadronShowerProfile.

Definition at line 4 of file GflashProtonShowerProfile.cc.

References GflashHadronShowerProfile::depthScale(), GflashHadronShowerProfile::energyScale, create_public_lumi_plots::exp, GflashHadronShowerProfile::fTanh(), GflashShowino::getEnergy(), GflashHadronShowerProfile::getFluctuationVector(), GflashShowino::getPositionAtShower(), GflashShowino::getShowerType(), i, j, Gflash::kENCA, Gflash::kESPM, Gflash::kHB, Gflash::kHE, GflashHadronShowerProfile::lateralPar, Gflash::LengthCrystalEB, Gflash::LengthCrystalEE, create_public_lumi_plots::log, GflashHadronShowerProfile::longEcal, GflashHadronShowerProfile::longHcal, max(), Gflash::NPar, Gflash::Nrpar, position, Gflash::pro_correl_hadem, Gflash::pro_emscale, Gflash::pro_hadscale, Gflash::pro_par, Gflash::pro_rho, Gflash::pro_rpar, diffTwoXMLs::r1, diffTwoXMLs::r2, Gflash::RFrontCrystalEB, Gflash::Rmin, mathSSE::sqrt(), GflashHadronShowerProfile::theShowino, Gflash::ZFrontCrystalEE, and Gflash::Zmin.

{
  double einc = theShowino->getEnergy();
  Gflash3Vector position = theShowino->getPositionAtShower();
  int showerType = theShowino->getShowerType();

  // energy scale
  double energyMeanHcal = 0.0;
  double energySigmaHcal = 0.0;
  
  double r1 = 0.0;
  double r2 = 0.0;

  if(showerType == 0 || showerType == 1 || showerType == 4 || showerType == 5) {

    //@@@ energy dependent energyRho based on tuning with testbeam data
    double energyRho =  fTanh(einc,Gflash::pro_correl_hadem); 

    if(showerType == 0 || showerType == 1) {
      do {
    r1 = CLHEP::RandGaussQ::shoot();

    energyScale[Gflash::kESPM] = einc*(fTanh(einc,Gflash::pro_emscale[0]) + fTanh(einc,Gflash::pro_emscale[1])*r1);

    //LogNormal mean and sigma of Hcal energy
    energyMeanHcal  = (fTanh(einc,Gflash::pro_hadscale[0]) +
               fTanh(einc,Gflash::pro_hadscale[1])*depthScale(position.getRho(),Gflash::RFrontCrystalEB,Gflash::LengthCrystalEB));
    energySigmaHcal = (fTanh(einc,Gflash::pro_hadscale[2]) +  
               fTanh(einc,Gflash::pro_hadscale[3])*depthScale(position.getRho(),Gflash::RFrontCrystalEB,Gflash::LengthCrystalEB));
         
    r2 = CLHEP::RandGaussQ::shoot();
    energyScale[Gflash::kHB] = exp(energyMeanHcal+energySigmaHcal*(energyRho*r1 + sqrt(1.0- energyRho*energyRho)*r2 ));
      } while ( energyScale[Gflash::kESPM] < 0 || energyScale[Gflash::kHB] > einc*1.5 ); 
    }
    else {
      do {
    r1 = CLHEP::RandGaussQ::shoot();
    energyScale[Gflash::kENCA] = einc*(fTanh(einc,Gflash::pro_emscale[0]) + fTanh(einc,Gflash::pro_emscale[1])*r1);
    
    //@@@extend depthScale for HE
    energyMeanHcal  = (fTanh(einc,Gflash::pro_hadscale[0]) + 
               fTanh(einc,Gflash::pro_hadscale[1])*depthScale(std::fabs(position.getZ()),Gflash::ZFrontCrystalEE,Gflash::LengthCrystalEE));
    energySigmaHcal = (fTanh(einc,Gflash::pro_hadscale[2]) +
               fTanh(einc,Gflash::pro_hadscale[3])*depthScale(std::fabs(position.getZ()),Gflash::ZFrontCrystalEE,Gflash::LengthCrystalEE));
    r2 = CLHEP::RandGaussQ::shoot();
    energyScale[Gflash::kHE] = exp(energyMeanHcal+energySigmaHcal*(energyRho*r1 + sqrt(1.0- energyRho*energyRho)*r2 ));
      } while ( energyScale[Gflash::kENCA] < 0 || energyScale[Gflash::kHE] > einc*1.5 ); 
    }
  }
  else if(showerType == 2 || showerType == 3 || showerType == 6 || showerType == 7) { 
    //Hcal response for mip-like pions (mip)
    
    energyMeanHcal  = fTanh(einc,Gflash::pro_hadscale[4]);
    energySigmaHcal = fTanh(einc,Gflash::pro_hadscale[5]);

    double gap_corr = einc*fTanh(einc,Gflash::pro_hadscale[6]);

    if(showerType == 2 || showerType == 3) {
      energyScale[Gflash::kESPM] = 0.0;

      do {
    r1 = CLHEP::RandGaussQ::shoot();
    energyScale[Gflash::kHB] = exp(energyMeanHcal+energySigmaHcal*r1);
      } while ( energyScale[Gflash::kHB] > einc*1.5 );

      if(showerType == 2) {
        energyScale[Gflash::kHE] = std::max(0.0,energyScale[Gflash::kHB]
                                 - gap_corr*depthScale(position.getRho(),Gflash::Rmin[Gflash::kHB],28.));
      }
    }
    else if(showerType == 6 || showerType == 7 ) {
      energyScale[Gflash::kENCA] = 0.0;

      do {
    r1 = CLHEP::RandGaussQ::shoot();
        energyMeanHcal += std::log(1.0-fTanh(einc,Gflash::pro_hadscale[7]));
        energyScale[Gflash::kHE] = exp(energyMeanHcal+energySigmaHcal*r1);
      } while ( energyScale[Gflash::kHE] > einc*1.5 );

      if(showerType == 6) {
        energyScale[Gflash::kHE] = std::max(0.0,energyScale[Gflash::kHE]
                                 - gap_corr*depthScale(std::fabs(position.getZ()),Gflash::Zmin[Gflash::kHE],66.));
      }
    }
  }

  // parameters for the longitudinal profiles
  //@@@check longitudinal profiles of endcaps for possible variations

  double *rhoHcal = new double [2*Gflash::NPar];
  double *correlationVectorHcal = new double [Gflash::NPar*(Gflash::NPar+1)/2];

  //@@@until we have a separate parameterization for Endcap 

  if(showerType>3) {
    showerType -= 4;
  }
  //no separate parameterization before crystal
  if(showerType==0) showerType = 1; 

  //Hcal parameters are always needed regardless of showerType
  for(int i = 0 ; i < 2*Gflash::NPar ; i++ ) {
    rhoHcal[i] = fTanh(einc,Gflash::pro_rho[i + showerType*2*Gflash::NPar]);
  }

  getFluctuationVector(rhoHcal,correlationVectorHcal);

  double normalZ[Gflash::NPar];
  for (int i = 0; i < Gflash::NPar ; i++) normalZ[i] = CLHEP::RandGaussQ::shoot();
  
  for(int i = 0 ; i < Gflash::NPar ; i++) {
    double correlationSum = 0.0;

    for(int j = 0 ; j < i+1 ; j++) {
      correlationSum += correlationVectorHcal[i*(i+1)/2+j]*normalZ[j];
    }
    longHcal[i] = fTanh(einc,Gflash::pro_par[i+showerType*Gflash::NPar]) +
          fTanh(einc,Gflash::pro_par[i+(4+showerType)*Gflash::NPar])*correlationSum;
  }
  delete [] rhoHcal;
  delete [] correlationVectorHcal;

  // lateral parameters for Hcal

  for (int i = 0 ; i < Gflash::Nrpar ; i++) {
    lateralPar[Gflash::kHB][i] = fTanh(einc,Gflash::pro_rpar[i+showerType*Gflash::Nrpar]);
    lateralPar[Gflash::kHE][i] = lateralPar[Gflash::kHB][i];
  }

  //Ecal parameters are needed if and only if the shower starts inside the crystal

  if(showerType == 1) {

    double *rhoEcal = new double [2*Gflash::NPar];
    double *correlationVectorEcal = new double [Gflash::NPar*(Gflash::NPar+1)/2];
    for(int i = 0 ; i < 2*Gflash::NPar ; i++ ) rhoEcal[i] = fTanh(einc,Gflash::pro_rho[i]);

    getFluctuationVector(rhoEcal,correlationVectorEcal);

    for(int i = 0 ; i < Gflash::NPar ; i++) normalZ[i] = CLHEP::RandGaussQ::shoot();
    for(int i = 0 ; i < Gflash::NPar ; i++) {
      double correlationSum = 0.0;

      for(int j = 0 ; j < i+1 ; j++) {
    correlationSum += correlationVectorEcal[i*(i+1)/2+j]*normalZ[j];
      }
      longEcal[i] = fTanh(einc,Gflash::pro_par[i]) +
                0.5*fTanh(einc,Gflash::pro_par[i+4*Gflash::NPar])*correlationSum;
    }

    delete [] rhoEcal;
    delete [] correlationVectorEcal;

    // lateral parameters for Ecal

    for (int i = 0 ; i < Gflash::Nrpar ; i++) {
      lateralPar[Gflash::kESPM][i] = fTanh(einc,Gflash::pro_rpar[i]);
      lateralPar[Gflash::kENCA][i] = lateralPar[Gflash::kESPM][i];
    }
  }

}