CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Functions
DeDxTools Namespace Reference

Functions

void buildDiscrimMap (edm::Run const &run, const edm::EventSetup &iSetup, std::string Reccord, std::string ProbabilityMode, TH3F *&Prob_ChargePath)
 
int getCharge (const SiStripCluster *cluster, int &nSatStrip, const GeomDetUnit &detUnit, const std::vector< std::vector< float > > &calibGains, const unsigned int &m_off)
 
bool IsFarFromBorder (const TrajectoryStateOnSurface &trajState, const GeomDetUnit *it)
 
bool IsSpanningOver2APV (unsigned int FirstStrip, unsigned int ClusterSize)
 
void makeCalibrationMap (const std::string &m_calibrationPath, const TrackerGeometry &tkGeom, std::vector< std::vector< float > > &calibGains, const unsigned int &m_off)
 
bool shapeSelection (const SiStripCluster &ampls)
 

Function Documentation

void DeDxTools::buildDiscrimMap ( edm::Run const &  run,
const edm::EventSetup iSetup,
std::string  Reccord,
std::string  ProbabilityMode,
TH3F *&  Prob_ChargePath 
)

Definition at line 213 of file DeDxTools.cc.

References edm::hlt::Exception, edm::EventSetup::get(), i, j, roll_playback::k, prof2calltree::l, tmp, SiStripMonitorClusterAlca_cfi::xmax, SiStripMonitorClusterAlca_cfi::xmin, SiStripMonitorClusterAlca_cfi::ymax, SiStripMonitorClusterAlca_cfi::ymin, SiStripMonitorClusterAlca_cfi::zmax, and SiStripMonitorClusterAlca_cfi::zmin.

Referenced by ASmirnovDeDxDiscriminator::beginRun(), SmirnovDeDxDiscriminator::beginRun(), ProductDeDxDiscriminator::beginRun(), and BTagLikeDeDxDiscriminator::beginRun().

214 {
216  if( strcmp(Reccord.c_str(),"SiStripDeDxMip_3D_Rcd")==0){
217  iSetup.get<SiStripDeDxMip_3D_Rcd>() .get(deDxMapHandle);
218  }else if(strcmp(Reccord.c_str(),"SiStripDeDxPion_3D_Rcd")==0){
219  iSetup.get<SiStripDeDxPion_3D_Rcd>().get(deDxMapHandle);
220  }else if(strcmp(Reccord.c_str(),"SiStripDeDxKaon_3D_Rcd")==0){
221  iSetup.get<SiStripDeDxKaon_3D_Rcd>().get(deDxMapHandle);
222  }else if(strcmp(Reccord.c_str(),"SiStripDeDxProton_3D_Rcd")==0){
223  iSetup.get<SiStripDeDxProton_3D_Rcd>().get(deDxMapHandle);
224  }else if(strcmp(Reccord.c_str(),"SiStripDeDxElectron_3D_Rcd")==0){
225  iSetup.get<SiStripDeDxElectron_3D_Rcd>().get(deDxMapHandle);
226  }else{
227  throw cms::Exception("WrongReccord for dEdx") << "The reccord : " << Reccord << "is unknown\n";
228  }
229 
230  float xmin = deDxMapHandle->rangeX().min;
231  float xmax = deDxMapHandle->rangeX().max;
232  float ymin = deDxMapHandle->rangeY().min;
233  float ymax = deDxMapHandle->rangeY().max;
234  float zmin = deDxMapHandle->rangeZ().min;
235  float zmax = deDxMapHandle->rangeZ().max;
236 
237  if(Prob_ChargePath)delete Prob_ChargePath;
238  Prob_ChargePath = new TH3F ("Prob_ChargePath" , "Prob_ChargePath" , deDxMapHandle->numberOfBinsX(), xmin, xmax, deDxMapHandle->numberOfBinsY() , ymin, ymax, deDxMapHandle->numberOfBinsZ(), zmin, zmax);
239 
240  if(strcmp(ProbabilityMode.c_str(),"Accumulation")==0){
241  for(int i=0;i<=Prob_ChargePath->GetXaxis()->GetNbins()+1;i++){
242  for(int j=0;j<=Prob_ChargePath->GetYaxis()->GetNbins()+1;j++){
243  float Ni = 0;
244  for(int k=0;k<=Prob_ChargePath->GetZaxis()->GetNbins()+1;k++){ Ni+=deDxMapHandle->binContent(i,j,k);}
245  for(int k=0;k<=Prob_ChargePath->GetZaxis()->GetNbins()+1;k++){
246  float tmp = 0;
247  for(int l=0;l<=k;l++){ tmp+=deDxMapHandle->binContent(i,j,l);}
248  if(Ni>0){
249  Prob_ChargePath->SetBinContent (i, j, k, tmp/Ni);
250  }else{
251  Prob_ChargePath->SetBinContent (i, j, k, 0);
252  }
253  }
254  }
255  }
256  }else if(strcmp(ProbabilityMode.c_str(),"Integral")==0){
257  for(int i=0;i<=Prob_ChargePath->GetXaxis()->GetNbins()+1;i++){
258  for(int j=0;j<=Prob_ChargePath->GetYaxis()->GetNbins()+1;j++){
259  float Ni = 0;
260  for(int k=0;k<=Prob_ChargePath->GetZaxis()->GetNbins()+1;k++){ Ni+=deDxMapHandle->binContent(i,j,k);}
261  for(int k=0;k<=Prob_ChargePath->GetZaxis()->GetNbins()+1;k++){
262  float tmp = deDxMapHandle->binContent(i,j,k);
263  if(Ni>0){
264  Prob_ChargePath->SetBinContent (i, j, k, tmp/Ni);
265  }else{
266  Prob_ChargePath->SetBinContent (i, j, k, 0);
267  }
268  }
269  }
270  }
271  }else{
272  throw cms::Exception("WrongConfig for dEdx") << "The ProbabilityMode: " << ProbabilityMode << "is unknown\n";
273  }
274 }
int i
Definition: DBlmapReader.cc:9
int j
Definition: DBlmapReader.cc:9
const T & get() const
Definition: EventSetup.h:55
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
int DeDxTools::getCharge ( const SiStripCluster cluster,
int &  nSatStrip,
const GeomDetUnit detUnit,
const std::vector< std::vector< float > > &  calibGains,
const unsigned int &  m_off 
)

Definition at line 161 of file DeDxTools.cc.

References SiStripCluster::amplitudes(), SiStripCluster::firstStrip(), i, and GeomDet::index().

Referenced by DeDxDiscriminatorLearner::processHit(), HSCPDeDxInfoProducer::processHit(), DeDxEstimatorProducer::processHit(), and TrackFitter::run().

162 {
163  const auto & Ampls = cluster->amplitudes();
164 
165  nSatStrip = 0;
166  int charge = 0;
167  for(unsigned int i=0;i<Ampls.size();i++){
168  int calibratedCharge = Ampls[i];
169 
170  if(calibGains.size()!=0){
171  auto & gains = calibGains[detUnit.index()-m_off];
172  calibratedCharge = (int)(calibratedCharge / gains[(cluster->firstStrip()+i)/128] );
173  if(calibratedCharge>=1024){
174  calibratedCharge = 255;
175  }else if(calibratedCharge>=255){
176  calibratedCharge = 254;
177  }
178  }
179 
180  charge+=calibratedCharge;
181  if(calibratedCharge>=254)nSatStrip++;
182  }
183 
184  return charge;
185 }
int i
Definition: DBlmapReader.cc:9
uint16_t firstStrip() const
int index() const
Definition: GeomDet.h:97
const std::vector< uint8_t > & amplitudes() const
bool DeDxTools::IsFarFromBorder ( const TrajectoryStateOnSurface trajState,
const GeomDetUnit it 
)

Definition at line 303 of file DeDxTools.cc.

References Surface::bounds(), Bounds::length(), TrajectoryStateOnSurface::localError(), TrajectoryStateOnSurface::localPosition(), Parameters::parameters, LocalTrajectoryError::positionError(), GeomDet::surface(), and LocalError::yy().

Referenced by DeDxDiscriminatorLearner::processHit().

304 {
305  if (dynamic_cast<const StripGeomDetUnit*>(it)==0 && dynamic_cast<const PixelGeomDetUnit*>(it)==0) {
306  edm::LogInfo("DeDxTools::IsFarFromBorder") << "this detID doesn't seem to belong to the Tracker" << std::endl;
307  return false;
308  }
309 
310  LocalPoint HitLocalPos = trajState.localPosition();
311  LocalError HitLocalError = trajState.localError().positionError() ;
312 
313  const BoundPlane plane = it->surface();
314  const TrapezoidalPlaneBounds* trapezoidalBounds( dynamic_cast<const TrapezoidalPlaneBounds*>(&(plane.bounds())));
315  const RectangularPlaneBounds* rectangularBounds( dynamic_cast<const RectangularPlaneBounds*>(&(plane.bounds())));
316 
317  double DistFromBorder = 1.0;
318  //double HalfWidth = it->surface().bounds().width() /2.0;
319  double HalfLength = it->surface().bounds().length() /2.0;
320 
321  if(trapezoidalBounds)
322  {
323  std::array<const float, 4> const & parameters = (*trapezoidalBounds).parameters();
324  HalfLength = parameters[3];
325  //double t = (HalfLength + HitLocalPos.y()) / (2*HalfLength) ;
326  //HalfWidth = parameters[0] + (parameters[1]-parameters[0]) * t;
327  }else if(rectangularBounds){
328  //HalfWidth = it->surface().bounds().width() /2.0;
329  HalfLength = it->surface().bounds().length() /2.0;
330  }else{return false;}
331 
332 // if (fabs(HitLocalPos.x())+HitLocalError.xx() >= (HalfWidth - DistFromBorder) ) return false;//Don't think is really necessary
333  if (fabs(HitLocalPos.y())+HitLocalError.yy() >= (HalfLength - DistFromBorder) ) return false;
334 
335  return true;
336 }
dictionary parameters
Definition: Parameters.py:2
virtual float length() const =0
T y() const
Definition: PV3DBase.h:63
const Bounds & bounds() const
Definition: Surface.h:128
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:40
LocalError positionError() const
float yy() const
Definition: LocalError.h:26
const LocalTrajectoryError & localError() const
bool DeDxTools::IsSpanningOver2APV ( unsigned int  FirstStrip,
unsigned int  ClusterSize 
)

Definition at line 278 of file DeDxTools.cc.

Referenced by DeDxDiscriminatorLearner::processHit().

279 {
280  if(FirstStrip==0 ) return true;
281  if(FirstStrip==128 ) return true;
282  if(FirstStrip==256 ) return true;
283  if(FirstStrip==384 ) return true;
284  if(FirstStrip==512 ) return true;
285  if(FirstStrip==640 ) return true;
286 
287  if(FirstStrip<=127 && FirstStrip+ClusterSize>127) return true;
288  if(FirstStrip<=255 && FirstStrip+ClusterSize>255) return true;
289  if(FirstStrip<=383 && FirstStrip+ClusterSize>383) return true;
290  if(FirstStrip<=511 && FirstStrip+ClusterSize>511) return true;
291  if(FirstStrip<=639 && FirstStrip+ClusterSize>639) return true;
292 
293  if(FirstStrip+ClusterSize==127 ) return true;
294  if(FirstStrip+ClusterSize==255 ) return true;
295  if(FirstStrip+ClusterSize==383 ) return true;
296  if(FirstStrip+ClusterSize==511 ) return true;
297  if(FirstStrip+ClusterSize==639 ) return true;
298  if(FirstStrip+ClusterSize==767 ) return true;
299 
300  return false;
301 }
void DeDxTools::makeCalibrationMap ( const std::string &  m_calibrationPath,
const TrackerGeometry tkGeom,
std::vector< std::vector< float > > &  calibGains,
const unsigned int &  m_off 
)

Definition at line 188 of file DeDxTools.cc.

References TrackerGeometry::detUnits(), TrackerGeometry::idToDetUnit(), and GeomDet::index().

Referenced by DeDxDiscriminatorLearner::algoBeginJob(), HSCPDeDxInfoProducer::beginRun(), and DeDxEstimatorProducer::beginRun().

189 {
190  auto const & dus = tkGeom.detUnits();
191  calibGains.resize(dus.size()-m_off);
192 
193  TChain* t1 = new TChain("SiStripCalib/APVGain");
194  t1->Add(m_calibrationPath.c_str());
195 
196  unsigned int tree_DetId;
197  unsigned char tree_APVId;
198  double tree_Gain;
199  t1->SetBranchAddress("DetId" ,&tree_DetId );
200  t1->SetBranchAddress("APVId" ,&tree_APVId );
201  t1->SetBranchAddress("Gain" ,&tree_Gain );
202 
203  for (unsigned int ientry = 0; ientry < t1->GetEntries(); ientry++) {
204  t1->GetEntry(ientry);
205  auto& gains = calibGains[tkGeom.idToDetUnit(DetId(tree_DetId))->index()-m_off];
206  if(gains.size()<(size_t)(tree_APVId+1)){gains.resize(tree_APVId+1);}
207  gains[tree_APVId] = tree_Gain;
208  }
209  t1->Delete();
210 }
virtual const TrackerGeomDet * idToDetUnit(DetId) const
Return the pointer to the GeomDetUnit corresponding to a given DetId.
virtual const DetUnitContainer & detUnits() const
Returm a vector of all GeomDetUnit.
int index() const
Definition: GeomDet.h:97
Definition: DetId.h:18
bool DeDxTools::shapeSelection ( const SiStripCluster ampls)

Definition at line 10 of file DeDxTools.cc.

References SiStripCluster::amplitudes(), and i.

Referenced by DeDxDiscriminatorLearner::processHit(), and DeDxEstimatorProducer::processHit().

11 {
12  // ---------------- COMPTAGE DU NOMBRE DE MAXIMA --------------------------
13  //----------------------------------------------------------------------------
14 // printf("ShapeTest \n");
15  auto const & ampls = clus.amplitudes();
16  Int_t NofMax=0; Int_t recur255=1; Int_t recur254=1;
17  bool MaxOnStart=false;bool MaxInMiddle=false, MaxOnEnd =false;
18  Int_t MaxPos=0;
19  // Début avec max
20  if(ampls.size()!=1 && ((ampls[0]>ampls[1])
21  || (ampls.size()>2 && ampls[0]==ampls[1] && ampls[1]>ampls[2] && ampls[0]!=254 && ampls[0]!=255)
22  || (ampls.size()==2 && ampls[0]==ampls[1] && ampls[0]!=254 && ampls[0]!=255)) ){
23  NofMax=NofMax+1; MaxOnStart=true; }
24 
25  // Maximum entouré
26  if(ampls.size()>2){
27  for (unsigned int i =1; i < ampls.size()-1U; i++) {
28  if( (ampls[i]>ampls[i-1] && ampls[i]>ampls[i+1])
29  || (ampls.size()>3 && i>0 && i<ampls.size()-2U && ampls[i]==ampls[i+1] && ampls[i]>ampls[i-1] && ampls[i]>ampls[i+2] && ampls[i]!=254 && ampls[i]!=255) ){
30  NofMax=NofMax+1; MaxInMiddle=true; MaxPos=i;
31  }
32  if(ampls[i]==255 && ampls[i]==ampls[i-1]) {
33  recur255=recur255+1;
34  MaxPos=i-(recur255/2);
35  if(ampls[i]>ampls[i+1]){NofMax=NofMax+1;MaxInMiddle=true;}
36  }
37  if(ampls[i]==254 && ampls[i]==ampls[i-1]) {
38  recur254=recur254+1;
39  MaxPos=i-(recur254/2);
40  if(ampls[i]>ampls[i+1]){NofMax=NofMax+1;MaxInMiddle=true;}
41  }
42  }
43  }
44  // Fin avec un max
45  if(ampls.size()>1){
46  if(ampls[ampls.size()-1]>ampls[ampls.size()-2]
47  || (ampls.size()>2 && ampls[ampls.size()-1]==ampls[ampls.size()-2] && ampls[ampls.size()-2]>ampls[ampls.size()-3] )
48  || ampls[ampls.size()-1]==255){
49  NofMax=NofMax+1; MaxOnEnd=true; }
50  }
51  // Si une seule strip touchée
52  if(ampls.size()==1){ NofMax=1;}
53 
54 
55 
56 
57  // --- SELECTION EN FONCTION DE LA FORME POUR LES MAXIMA UNIQUES ---------
58  //------------------------------------------------------------------------
59  /*
60  ____
61  | |____
62  ____| | |
63  | | | |____
64  ____| | | | |
65  | | | | | |____
66  __|____|____|____|____|____|____|__
67  C_Mnn C_Mn C_M C_D C_Dn C_Dnn
68  */
69 // bool shapetest=true;
70  bool shapecdtn=false;
71 
72 // Float_t C_M; Float_t C_D; Float_t C_Mn; Float_t C_Dn; Float_t C_Mnn; Float_t C_Dnn;
73  Float_t C_M=0.0; Float_t C_D=0.0; Float_t C_Mn=10000; Float_t C_Dn=10000; Float_t C_Mnn=10000; Float_t C_Dnn=10000;
74  Int_t CDPos;
75  Float_t coeff1=1.7; Float_t coeff2=2.0;
76  Float_t coeffn=0.10; Float_t coeffnn=0.02; Float_t noise=4.0;
77 
78  if(NofMax==1){
79 
80  if(MaxOnStart==true){
81  C_M=(Float_t)ampls[0]; C_D=(Float_t)ampls[1];
82  if(ampls.size()<3) shapecdtn=true ;
83  else if(ampls.size()==3){C_Dn=(Float_t)ampls[2] ; if(C_Dn<=coeff1*coeffn*C_D+coeff2*coeffnn*C_M+2*noise || C_D==255) shapecdtn=true;}
84  else if(ampls.size()>3){ C_Dn=(Float_t)ampls[2]; C_Dnn=(Float_t)ampls[3] ;
85  if((C_Dn<=coeff1*coeffn*C_D+coeff2*coeffnn*C_M+2*noise || C_D==255)
86  && C_Dnn<=coeff1*coeffn*C_Dn+coeff2*coeffnn*C_D+2*noise){
87  shapecdtn=true;}
88  }
89  }
90 
91  if(MaxOnEnd==true){
92  C_M=(Float_t)ampls[ampls.size()-1]; C_D=(Float_t)ampls[ampls.size()-2];
93  if(ampls.size()<3) shapecdtn=true ;
94  else if(ampls.size()==3){C_Dn=(Float_t)ampls[0] ; if(C_Dn<=coeff1*coeffn*C_D+coeff2*coeffnn*C_M+2*noise || C_D==255) shapecdtn=true;}
95  else if(ampls.size()>3){C_Dn=(Float_t)ampls[ampls.size()-3] ; C_Dnn=(Float_t)ampls[ampls.size()-4] ;
96  if((C_Dn<=coeff1*coeffn*C_D+coeff2*coeffnn*C_M+2*noise || C_D==255)
97  && C_Dnn<=coeff1*coeffn*C_Dn+coeff2*coeffnn*C_D+2*noise){
98  shapecdtn=true;}
99  }
100  }
101 
102  if(MaxInMiddle==true){
103  C_M=(Float_t)ampls[MaxPos];
104  int LeftOfMaxPos=MaxPos-1;if(LeftOfMaxPos<=0)LeftOfMaxPos=0;
105  int RightOfMaxPos=MaxPos+1;if(RightOfMaxPos>=(int)ampls.size())RightOfMaxPos=ampls.size()-1;
106  //int after = RightOfMaxPos; int before = LeftOfMaxPos; if (after>=(int)ampls.size() || before<0) std::cout<<"invalid read MaxPos:"<<MaxPos <<"size:"<<ampls.size() <<std::endl;
107  if(ampls[LeftOfMaxPos]<ampls[RightOfMaxPos]){ C_D=(Float_t)ampls[RightOfMaxPos]; C_Mn=(Float_t)ampls[LeftOfMaxPos];CDPos=RightOfMaxPos;} else{ C_D=(Float_t)ampls[LeftOfMaxPos]; C_Mn=(Float_t)ampls[RightOfMaxPos];CDPos=LeftOfMaxPos;}
108  if(C_Mn<coeff1*coeffn*C_M+coeff2*coeffnn*C_D+2*noise || C_M==255){
109  if(ampls.size()==3) shapecdtn=true ;
110  else if(ampls.size()>3){
111  if(CDPos>MaxPos){
112  if(ampls.size()-CDPos-1==0){
113  C_Dn=0; C_Dnn=0;
114  }
115  if(ampls.size()-CDPos-1==1){
116  C_Dn=(Float_t)ampls[CDPos+1];
117  C_Dnn=0;
118  }
119  if(ampls.size()-CDPos-1>1){
120  C_Dn=(Float_t)ampls[CDPos+1];
121  C_Dnn=(Float_t)ampls[CDPos+2];
122  }
123  if(MaxPos>=2){
124  C_Mnn=(Float_t)ampls[MaxPos-2];
125  }
126  else if(MaxPos<2) C_Mnn=0;
127  }
128  if(CDPos<MaxPos){
129  if(CDPos==0){
130  C_Dn=0; C_Dnn=0;
131  }
132  if(CDPos==1){
133  C_Dn=(Float_t)ampls[0];
134  C_Dnn=0;
135  }
136  if(CDPos>1){
137  C_Dn=(Float_t)ampls[CDPos-1];
138  C_Dnn=(Float_t)ampls[CDPos-2];
139  }
140  if(ampls.size()-LeftOfMaxPos>1 && MaxPos+2<(int)(ampls.size())-1){
141  C_Mnn=(Float_t)ampls[MaxPos+2];
142  }else C_Mnn=0;
143  }
144  if((C_Dn<=coeff1*coeffn*C_D+coeff2*coeffnn*C_M+2*noise || C_D==255)
145  && C_Mnn<=coeff1*coeffn*C_Mn+coeff2*coeffnn*C_M+2*noise
146  && C_Dnn<=coeff1*coeffn*C_Dn+coeff2*coeffnn*C_D+2*noise) {
147  shapecdtn=true;
148  }
149 
150  }
151  }
152  }
153  }
154  if(ampls.size()==1){shapecdtn=true;}
155 
156  return shapecdtn;
157 }
int i
Definition: DBlmapReader.cc:9
const std::vector< uint8_t > & amplitudes() const