#include <MSLayersKeeperX0Averaged.h>

Inheritance diagram for MSLayersKeeperX0Averaged:

Public Member Functions
MSLayer	layer (const DetLayer *layer) const override

const MSLayersAtAngle &	layers (float cotTheta) const override

	MSLayersKeeperX0Averaged (const GeometricSearchTracker &tracker, const MagneticField &bfield)

	~MSLayersKeeperX0Averaged () override=default

Public Member Functions inherited from MSLayersKeeper
virtual	~MSLayersKeeper ()

Private Attributes
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 10 of file MSLayersKeeperX0Averaged.h.

Constructor & Destructor Documentation

◆ MSLayersKeeperX0Averaged()

MSLayersKeeperX0Averaged::MSLayersKeeperX0Averaged	(	const GeometricSearchTracker &	tracker,
		const MagneticField &	bfield
	)

Definition at line 7 of file MSLayersKeeperX0Averaged.cc.

References GeomDetEnumerators::barrel, GeomDetEnumerators::endcap, PVValHelper::eta, spr::find(), MSLayersAtAngle::findLayer(), relativeConstraints::geom, MSLayersKeeper::getDataX0(), mps_fire::i, ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::it, hgcalTBTopologyTester_cfi::layers, MSLayersKeeperX0AtEta::layers(), LaserClient_cfi::nbins, MSLayersKeeper::setDataX0(), funct::sin(), slope, mathSSE::sqrt(), MSLayer::sumX0D(), MSLayersKeeperX0AtEta::theHalfNBins, theLayersData, MSLayersKeeperX0AtEta::theLayersData, PbPb_ZMuSkimMuonDPG_cff::tracker, and MSLayer::DataX0::x0.

                                                                                                                      {
   //  cout << "HERE INITIALISATION! MSLayersKeeperX0Averaged"<<endl;
   MSLayersKeeperX0AtEta layersX0Eta(tracker, bfield);
   MultipleScatteringGeometry geom(tracker);
   vector<MSLayer> allLayers = geom.detLayers();
   vector<MSLayer>::iterator it;
   for (int i = -1; i <= 1; i++) {
     float eta = i * (-1.8);
     vector<MSLayer> tmpLayers = geom.otherLayers(eta);
     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;
 }