#include <MSLayersKeeperX0Averaged.h>

Inheritance diagram for MSLayersKeeperX0Averaged:

Public Member Functions
virtual void	init (const edm::EventSetup &iSetup)

virtual MSLayer	layer (const DetLayer *layer) const

virtual const MSLayersAtAngle &	layers (float cotTheta) const

	MSLayersKeeperX0Averaged ()

virtual	~MSLayersKeeperX0Averaged ()

Public Member Functions inherited from MSLayersKeeper
virtual	~MSLayersKeeper ()

Private Attributes
bool	isInitialised

MSLayersAtAngle	theLayersData

Additional Inherited Members
Protected Types inherited from MSLayersKeeper
typedef MSLayer::DataX0	DataX0

Static Protected Member Functions inherited from MSLayersKeeper
static const DataX0 &	getDataX0 (const MSLayer &l)

static void	setDataX0 (MSLayer &l, const DataX0 &x0Data)

Detailed Description

Definition at line 6 of file MSLayersKeeperX0Averaged.h.

Constructor & Destructor Documentation

MSLayersKeeperX0Averaged::MSLayersKeeperX0Averaged ( )

inline

Definition at line 8 of file MSLayersKeeperX0Averaged.h.

8 : isInitialised(false) { }

MSLayersKeeperX0Averaged::isInitialised

bool isInitialised

Definition: MSLayersKeeperX0Averaged.h:17

virtual MSLayersKeeperX0Averaged::~MSLayersKeeperX0Averaged ( )

inlinevirtual

Definition at line 9 of file MSLayersKeeperX0Averaged.h.

9 { }

Member Function Documentation

void MSLayersKeeperX0Averaged::init ( const edm::EventSetup & iSetup )

virtual

Reimplemented from MSLayersKeeper.

Definition at line 7 of file MSLayersKeeperX0Averaged.cc.

References GeomDetEnumerators::barrel, MultipleScatteringGeometry::detLayers(), GeomDetEnumerators::endcap, eta(), spr::find(), MSLayersAtAngle::findLayer(), relativeConstraints::geom, MSLayersKeeper::getDataX0(), i, MSLayersKeeperX0AtEta::init(), isInitialised, MSLayersKeeperX0AtEta::layers(), LayerTriplets::layers(), pileupCalc::nbins, MultipleScatteringGeometry::otherLayers(), MSLayersKeeper::setDataX0(), funct::sin(), slope, mathSSE::sqrt(), MSLayer::sumX0D(), MSLayersKeeperX0AtEta::theHalfNBins, theLayersData, MSLayersKeeperX0AtEta::theLayersData, and MSLayer::DataX0::x0.

 {
   if (isInitialised) return;
   isInitialised = true;
 //  cout << "HERE INITIALISATION! MSLayersKeeperX0Averaged"<<endl;
   MSLayersKeeperX0AtEta layersX0Eta;
   layersX0Eta.init(iSetup);
   MultipleScatteringGeometry geom(iSetup);
   vector<MSLayer> allLayers = geom.detLayers(iSetup);
   vector<MSLayer>::iterator it;
   for (int i=-1;i<=1;i++) {
     float eta = i*(-1.8);
     vector<MSLayer> tmpLayers = geom.otherLayers(eta,iSetup);
     vector<MSLayer>::const_iterator ic;
     for (ic = tmpLayers.begin(); ic != tmpLayers.end(); ic++) {
       it = find(allLayers.begin(), allLayers.end(), *ic);  
       if (it == allLayers.end()) allLayers.push_back(*ic);
     }
   }
 
   
 
 for (it = allLayers.begin(); it != allLayers.end(); it++) {
     float cotTheta = (it->face()==GeomDetEnumerators::barrel) ?
         it->range().mean()/it->position()
       : it->position()/it->range().mean();
 
     int nbins = 0;
     float sumX0 = 0.;
     for (int ibin = 0; ibin < 2*layersX0Eta.theHalfNBins; ibin++) { 
       const MSLayersAtAngle & layers = layersX0Eta.theLayersData[ibin];
       const MSLayer * aLayer = layers.findLayer(*it);
       if (aLayer) { nbins++; sumX0 += getDataX0(*aLayer).x0; }
     }
     if ( nbins==0) nbins=1;
 
     float hrange= (it->range().max()-it->range().min())/2.;
     DataX0 dataX0;
     if (it->face()==GeomDetEnumerators::endcap) {
       float cot1 = it->position()/(it->range().mean()-hrange/2); 
       float cot2 = it->position()/(it->range().mean()+hrange/2); 
       const MSLayer * aLayer1 = layersX0Eta.layers(cot1).findLayer(*it);
       const MSLayer * aLayer2 = layersX0Eta.layers(cot1).findLayer(*it);
       float sum1 = aLayer1 ? aLayer1->sumX0D(cot1) : 0.;
       float sum2 = aLayer2 ? aLayer2->sumX0D(cot2) : 0.;
       float slope = (sum2-sum1)/(1/cot2-1/cot1);
       float sumX0D = sum1 + slope*(1/cotTheta-1/cot1);
       dataX0 = DataX0(sumX0/nbins, sumX0D, cotTheta);
       dataX0.setForwardSumX0DSlope(slope);
     } else {
       float sumX0D = 0;
       int nb=10;
       for (int i=0; i<nb; i++) {
         float cot = (it->range().mean()+(2*i+1-nb)*hrange/nb)/it->position();
         float sin = 1/sqrt(1+cot*cot);
         const MSLayer * aLayer = layersX0Eta.layers(cot).findLayer(*it);
         if (aLayer) sumX0D += aLayer->sumX0D(cot) * sqrt(sin);
       } 
       dataX0 = DataX0(sumX0/nbins, sumX0D/nb, 0);
     }
     setDataX0(*it, dataX0);
   }
   theLayersData = MSLayersAtAngle(allLayers);
 //  cout << "MSLayersKeeperX0Averaged - LAYERS:"<<endl;
 //  theLayersData.print();
 //  cout << "END OF LAYERS"<<endl;
 }