00001 #include "FastShowerUtils/Normalisers/HadEarlyEcalCoreFracBar.h"
00002
00003 #include "FastShowerUtils/ParticleParameters.h"
00004
00005 #include "CLHEP/Random/RandFlat.h"
00006
00007 #include <algorithm>
00008 #include <cmath>
00009
00010 namespace FastShower{
00011
00012 double HadEarlyEcalCoreFracBar::value(const ParticleParameters& pp) const{
00013
00014 double energy = pp.energy();
00015
00016 double p1 = -0.0537748*std::exp(-0.196332*energy) + 1.34224;
00017 double p2 = 0.100296*std::exp(-0.108459*energy) + 0.183062;
00018 double p3 = 3.78374*std::exp(-0.04592*energy) - 4.39248;
00019 double p4 = -0.9031114*std::exp(-0.00512578*energy) + 1.99314;
00020 double p5 = 2.7764447 + 0.0041359*energy;
00021
00022 double result;
00023 double rndm = RandFlat::shoot();
00024 double yref = std::exp(p4*std::log(p5/(p3+std::exp((p1-0.)/p2))));
00025 if ((rndm>=yref) && (rndm<=1.0)) {
00026 result = p1-p2*std::log(p5*std::exp(-std::log(rndm)/p4)-p3);
00027 }else if ((rndm>0.) && (rndm<yref)) {
00028 result = RandFlat::shoot()*0.02;
00029 }else{
00030 result = RandFlat::shoot();
00031 }
00032 result = max(result, 0.0);
00033 result = min(result, 1.0);
00034
00035 return result;
00036 }
00037
00038 IFnOfParticleParameters* HadEarlyEcalCoreFracBar::clone() const {
00039 return new HadEarlyEcalCoreFracBar(*this);
00040 }
00041 }
00042
00043
00044