CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
SiPixelGaussianSmearingRecHitConverterAlgorithm Class Reference

#include <SiPixelGaussianSmearingRecHitConverterAlgorithm.h>

Public Member Functions

LocalError getError ()
 
double getErrorX ()
 
double getErrorY ()
 
double getErrorZ ()
 
unsigned int getPixelMultiplicityAlpha ()
 
unsigned int getPixelMultiplicityBeta ()
 
Local3DPoint getPosition ()
 
double getPositionX ()
 
double getPositionY ()
 
double getPositionZ ()
 
 SiPixelGaussianSmearingRecHitConverterAlgorithm (const edm::ParameterSet &pset, GeomDetType::SubDetector pixelPart, const RandomEngine *engine)
 
void smearHit (const PSimHit &simHit, const PixelGeomDetUnit *detUnit, const double boundX, const double boundY)
 
virtual ~SiPixelGaussianSmearingRecHitConverterAlgorithm ()
 

Private Member Functions

void initializeBarrel ()
 
void initializeForward ()
 
bool isFlipped (const PixelGeomDetUnit *theDet) const
 

Private Attributes

bool isForward
 
edm::ParameterSet pset_
 
const RandomEnginerandom
 
double rescotAlpha_binMin
 
unsigned int rescotAlpha_binN
 
double rescotAlpha_binWidth
 
double rescotBeta_binMin
 
unsigned int rescotBeta_binN
 
double rescotBeta_binWidth
 
int resqbin_binMin
 
unsigned int resqbin_binN
 
int resqbin_binWidth
 
int tempId
 
SiPixelTemplate templ
 
unsigned int theClslenx
 
unsigned int theClsleny
 
LocalError theError
 
double theErrorX
 
double theErrorY
 
double theErrorZ
 
unsigned int theLayer
 
GeomDetType::SubDetector thePixelPart
 
TFile * thePixelResolutionFile1
 
TFile * thePixelResolutionFile2
 
std::string thePixelResolutionFileName1
 
std::string thePixelResolutionFileName2
 
Local3DPoint thePosition
 
double thePositionX
 
double thePositionY
 
double thePositionZ
 
std::map< unsigned, const
SimpleHistogramGenerator * > 
theXHistos
 
std::map< unsigned, const
SimpleHistogramGenerator * > 
theYHistos
 
bool useCMSSWPixelParameterization
 

Detailed Description

Definition at line 36 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Constructor & Destructor Documentation

SiPixelGaussianSmearingRecHitConverterAlgorithm::SiPixelGaussianSmearingRecHitConverterAlgorithm ( const edm::ParameterSet pset,
GeomDetType::SubDetector  pixelPart,
const RandomEngine engine 
)
explicit

Definition at line 41 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.cc.

References gather_cfg::cout, edm::hlt::Exception, edm::ParameterSet::getParameter(), initializeBarrel(), initializeForward(), isForward, GeomDetEnumerators::PixelBarrel, GeomDetEnumerators::PixelEndcap, pset_, SiPixelTemplate::pushfile(), tempId, templ, thePixelPart, thePixelResolutionFile1, thePixelResolutionFile2, thePixelResolutionFileName1, thePixelResolutionFileName2, theXHistos, theYHistos, and useCMSSWPixelParameterization.

45 :
46  pset_(pset),
47  thePixelPart(pixelPart),
48  random(engine)
49 {
50  // Switch between old (ORCA) and new (CMSSW) pixel parameterization
51  useCMSSWPixelParameterization = pset.getParameter<bool>("UseCMSSWPixelParametrization");
52 
54  isForward = false;
55  thePixelResolutionFileName1 = pset_.getParameter<string>( "NewPixelBarrelResolutionFile1" );
57  TFile( edm::FileInPath( thePixelResolutionFileName1 ).fullPath().c_str() ,"READ");
58  thePixelResolutionFileName2 = pset_.getParameter<string>( "NewPixelBarrelResolutionFile2" );
60  TFile( edm::FileInPath( thePixelResolutionFileName2 ).fullPath().c_str() ,"READ");
62  tempId = pset_.getParameter<int> ( "templateIdBarrel" );
63  if( ! templ.pushfile(tempId) )
64  throw cms::Exception("SiPixelGaussianSmearingRecHitConverterAlgorithm:")
65  <<"SiPixel Barrel Template Not Loaded Correctly!"<<endl;
66 #ifdef FAMOS_DEBUG
67  cout<<"The Barrel map size is "<<theXHistos.size()<<" and "<<theYHistos.size()<<endl;
68 #endif
69  }
70  else
72  isForward = true;
73  thePixelResolutionFileName1 = pset_.getParameter<string>( "NewPixelForwardResolutionFile" );
74  thePixelResolutionFile1 = new
75  TFile( edm::FileInPath( thePixelResolutionFileName1 ).fullPath().c_str() ,"READ");
77  tempId = pset_.getParameter<int> ( "templateIdForward" );
78  if( ! templ.pushfile(tempId) )
79  throw cms::Exception("SiPixelGaussianSmearingRecHitConverterAlgorithm:")
80  <<"SiPixel Forward Template Not Loaded Correctly!"<<endl;
81 #ifdef FAMOS_DEBUG
82  cout<<"The Forward map size is "<<theXHistos.size()<<" and "<<theYHistos.size()<<endl;
83 #endif
84  }
85  else
86  throw cms::Exception("SiPixelGaussianSmearingRecHitConverterAlgorithm :")
87  <<"Not a pixel detector"<<endl;
88 }
T getParameter(std::string const &) const
std::map< unsigned, const SimpleHistogramGenerator * > theXHistos
bool pushfile(int filenum)
std::map< unsigned, const SimpleHistogramGenerator * > theYHistos
tuple cout
Definition: gather_cfg.py:121
SiPixelGaussianSmearingRecHitConverterAlgorithm::~SiPixelGaussianSmearingRecHitConverterAlgorithm ( )
virtual

Definition at line 90 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.cc.

References theXHistos, and theYHistos.

91 {
92 
93  std::map<unsigned,const SimpleHistogramGenerator*>::const_iterator it;
94  for ( it=theXHistos.begin(); it!=theXHistos.end(); ++it )
95  delete it->second;
96  for ( it=theYHistos.begin(); it!=theYHistos.end(); ++it )
97  delete it->second;
98 
99  theXHistos.clear();
100  theYHistos.clear();
101 
102 }
std::map< unsigned, const SimpleHistogramGenerator * > theXHistos
std::map< unsigned, const SimpleHistogramGenerator * > theYHistos

Member Function Documentation

LocalError SiPixelGaussianSmearingRecHitConverterAlgorithm::getError ( )
inline
double SiPixelGaussianSmearingRecHitConverterAlgorithm::getErrorX ( )
inline
double SiPixelGaussianSmearingRecHitConverterAlgorithm::getErrorY ( )
inline
double SiPixelGaussianSmearingRecHitConverterAlgorithm::getErrorZ ( )
inline
unsigned int SiPixelGaussianSmearingRecHitConverterAlgorithm::getPixelMultiplicityAlpha ( )
inline
unsigned int SiPixelGaussianSmearingRecHitConverterAlgorithm::getPixelMultiplicityBeta ( )
inline
Local3DPoint SiPixelGaussianSmearingRecHitConverterAlgorithm::getPosition ( )
inline
double SiPixelGaussianSmearingRecHitConverterAlgorithm::getPositionX ( )
inline
double SiPixelGaussianSmearingRecHitConverterAlgorithm::getPositionY ( )
inline
double SiPixelGaussianSmearingRecHitConverterAlgorithm::getPositionZ ( )
inline
void SiPixelGaussianSmearingRecHitConverterAlgorithm::initializeBarrel ( )
private

Definition at line 487 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.cc.

References random, rescotAlpha_binMin, rescotAlpha_binN, rescotAlpha_binWidth, rescotBeta_binMin, rescotBeta_binN, rescotBeta_binWidth, resqbin_binMin, resqbin_binN, resqbin_binWidth, thePixelResolutionFile1, thePixelResolutionFile2, theXHistos, and theYHistos.

Referenced by SiPixelGaussianSmearingRecHitConverterAlgorithm().

488 {
489  //Hard coded at the moment, can easily be changed to be configurable
490  rescotAlpha_binMin = -0.2;
491  rescotAlpha_binWidth = 0.08 ;
492  rescotAlpha_binN = 5;
493  rescotBeta_binMin = -5.5;
494  rescotBeta_binWidth = 1.0;
495  rescotBeta_binN = 11;
496  resqbin_binMin = 0;
497  resqbin_binWidth = 1;
498  resqbin_binN = 4;
499 
500  // Initialize the barrel histos once and for all, and prepare the random generation
501  for ( unsigned cotalphaHistBin=1; cotalphaHistBin<=rescotAlpha_binN; ++cotalphaHistBin )
502  for ( unsigned cotbetaHistBin=1; cotbetaHistBin<=rescotBeta_binN; ++cotbetaHistBin ) {
503  unsigned int singleBigPixelHistN = 1 * 100000
504  + cotalphaHistBin * 100
505  + cotbetaHistBin ;
506  theXHistos[singleBigPixelHistN] = new SimpleHistogramGenerator(
507  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustBPIX/hx%u" , singleBigPixelHistN ) ),
508  random);
509  theYHistos[singleBigPixelHistN] = new SimpleHistogramGenerator(
510  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustBPIX/hy%u" , singleBigPixelHistN ) ),
511  random);
512  unsigned int singlePixelHistN = 1 * 100000 + 1 * 1000
513  + cotalphaHistBin * 100
514  + cotbetaHistBin ;
515  theXHistos[singlePixelHistN] = new SimpleHistogramGenerator(
516  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustBPIX/hx%u" , singlePixelHistN ) ),
517  random);
518  theYHistos[singlePixelHistN] = new SimpleHistogramGenerator(
519  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustBPIX/hy%u" , singlePixelHistN ) ),
520  random);
521  for( unsigned qbinBin=1; qbinBin<=resqbin_binN; ++qbinBin ) {
522  unsigned int edgePixelHistN = cotalphaHistBin * 1000
523  + cotbetaHistBin * 10
524  + qbinBin;
525  theXHistos[edgePixelHistN] = new SimpleHistogramGenerator(
526  (TH1F*) thePixelResolutionFile2->Get( Form( "DQMData/clustBPIX/hx0%u" ,edgePixelHistN ) ), random );
527  theYHistos[edgePixelHistN] = new SimpleHistogramGenerator(
528  (TH1F*) thePixelResolutionFile2->Get( Form( "DQMData/clustBPIX/hy0%u" ,edgePixelHistN ) ), random );
529  unsigned int multiPixelBigHistN = 1 * 1000000 + 1 * 100000
530  + cotalphaHistBin * 1000
531  + cotbetaHistBin * 10
532  + qbinBin;
533  theXHistos[multiPixelBigHistN] = new SimpleHistogramGenerator(
534  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustBPIX/hx%u" ,multiPixelBigHistN ) ), random );
535  theYHistos[multiPixelBigHistN] = new SimpleHistogramGenerator(
536  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustBPIX/hy%u" ,multiPixelBigHistN ) ), random );
537  unsigned int multiPixelHistN = 1 * 1000000 + 1 * 100000 + 1 * 10000
538  + cotalphaHistBin * 1000
539  + cotbetaHistBin * 10
540  + qbinBin;
541  theXHistos[multiPixelHistN] = new SimpleHistogramGenerator(
542  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustBPIX/hx%u" , multiPixelHistN ) ),
543  random);
544  theYHistos[multiPixelHistN] = new SimpleHistogramGenerator(
545  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustBPIX/hy%u" , multiPixelHistN ) ),
546  random);
547  } //end for qbinBin
548  }//end for cotalphaHistBin, cotbetaHistBin
549 }
std::map< unsigned, const SimpleHistogramGenerator * > theXHistos
std::map< unsigned, const SimpleHistogramGenerator * > theYHistos
void SiPixelGaussianSmearingRecHitConverterAlgorithm::initializeForward ( )
private

Definition at line 551 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.cc.

References random, rescotAlpha_binMin, rescotAlpha_binN, rescotAlpha_binWidth, rescotBeta_binMin, rescotBeta_binN, rescotBeta_binWidth, resqbin_binMin, resqbin_binN, resqbin_binWidth, thePixelResolutionFile1, theXHistos, and theYHistos.

Referenced by SiPixelGaussianSmearingRecHitConverterAlgorithm().

552 {
553  //Hard coded at the moment, can easily be changed to be configurable
554  rescotAlpha_binMin = 0.1;
555  rescotAlpha_binWidth = 0.1 ;
556  rescotAlpha_binN = 4;
557  rescotBeta_binMin = 0.;
558  rescotBeta_binWidth = 0.15;
559  rescotBeta_binN = 4;
560  resqbin_binMin = 0;
561  resqbin_binWidth = 1;
562  resqbin_binN = 4;
563 
564  // Initialize the forward histos once and for all, and prepare the random generation
565  for ( unsigned cotalphaHistBin=1; cotalphaHistBin<=rescotAlpha_binN; ++cotalphaHistBin )
566  for ( unsigned cotbetaHistBin=1; cotbetaHistBin<=rescotBeta_binN; ++cotbetaHistBin )
567  for( unsigned qbinBin=1; qbinBin<=resqbin_binN; ++qbinBin ) {
568  unsigned int edgePixelHistN = cotalphaHistBin * 1000 + cotbetaHistBin * 10 + qbinBin;
569  theXHistos[edgePixelHistN] = new SimpleHistogramGenerator(
570  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustFPIX/fhx0%u" ,edgePixelHistN ) ), random );
571  theYHistos[edgePixelHistN] = new SimpleHistogramGenerator(
572  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustFPIX/fhy0%u" ,edgePixelHistN ) ), random );
573  unsigned int PixelHistN = 100000 + 10000 + cotalphaHistBin * 1000 + cotbetaHistBin * 10 + qbinBin;
574  theXHistos[PixelHistN] = new SimpleHistogramGenerator(
575  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustFPIX/fhx%u" ,PixelHistN ) ), random );
576  theYHistos[PixelHistN] = new SimpleHistogramGenerator(
577  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustFPIX/fhy%u" ,PixelHistN ) ), random );
578  }//end cotalphaHistBin, cotbetaHistBin, qbinBin
579 
580  for ( unsigned cotalphaHistBin=1; cotalphaHistBin<=rescotAlpha_binN; ++cotalphaHistBin )
581  for ( unsigned cotbetaHistBin=1; cotbetaHistBin<=rescotBeta_binN; ++cotbetaHistBin )
582  {
583  unsigned int SingleBigPixelHistN = 100000 + cotalphaHistBin * 100 + cotbetaHistBin;
584  theXHistos[SingleBigPixelHistN] = new SimpleHistogramGenerator(
585  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustFPIX/fhx%u" ,SingleBigPixelHistN ) ), random );
586  theYHistos[SingleBigPixelHistN] = new SimpleHistogramGenerator(
587  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustFPIX/fhy%u" ,SingleBigPixelHistN ) ), random );
588  unsigned int SinglePixelHistN = 100000 + 1000 + cotalphaHistBin * 100 + cotbetaHistBin;
589  theXHistos[SinglePixelHistN] = new SimpleHistogramGenerator(
590  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustFPIX/fhx%u" ,SinglePixelHistN ) ), random );
591  theYHistos[SinglePixelHistN] = new SimpleHistogramGenerator(
592  (TH1F*) thePixelResolutionFile1->Get( Form( "DQMData/clustFPIX/fhy%u" ,SinglePixelHistN ) ), random );
593 
594  }
595 }
std::map< unsigned, const SimpleHistogramGenerator * > theXHistos
std::map< unsigned, const SimpleHistogramGenerator * > theYHistos
bool SiPixelGaussianSmearingRecHitConverterAlgorithm::isFlipped ( const PixelGeomDetUnit theDet) const
private

Definition at line 478 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.cc.

References PV3DBase< T, PVType, FrameType >::perp(), GeomDet::surface(), and Surface::toGlobal().

Referenced by smearHit().

478  {
479  // Check the relative position of the local +/- z in global coordinates.
480  float tmp1 = theDet->surface().toGlobal(Local3DPoint(0.,0.,0.)).perp();
481  float tmp2 = theDet->surface().toGlobal(Local3DPoint(0.,0.,1.)).perp();
482  // std::cout << " 1: " << tmp1 << " 2: " << tmp2 << std::endl;
483  if ( tmp2<tmp1 ) return true;
484  else return false;
485 }
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:78
T perp() const
Definition: PV3DBase.h:71
Point3DBase< float, LocalTag > Local3DPoint
Definition: LocalPoint.h:9
const BoundPlane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:35
void SiPixelGaussianSmearingRecHitConverterAlgorithm::smearHit ( const PSimHit simHit,
const PixelGeomDetUnit detUnit,
const double  boundX,
const double  boundY 
)

Definition at line 104 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.cc.

References BHX, BHY, BXSIZE, BYSIZE, gather_cfg::cout, alignCSCRings::e, prof2calltree::edge, RandomEngine::flatShoot(), i, SiPixelTemplate::interpolate(), isFlipped(), isForward, PSimHit::localDirection(), PSimHit::localPosition(), microntocm, PSimHit::momentumAtEntry(), PixelTopology::ncolumns(), PixelTopology::nrows(), PixelTopology::pitch(), SiPixelTemplate::qbin_dist(), random, rescotAlpha_binMin, rescotAlpha_binN, rescotAlpha_binWidth, rescotBeta_binMin, rescotBeta_binN, rescotBeta_binWidth, SiPixelTemplate::s50(), PixelGeomDetUnit::specificTopology(), SiPixelTemplate::temperrors(), tempId, templ, theClslenx, theClsleny, theError, theErrorX, theErrorY, theErrorZ, thePixelPart, thePosition, thePositionX, thePositionY, thePositionZ, theXHistos, theYHistos, Vector3DBase< T, FrameTag >::unit(), PV2DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::x(), xsize, SiPixelTemplate::xtemp(), PV2DBase< T, PVType, FrameType >::y(), PV3DBase< T, PVType, FrameType >::y(), ysize, SiPixelTemplate::ytemp(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by SiTrackerGaussianSmearingRecHitConverter::gaussianSmearing().

109 {
110 
111 #ifdef FAMOS_DEBUG
112  std::cout << " Pixel smearing in " << thePixelPart
113  << std::endl;
114 #endif
115  //
116  // at the beginning the position is the Local Point in the local pixel module reference frame
117  // same code as in PixelCPEBase
118  LocalVector localDir = simHit.momentumAtEntry().unit();
119  float locx = localDir.x();
120  float locy = localDir.y();
121  float locz = localDir.z();
122 
123  // alpha: angle with respect to local x axis in local (x,z) plane
124  float cotalpha = locx/locz;
125  if ( isFlipped( detUnit ) ) { // &&& check for FPIX !!!
126 #ifdef FAMOS_DEBUG
127  std::cout << " isFlipped " << std::endl;
128 #endif
129  }
130  // beta: angle with respect to local y axis in local (y,z) plane
131  float cotbeta = locy/locz;
132  float sign=1.;
133  if( isForward ) {
134  if( cotbeta < 0 ) sign=-1.;
135  cotbeta = sign*cotbeta;
136  }
137 
138 
139  //
140 #ifdef FAMOS_DEBUG
141  std::cout << " Local Direction " << simHit.localDirection()
142  << " cotalpha(x) = " << cotalpha
143  << " cotbeta(y) = " << cotbeta
144  << std::endl;
145 #endif
146 
147  const PixelTopology* theSpecificTopology = &(detUnit->specificTopology());
148  RectangularPixelTopology rectPixelTopology(theSpecificTopology->nrows(),
149  theSpecificTopology->ncolumns(),
150  theSpecificTopology->pitch().first,
151  theSpecificTopology->pitch().second);
152  const int nrows = theSpecificTopology->nrows();
153  const int ncolumns = theSpecificTopology->ncolumns();
154 
155  const Local3DPoint lp = simHit.localPosition();
156  //Transform local position to measurement position
157  const MeasurementPoint mp = rectPixelTopology.measurementPosition( lp );
158  float mpy = mp.y();
159  float mpx = mp.x();
160  //Get the center of the struck pixel in measurement position
161  float pixelCenterY = 0.5 + (int)mpy;
162  float pixelCenterX = 0.5 + (int)mpx;
163 #ifdef FAMOS_DEBUG
164  cout<<"Struck pixel center at pitch units x: "<<pixelCenterX<<" y: "<<pixelCenterY<<endl;
165 #endif
166 
167  const MeasurementPoint mpCenter(pixelCenterX, pixelCenterY);
168  //Transform the center of the struck pixel back into local position
169  const Local3DPoint lpCenter = rectPixelTopology.localPosition( mpCenter );
170 #ifdef FAMOS_DEBUG
171  cout<<"Struck point at cm x: "<<lp.x()<<" y: "<<lp.y()<<endl;
172  cout<<"Struck pixel center at cm x: "<<lpCenter.x()<<" y: "<<lpCenter.y()<<endl;
173  cout<<"The boundX is "<<boundX<<" boundY is "<<boundY<<endl;
174 #endif
175 
176  //Get the relative position of struck point to the center of the struck pixel
177  float xtrk = lp.x() - lpCenter.x();
178  float ytrk = lp.y() - lpCenter.y();
179  //Pixel Y, X pitch
180  const float ysize={0.015}, xsize={0.01};
181  //Variables for SiPixelTemplate input, see SiPixelTemplate reco
182  float yhit = 20. + 8.*(ytrk/ysize);
183  float xhit = 20. + 8.*(xtrk/xsize);
184  int ybin = (int)yhit;
185  int xbin = (int)xhit;
186  float yfrac= yhit - (float)ybin;
187  float xfrac= xhit - (float)xbin;
188  //Protect againt ybin, xbin being outside of range [0-39]
189  if( ybin < 0 ) ybin = 0;
190  if( ybin > 39 ) ybin = 39;
191  if( xbin < 0 ) xbin = 0;
192  if( xbin > 39 ) xbin = 39;
193 
194  //Variables for SiPixelTemplate output
195  //qBin -- normalized pixel charge deposition
196  float qbin_frac[4];
197  //Single pixel cluster projection possibility
198  float ny1_frac, ny2_frac, nx1_frac, nx2_frac;
199  bool singlex = false, singley = false;
200  templ.interpolate(tempId, cotalpha, cotbeta);
201  templ.qbin_dist(tempId, cotalpha, cotbeta, qbin_frac, ny1_frac, ny2_frac, nx1_frac, nx2_frac );
202  int nqbin;
203 
204  double xsizeProbability = random->flatShoot();
205  double ysizeProbability = random->flatShoot();
206  bool hitbigx = rectPixelTopology.isItBigPixelInX( (int)mpx );
207  bool hitbigy = rectPixelTopology.isItBigPixelInY( (int)mpy );
208 
209  if( hitbigx )
210  if( xsizeProbability < nx2_frac ) singlex = true;
211  else singlex = false;
212  else
213  if( xsizeProbability < nx1_frac ) singlex = true;
214  else singlex = false;
215 
216  if( hitbigy )
217  if( ysizeProbability < ny2_frac ) singley = true;
218  else singley = false;
219  else
220  if( ysizeProbability < ny1_frac ) singley = true;
221  else singley = false;
222 
223 
224 
225  // random multiplicity for alpha and beta
226  double qbinProbability = random->flatShoot();
227  for(int i = 0; i<4; ++i) {
228  nqbin = i;
229  if(qbinProbability < qbin_frac[i]) break;
230  }
231 
232  //Store interpolated pixel cluster profile
233  //BYSIZE, BXSIZE, const definition from SiPixelTemplate
234  float ytempl[41][BYSIZE] = {{0}}, xtempl[41][BXSIZE] = {{0}} ;
235  templ.ytemp(0, 40, ytempl);
236  templ.xtemp(0, 40, xtempl);
237 
238  std::vector<double> ytemp(BYSIZE);
239  for( int i=0; i<BYSIZE; ++i) {
240  ytemp[i]=(1.-yfrac)*ytempl[ybin][i]+yfrac*ytempl[ybin+1][i];
241  }
242 
243  std::vector<double> xtemp(BXSIZE);
244  for(int i=0; i<BXSIZE; ++i) {
245  xtemp[i]=(1.-xfrac)*xtempl[xbin][i]+xfrac*xtempl[xbin+1][i];
246  }
247 
248  //Pixel readout threshold
249  const float qThreshold = templ.s50()*2.0;
250 
251  //Cut away pixels below readout threshold
252  //For cluster lengths calculation
253  int offsetX1=0, offsetX2=0, offsetY1=0, offsetY2=0;
254  int firstY, lastY, firstX, lastX;
255  for( firstY = 0; firstY < BYSIZE; ++firstY ) {
256  bool yCluster = ytemp[firstY] > qThreshold ;
257  if( yCluster )
258  {
259  offsetY1 = BHY -firstY;
260  break;
261  }
262  }
263  for( lastY = firstY; lastY < BYSIZE; ++lastY )
264  {
265  bool yCluster = ytemp[lastY] > qThreshold ;
266  if( !yCluster )
267  {
268  lastY = lastY - 1;
269  offsetY2 = lastY - BHY;
270  break;
271  }
272  }
273 
274  for( firstX = 0; firstX < BXSIZE; ++firstX ) {
275  bool xCluster = xtemp[firstX] > qThreshold ;
276  if( xCluster ) {
277  offsetX1 = BHX - firstX;
278  break;
279  }
280  }
281  for( lastX = firstX; lastX < BXSIZE; ++ lastX ) {
282  bool xCluster = xtemp[lastX] > qThreshold ;
283  if( !xCluster ) {
284  lastX = lastX - 1;
285  offsetX2 = lastX - BHX;
286  break;
287  }
288  }
289 
290  bool edge, edgex, edgey;
291  // bool bigx, bigy;
292  unsigned int clslenx = offsetX1 + offsetX2 + 1;
293  unsigned int clsleny = offsetY1 + offsetY2 + 1;
294 
295  theClslenx = clslenx;
296  theClsleny = clsleny;
297 
298  int firstPixelInX = (int)mpx - offsetX1 ;
299  int firstPixelInY = (int)mpy - offsetY1 ;
300  int lastPixelInX = (int)mpx + offsetX2 ;
301  int lastPixelInY = (int)mpy + offsetY2 ;
302  firstPixelInX = (firstPixelInX >= 0) ? firstPixelInX : 0 ;
303  firstPixelInY = (firstPixelInY >= 0) ? firstPixelInY : 0 ;
304  lastPixelInX = (lastPixelInX < nrows ) ? lastPixelInX : nrows-1 ;
305  lastPixelInY = (lastPixelInY < ncolumns ) ? lastPixelInY : ncolumns-1;
306 
307  edgex = rectPixelTopology.isItEdgePixelInX( firstPixelInX ) || rectPixelTopology.isItEdgePixelInX( lastPixelInX );
308  edgey = rectPixelTopology.isItEdgePixelInY( firstPixelInY ) || rectPixelTopology.isItEdgePixelInY( lastPixelInY );
309  edge = edgex || edgey;
310 
311  // bigx = rectPixelTopology.isItBigPixelInX( firstPixelInX ) || rectPixelTopology.isItBigPixelInX( lastPixelInX );
312  // bigy = rectPixelTopology.isItBigPixelInY( firstPixelInY ) || rectPixelTopology.isItBigPixelInY( lastPixelInY );
313  bool hasBigPixelInX = rectPixelTopology.containsBigPixelInX( firstPixelInX, lastPixelInX );
314  bool hasBigPixelInY = rectPixelTopology.containsBigPixelInY( firstPixelInY, lastPixelInY );
315 
316  //Variables for SiPixelTemplate pixel hit error output
317  float sigmay, sigmax, sy1, sy2, sx1, sx2;
318  templ.temperrors(tempId, cotalpha, cotbeta, nqbin, sigmay, sigmax, sy1, sy2, sx1, sx2 );
319 
320  // define private mebers --> Errors
321  if( edge ) {
322  if( edgex && !edgey ) {
323  theErrorX = 23.0*microntocm;
324  theErrorY = 39.0*microntocm;
325  }
326  else if( !edgex && edgey ) {
327  theErrorX = 24.0*microntocm;
328  theErrorY = 96.0*microntocm;
329  }
330  else
331  {
332  theErrorX = 31.0*microntocm;
333  theErrorY = 90.0*microntocm;
334  }
335 
336  }
337  else {
338  if( singlex )
339  if ( hitbigx )
340  theErrorX = sx2*microntocm;
341  else
342  theErrorX = sx1*microntocm;
343  else theErrorX = sigmax*microntocm;
344  if( singley )
345  if( hitbigy )
346  theErrorY = sy2*microntocm;
347  else
348  theErrorY = sy1*microntocm;
349  else theErrorY = sigmay*microntocm;
350  }
351  theErrorZ = 1e-8; // 1 um means zero
353  // Local Error is 2D: (xx,xy,yy), square of sigma in first an third position
354  // as for resolution matrix
355  //
356 #ifdef FAMOS_DEBUG
357  std::cout << " Pixel Errors "
358  << "\talpha(x) = " << theErrorX
359  << "\tbeta(y) = " << theErrorY
360  << std::endl;
361 #endif
362  // Generate position
363  // get resolution histograms
364  int cotalphaHistBin = (int)( ( cotalpha - rescotAlpha_binMin ) / rescotAlpha_binWidth + 1 );
365  int cotbetaHistBin = (int)( ( cotbeta - rescotBeta_binMin ) / rescotBeta_binWidth + 1 );
366  // protection against out-of-range (undeflows and overflows)
367  if( cotalphaHistBin < 1 ) cotalphaHistBin = 1;
368  if( cotbetaHistBin < 1 ) cotbetaHistBin = 1;
369  if( cotalphaHistBin > (int)rescotAlpha_binN ) cotalphaHistBin = (int)rescotAlpha_binN;
370  if( cotbetaHistBin > (int)rescotBeta_binN ) cotbetaHistBin = (int)rescotBeta_binN;
371  //
372  unsigned int theXHistN;
373  unsigned int theYHistN;
374  if( !isForward ) {
375  if(edge)
376  {
377  theXHistN = cotalphaHistBin * 1000 + cotbetaHistBin * 10 + (nqbin+1);
378  theYHistN = theXHistN;
379  }
380  else
381  {
382  if(singlex)
383  {
384  if(hitbigx) theXHistN = 1 * 100000 + cotalphaHistBin * 100 + cotbetaHistBin ;
385  else theXHistN = 1 * 100000 + 1 * 1000 + cotalphaHistBin * 100 + cotbetaHistBin ;
386  }
387  else
388  {
389  if(hasBigPixelInX) theXHistN = 1 * 1000000 + 1 * 100000 + cotalphaHistBin * 1000 + cotbetaHistBin * 10 + (nqbin+1);
390  else theXHistN = 1 * 1000000 + 1 * 100000 + 1 * 10000 + cotalphaHistBin * 1000 + cotbetaHistBin * 10 + (nqbin+1);
391  }
392  if(singley)
393  {
394  if(hitbigy) theYHistN = 1 * 100000 + cotalphaHistBin * 100 + cotbetaHistBin ;
395  else theYHistN = 1 * 100000 + 1 * 1000 + cotalphaHistBin * 100 + cotbetaHistBin ;
396  }
397  else
398  {
399  if(hasBigPixelInY) theYHistN = 1 * 1000000 + 1 * 100000 + cotalphaHistBin * 1000 + cotbetaHistBin * 10 + (nqbin+1);
400  else theYHistN = 1 * 1000000 + 1 * 100000 + 1 * 10000 + cotalphaHistBin * 1000 + cotbetaHistBin * 10 + (nqbin+1);
401  }
402  }
403  }
404  else
405  {
406  if(edge)
407  {
408  theXHistN = cotalphaHistBin * 1000 + cotbetaHistBin * 10 + (nqbin+1);
409  theYHistN = theXHistN;
410  }
411  else
412  {
413  if( singlex )
414  if( hitbigx )
415  theXHistN = 100000 + cotalphaHistBin * 100 + cotbetaHistBin;
416  else
417  theXHistN = 100000 + 1000 + cotalphaHistBin * 100 + cotbetaHistBin;
418  else
419  theXHistN = 100000 + 10000 + cotalphaHistBin * 1000 + cotbetaHistBin * 10 + (nqbin+1);
420  if( singley )
421  if( hitbigy )
422  theYHistN = 100000 + cotalphaHistBin * 100 + cotbetaHistBin;
423  else
424  theYHistN = 100000 + 1000 + cotalphaHistBin * 100 + cotbetaHistBin;
425  else
426  theYHistN = 100000 + 10000 + cotalphaHistBin * 1000 + cotbetaHistBin * 10 + (nqbin+1);
427  }
428  }
429  unsigned int counter = 0;
430  do {
431  //
432  // Smear the hit Position
433  thePositionX = theXHistos[theXHistN]->generate();
434  thePositionY = theYHistos[theYHistN]->generate();
435  if( isForward ) thePositionY *= sign;
436  thePositionZ = 0.0; // set at the centre of the active area
437  //protect from empty resolution histograms
438  //if( thePositionX > 0.0799 ) thePositionX = 0;
439  //if( thePositionY > 0.0799 ) thePositionY = 0;
440  thePosition =
441  Local3DPoint(simHit.localPosition().x() + thePositionX ,
442  simHit.localPosition().y() + thePositionY ,
443  simHit.localPosition().z() + thePositionZ );
444 #ifdef FAMOS_DEBUG
445  std::cout << " Detector bounds: "
446  << "\t\tx = " << boundX
447  << "\ty = " << boundY
448  << std::endl;
449  std::cout << " Generated local position "
450  << "\tx = " << thePosition.x()
451  << "\ty = " << thePosition.y()
452  << std::endl;
453 #endif
454  counter++;
455  if(counter > 20) {
456  thePosition = Local3DPoint(simHit.localPosition().x(),
457  simHit.localPosition().y(),
458  simHit.localPosition().z());
459  break;
460  }
461  } while(fabs(thePosition.x()) > boundX || fabs(thePosition.y()) > boundY);
462 
463 }
int i
Definition: DBlmapReader.cc:9
#define BXSIZE
T y() const
Definition: PV2DBase.h:45
#define BYSIZE
virtual int ncolumns() const =0
LocalVector momentumAtEntry() const
The momentum of the track that produced the hit, at entry point.
Definition: PSimHit.h:47
T y() const
Definition: PV3DBase.h:62
virtual int nrows() const =0
const Int_t ysize
std::map< unsigned, const SimpleHistogramGenerator * > theXHistos
bool interpolate(int id, float cotalpha, float cotbeta, float locBz)
void ytemp(int fybin, int lybin, float ytemplate[41][21+4])
Local3DPoint localPosition() const
Definition: PSimHit.h:44
#define BHX
T z() const
Definition: PV3DBase.h:63
void temperrors(int id, float cotalpha, float cotbeta, int qBin, float &sigmay, float &sigmax, float &sy1, float &sy2, float &sx1, float &sx2)
float s50()
1/2 of the pixel threshold signal in electrons
virtual std::pair< float, float > pitch() const =0
LocalVector localDirection() const
Obsolete. Same as momentumAtEntry().unit(), for backward compatibility.
Definition: PSimHit.h:52
Point3DBase< float, LocalTag > Local3DPoint
Definition: LocalPoint.h:9
Vector3DBase unit() const
Definition: Vector3DBase.h:57
void xtemp(int fxbin, int lxbin, float xtemplate[41][13+4])
void qbin_dist(int id, float cotalpha, float cotbeta, float qbin_frac[4], float &ny1_frac, float &ny2_frac, float &nx1_frac, float &nx2_frac)
std::map< unsigned, const SimpleHistogramGenerator * > theYHistos
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
double flatShoot(double xmin=0.0, double xmax=1.0) const
Definition: RandomEngine.h:30
#define BHY
tuple cout
Definition: gather_cfg.py:121
const Int_t xsize
T x() const
Definition: PV2DBase.h:44
T x() const
Definition: PV3DBase.h:61

Member Data Documentation

bool SiPixelGaussianSmearingRecHitConverterAlgorithm::isForward
private
edm::ParameterSet SiPixelGaussianSmearingRecHitConverterAlgorithm::pset_
private
const RandomEngine* SiPixelGaussianSmearingRecHitConverterAlgorithm::random
private
double SiPixelGaussianSmearingRecHitConverterAlgorithm::rescotAlpha_binMin
private
unsigned int SiPixelGaussianSmearingRecHitConverterAlgorithm::rescotAlpha_binN
private
double SiPixelGaussianSmearingRecHitConverterAlgorithm::rescotAlpha_binWidth
private
double SiPixelGaussianSmearingRecHitConverterAlgorithm::rescotBeta_binMin
private
unsigned int SiPixelGaussianSmearingRecHitConverterAlgorithm::rescotBeta_binN
private
double SiPixelGaussianSmearingRecHitConverterAlgorithm::rescotBeta_binWidth
private
int SiPixelGaussianSmearingRecHitConverterAlgorithm::resqbin_binMin
private
unsigned int SiPixelGaussianSmearingRecHitConverterAlgorithm::resqbin_binN
private
int SiPixelGaussianSmearingRecHitConverterAlgorithm::resqbin_binWidth
private
int SiPixelGaussianSmearingRecHitConverterAlgorithm::tempId
private
SiPixelTemplate SiPixelGaussianSmearingRecHitConverterAlgorithm::templ
private
unsigned int SiPixelGaussianSmearingRecHitConverterAlgorithm::theClslenx
private
unsigned int SiPixelGaussianSmearingRecHitConverterAlgorithm::theClsleny
private
LocalError SiPixelGaussianSmearingRecHitConverterAlgorithm::theError
private

Definition at line 104 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Referenced by getError(), and smearHit().

double SiPixelGaussianSmearingRecHitConverterAlgorithm::theErrorX
private

Definition at line 105 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Referenced by getErrorX(), and smearHit().

double SiPixelGaussianSmearingRecHitConverterAlgorithm::theErrorY
private

Definition at line 106 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Referenced by getErrorY(), and smearHit().

double SiPixelGaussianSmearingRecHitConverterAlgorithm::theErrorZ
private

Definition at line 107 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Referenced by getErrorZ(), and smearHit().

unsigned int SiPixelGaussianSmearingRecHitConverterAlgorithm::theLayer
private
GeomDetType::SubDetector SiPixelGaussianSmearingRecHitConverterAlgorithm::thePixelPart
private
TFile* SiPixelGaussianSmearingRecHitConverterAlgorithm::thePixelResolutionFile1
private
TFile* SiPixelGaussianSmearingRecHitConverterAlgorithm::thePixelResolutionFile2
private
std::string SiPixelGaussianSmearingRecHitConverterAlgorithm::thePixelResolutionFileName1
private
std::string SiPixelGaussianSmearingRecHitConverterAlgorithm::thePixelResolutionFileName2
private
Local3DPoint SiPixelGaussianSmearingRecHitConverterAlgorithm::thePosition
private

Definition at line 100 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Referenced by getPosition(), and smearHit().

double SiPixelGaussianSmearingRecHitConverterAlgorithm::thePositionX
private

Definition at line 101 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Referenced by getPositionX(), and smearHit().

double SiPixelGaussianSmearingRecHitConverterAlgorithm::thePositionY
private

Definition at line 102 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Referenced by getPositionY(), and smearHit().

double SiPixelGaussianSmearingRecHitConverterAlgorithm::thePositionZ
private

Definition at line 103 of file SiPixelGaussianSmearingRecHitConverterAlgorithm.h.

Referenced by getPositionZ(), and smearHit().

std::map<unsigned,const SimpleHistogramGenerator*> SiPixelGaussianSmearingRecHitConverterAlgorithm::theXHistos
private
std::map<unsigned,const SimpleHistogramGenerator*> SiPixelGaussianSmearingRecHitConverterAlgorithm::theYHistos
private
bool SiPixelGaussianSmearingRecHitConverterAlgorithm::useCMSSWPixelParameterization
private