test
CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Classes | Functions
DeDxTools Namespace Reference

Classes

struct  RawHits
 

Functions

double genericAverage (const reco::DeDxHitCollection &, float expo=1.)
 
const SiStripCluster * GetCluster (const TrackerSingleRecHit *hit)
 
const SiStripCluster * GetCluster (const TrackerSingleRecHit &hit)
 
bool shapeSelection (const std::vector< uint8_t > &ampls)
 
void trajectoryRawHits (const edm::Ref< std::vector< Trajectory > > &trajectory, vector< RawHits > &hits, bool usePixel, bool useStrip)
 
void trajectoryRawHits (const edm::Ref< std::vector< Trajectory > > &trajectory, std::vector< RawHits > &hits, bool usePixel, bool useStrip)
 

Function Documentation

double DeDxTools::genericAverage ( const reco::DeDxHitCollection &  hits,
float  expo = 1. 
)

Definition at line 121 of file DeDxTools.cc.

References DeDxDiscriminatorTools::charge(), i, n, funct::pow(), and query::result.

Referenced by GenericAverageDeDxEstimator::dedx().

122 {
123  double result=0;
124  size_t n = hits.size();
125  for(size_t i = 0; i< n; i ++)
126  {
127  result+=pow(hits[i].charge(),expo);
128  }
129  return (n>0)?pow(result/n,1./expo):0.;
130 }
i
int i
Definition: DBlmapReader.cc:9
DeDxDiscriminatorTools::charge
double charge(const std::vector< uint8_t > &Ampls)
Definition: DeDxDiscriminatorTools.cc:43
query.result
tuple result
Definition: query.py:137
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
const SiStripCluster* DeDxTools::GetCluster ( const TrackerSingleRecHit *  hit)
inline

Definition at line 27 of file DeDxTools.h.

References TrackerSingleRecHit::stripCluster().

Referenced by DeDxEstimatorProducer::produce(), and DeDxDiscriminatorProducer::produce().

27 { return &hit->stripCluster();}
TrackerSingleRecHit::stripCluster
SiStripCluster const & stripCluster() const
Definition: TrackerSingleRecHit.h:64
const SiStripCluster* DeDxTools::GetCluster ( const TrackerSingleRecHit &  hit)
inline

Definition at line 28 of file DeDxTools.h.

References TrackerSingleRecHit::stripCluster().

28 {return &hit.stripCluster();}
TrackerSingleRecHit::stripCluster
SiStripCluster const & stripCluster() const
Definition: TrackerSingleRecHit.h:64
bool DeDxTools::shapeSelection ( const std::vector< uint8_t > &  ampls)

Definition at line 133 of file DeDxTools.cc.

References i.

Referenced by DeDxEstimatorProducer::produce(), and DeDxDiscriminatorProducer::produce().

134 {
135  // ---------------- COMPTAGE DU NOMBRE DE MAXIMA --------------------------
136  //----------------------------------------------------------------------------
137 // printf("ShapeTest \n");
138  Int_t NofMax=0; Int_t recur255=1; Int_t recur254=1;
139  bool MaxOnStart=false;bool MaxInMiddle=false, MaxOnEnd =false;
140  Int_t MaxPos=0;
141  // Début avec max
142  if(ampls.size()!=1 && ((ampls[0]>ampls[1])
143  || (ampls.size()>2 && ampls[0]==ampls[1] && ampls[1]>ampls[2] && ampls[0]!=254 && ampls[0]!=255)
144  || (ampls.size()==2 && ampls[0]==ampls[1] && ampls[0]!=254 && ampls[0]!=255)) ){
145  NofMax=NofMax+1; MaxOnStart=true; }
146 
147  // Maximum entouré
148  if(ampls.size()>2){
149  for (unsigned int i =1; i < ampls.size()-1; i++) {
150  if( (ampls[i]>ampls[i-1] && ampls[i]>ampls[i+1])
151  || (ampls.size()>3 && i>0 && i<ampls.size()-2 && ampls[i]==ampls[i+1] && ampls[i]>ampls[i-1] && ampls[i]>ampls[i+2] && ampls[i]!=254 && ampls[i]!=255) ){
152  NofMax=NofMax+1; MaxInMiddle=true; MaxPos=i;
153  }
154  if(ampls[i]==255 && ampls[i]==ampls[i-1]) {
155  recur255=recur255+1;
156  MaxPos=i-(recur255/2);
157  if(ampls[i]>ampls[i+1]){NofMax=NofMax+1;MaxInMiddle=true;}
158  }
159  if(ampls[i]==254 && ampls[i]==ampls[i-1]) {
160  recur254=recur254+1;
161  MaxPos=i-(recur254/2);
162  if(ampls[i]>ampls[i+1]){NofMax=NofMax+1;MaxInMiddle=true;}
163  }
164  }
165  }
166  // Fin avec un max
167  if(ampls.size()>1){
168  if(ampls[ampls.size()-1]>ampls[ampls.size()-2]
169  || (ampls.size()>2 && ampls[ampls.size()-1]==ampls[ampls.size()-2] && ampls[ampls.size()-2]>ampls[ampls.size()-3] )
170  || ampls[ampls.size()-1]==255){
171  NofMax=NofMax+1; MaxOnEnd=true; }
172  }
173  // Si une seule strip touchée
174  if(ampls.size()==1){ NofMax=1;}
175 
176 
177 
178 
179  // --- SELECTION EN FONCTION DE LA FORME POUR LES MAXIMA UNIQUES ---------
180  //------------------------------------------------------------------------
181  /*
182  ____
183  | |____
184  ____| | |
185  | | | |____
186  ____| | | | |
187  | | | | | |____
188  __|____|____|____|____|____|____|__
189  C_Mnn C_Mn C_M C_D C_Dn C_Dnn
190  */
191 // bool shapetest=true;
192  bool shapecdtn=false;
193 
194 // Float_t C_M; Float_t C_D; Float_t C_Mn; Float_t C_Dn; Float_t C_Mnn; Float_t C_Dnn;
195  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;
196  Int_t CDPos;
197  Float_t coeff1=1.7; Float_t coeff2=2.0;
198  Float_t coeffn=0.10; Float_t coeffnn=0.02; Float_t noise=4.0;
199 
200  if(NofMax==1){
201 
202  if(MaxOnStart==true){
203  C_M=(Float_t)ampls[0]; C_D=(Float_t)ampls[1];
204  if(ampls.size()<3) shapecdtn=true ;
205  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;}
206  else if(ampls.size()>3){ C_Dn=(Float_t)ampls[2]; C_Dnn=(Float_t)ampls[3] ;
207  if((C_Dn<=coeff1*coeffn*C_D+coeff2*coeffnn*C_M+2*noise || C_D==255)
208  && C_Dnn<=coeff1*coeffn*C_Dn+coeff2*coeffnn*C_D+2*noise){
209  shapecdtn=true;}
210  }
211  }
212 
213  if(MaxOnEnd==true){
214  C_M=(Float_t)ampls[ampls.size()-1]; C_D=(Float_t)ampls[ampls.size()-2];
215  if(ampls.size()<3) shapecdtn=true ;
216  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;}
217  else if(ampls.size()>3){C_Dn=(Float_t)ampls[ampls.size()-3] ; C_Dnn=(Float_t)ampls[ampls.size()-4] ;
218  if((C_Dn<=coeff1*coeffn*C_D+coeff2*coeffnn*C_M+2*noise || C_D==255)
219  && C_Dnn<=coeff1*coeffn*C_Dn+coeff2*coeffnn*C_D+2*noise){
220  shapecdtn=true;}
221  }
222  }
223 
224  if(MaxInMiddle==true){
225  C_M=(Float_t)ampls[MaxPos];
226  int LeftOfMaxPos=MaxPos-1;if(LeftOfMaxPos<=0)LeftOfMaxPos=0;
227  int RightOfMaxPos=MaxPos+1;if(RightOfMaxPos>=(int)ampls.size())RightOfMaxPos=ampls.size()-1;
228  //int after = RightOfMaxPos; int before = LeftOfMaxPos; if (after>=(int)ampls.size() || before<0) std::cout<<"invalid read MaxPos:"<<MaxPos <<"size:"<<ampls.size() <<std::endl;
229  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;}
230  if(C_Mn<coeff1*coeffn*C_M+coeff2*coeffnn*C_D+2*noise || C_M==255){
231  if(ampls.size()==3) shapecdtn=true ;
232  else if(ampls.size()>3){
233  if(CDPos>MaxPos){
234  if(ampls.size()-CDPos-1==0){
235  C_Dn=0; C_Dnn=0;
236  }
237  if(ampls.size()-CDPos-1==1){
238  C_Dn=(Float_t)ampls[CDPos+1];
239  C_Dnn=0;
240  }
241  if(ampls.size()-CDPos-1>1){
242  C_Dn=(Float_t)ampls[CDPos+1];
243  C_Dnn=(Float_t)ampls[CDPos+2];
244  }
245  if(MaxPos>=2){
246  C_Mnn=(Float_t)ampls[MaxPos-2];
247  }
248  else if(MaxPos<2) C_Mnn=0;
249  }
250  if(CDPos<MaxPos){
251  if(CDPos==0){
252  C_Dn=0; C_Dnn=0;
253  }
254  if(CDPos==1){
255  C_Dn=(Float_t)ampls[0];
256  C_Dnn=0;
257  }
258  if(CDPos>1){
259  C_Dn=(Float_t)ampls[CDPos-1];
260  C_Dnn=(Float_t)ampls[CDPos-2];
261  }
262  if(ampls.size()-LeftOfMaxPos>1 && MaxPos+2<(int)(ampls.size())-1){
263  C_Mnn=(Float_t)ampls[MaxPos+2];
264  }else C_Mnn=0;
265  }
266  if((C_Dn<=coeff1*coeffn*C_D+coeff2*coeffnn*C_M+2*noise || C_D==255)
267  && C_Mnn<=coeff1*coeffn*C_Mn+coeff2*coeffnn*C_M+2*noise
268  && C_Dnn<=coeff1*coeffn*C_Dn+coeff2*coeffnn*C_D+2*noise) {
269  shapecdtn=true;
270  }
271 
272  }
273  }
274  }
275  }
276  if(ampls.size()==1){shapecdtn=true;}
277 
278  return shapecdtn;
279 }
i
int i
Definition: DBlmapReader.cc:9
void DeDxTools::trajectoryRawHits ( const edm::Ref< std::vector< Trajectory > > &  trajectory,
vector< RawHits > &  hits,
bool  usePixel,
bool  useStrip 
)

Definition at line 10 of file DeDxTools.cc.

References DeDxTools::RawHits::angleCosine, DeDxTools::RawHits::charge, DeDxTools::RawHits::detId, i, TrajectoryStateOnSurface::isValid(), TrajectoryStateOnSurface::localDirection(), DeDxTools::RawHits::NSaturating, DeDxTools::RawHits::trajectoryMeasurement, and PV3DBase< T, PVType, FrameType >::z().

11  {
12 
13  // vector<RawHits> hits;
14 
15  const vector<TrajectoryMeasurement> & measurements = trajectory->measurements();
16  for(vector<TrajectoryMeasurement>::const_iterator it = measurements.begin(); it!=measurements.end(); it++){
17 
18  //FIXME: check that "updated" State is the best one (wrt state in the middle)
19  TrajectoryStateOnSurface trajState=it->updatedState();
20  if( !trajState.isValid()) continue;
21 
22  const TrackingRecHit * recHit=(*it->recHit()).hit();
23 
24  LocalVector trackDirection = trajState.localDirection();
25  double cosine = trackDirection.z()/trackDirection.mag();
26 
27  if(const SiStripMatchedRecHit2D* matchedHit=dynamic_cast<const SiStripMatchedRecHit2D*>(recHit)){
28  if(!useStrip) continue;
29 
30  RawHits mono,stereo;
31  mono.trajectoryMeasurement = &(*it);
32  stereo.trajectoryMeasurement = &(*it);
33  mono.angleCosine = cosine;
34  stereo.angleCosine = cosine;
35  const std::vector<uint8_t> & amplitudes = matchedHit->monoCluster().amplitudes();
36  mono.charge = accumulate(amplitudes.begin(), amplitudes.end(), 0);
37  mono.NSaturating =0;
38  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)mono.NSaturating++;}
39 
40  const std::vector<uint8_t> & amplitudesSt = matchedHit->stereoCluster().amplitudes();
41  stereo.charge = accumulate(amplitudesSt.begin(), amplitudesSt.end(), 0);
42  stereo.NSaturating =0;
43  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)stereo.NSaturating++;}
44 
45  mono.detId= matchedHit->monoId();
46  stereo.detId= matchedHit->stereoId();
47 
48  hits.push_back(mono);
49  hits.push_back(stereo);
50 
51  }else if(const ProjectedSiStripRecHit2D* projectedHit=dynamic_cast<const ProjectedSiStripRecHit2D*>(recHit)) {
52  if(!useStrip) continue;
53 
54  RawHits mono;
55 
56  mono.trajectoryMeasurement = &(*it);
57 
58  mono.angleCosine = cosine;
59  const std::vector<uint8_t> & amplitudes = projectedHit->cluster()->amplitudes();
60  mono.charge = accumulate(amplitudes.begin(), amplitudes.end(), 0);
61  mono.NSaturating =0;
62  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)mono.NSaturating++;}
63 
64  mono.detId= projectedHit->originalId();
65  hits.push_back(mono);
66 
67  }else if(const SiStripRecHit2D* singleHit=dynamic_cast<const SiStripRecHit2D*>(recHit)){
68  if(!useStrip) continue;
69 
70  RawHits mono;
71 
72  mono.trajectoryMeasurement = &(*it);
73 
74  mono.angleCosine = cosine;
75  const std::vector<uint8_t> & amplitudes = singleHit->cluster()->amplitudes();
76  mono.charge = accumulate(amplitudes.begin(), amplitudes.end(), 0);
77  mono.NSaturating =0;
78  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)mono.NSaturating++;}
79 
80  mono.detId= singleHit->geographicalId();
81  hits.push_back(mono);
82 
83  }else if(const SiStripRecHit1D* single1DHit=dynamic_cast<const SiStripRecHit1D*>(recHit)){
84  if(!useStrip) continue;
85 
86  RawHits mono;
87 
88  mono.trajectoryMeasurement = &(*it);
89 
90  mono.angleCosine = cosine;
91  const std::vector<uint8_t> & amplitudes = single1DHit->cluster()->amplitudes();
92  mono.charge = accumulate(amplitudes.begin(), amplitudes.end(), 0);
93  mono.NSaturating =0;
94  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)mono.NSaturating++;}
95 
96  mono.detId= single1DHit->geographicalId();
97  hits.push_back(mono);
98 
99 
100  }else if(const SiPixelRecHit* pixelHit=dynamic_cast<const SiPixelRecHit*>(recHit)){
101  if(!usePixel) continue;
102 
103  RawHits pixel;
104 
105  pixel.trajectoryMeasurement = &(*it);
106 
107  pixel.angleCosine = cosine;
108  pixel.charge = pixelHit->cluster()->charge();
109  pixel.NSaturating=-1;
110  pixel.detId= pixelHit->geographicalId();
111  hits.push_back(pixel);
112  }
113 
114  }
115  // return hits;
116  }
i
int i
Definition: DBlmapReader.cc:9
TrajectoryStateOnSurface::localDirection
LocalVector localDirection() const
Definition: TrajectoryStateOnSurface.h:122
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:67
PV3DBase::z
T z() const
Definition: PV3DBase.h:64
SiPixelRecHit
Our base class.
Definition: SiPixelRecHit.h:23
void DeDxTools::trajectoryRawHits ( const edm::Ref< std::vector< Trajectory > > &  trajectory,
std::vector< RawHits > &  hits,
bool  usePixel,
bool  useStrip 
)

Definition at line 10 of file DeDxTools.cc.

References DeDxTools::RawHits::angleCosine, DeDxTools::RawHits::charge, DeDxTools::RawHits::detId, i, TrajectoryStateOnSurface::isValid(), TrajectoryStateOnSurface::localDirection(), DeDxTools::RawHits::NSaturating, DeDxTools::RawHits::trajectoryMeasurement, and PV3DBase< T, PVType, FrameType >::z().

11  {
12 
13  // vector<RawHits> hits;
14 
15  const vector<TrajectoryMeasurement> & measurements = trajectory->measurements();
16  for(vector<TrajectoryMeasurement>::const_iterator it = measurements.begin(); it!=measurements.end(); it++){
17 
18  //FIXME: check that "updated" State is the best one (wrt state in the middle)
19  TrajectoryStateOnSurface trajState=it->updatedState();
20  if( !trajState.isValid()) continue;
21 
22  const TrackingRecHit * recHit=(*it->recHit()).hit();
23 
24  LocalVector trackDirection = trajState.localDirection();
25  double cosine = trackDirection.z()/trackDirection.mag();
26 
27  if(const SiStripMatchedRecHit2D* matchedHit=dynamic_cast<const SiStripMatchedRecHit2D*>(recHit)){
28  if(!useStrip) continue;
29 
30  RawHits mono,stereo;
31  mono.trajectoryMeasurement = &(*it);
32  stereo.trajectoryMeasurement = &(*it);
33  mono.angleCosine = cosine;
34  stereo.angleCosine = cosine;
35  const std::vector<uint8_t> & amplitudes = matchedHit->monoCluster().amplitudes();
36  mono.charge = accumulate(amplitudes.begin(), amplitudes.end(), 0);
37  mono.NSaturating =0;
38  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)mono.NSaturating++;}
39 
40  const std::vector<uint8_t> & amplitudesSt = matchedHit->stereoCluster().amplitudes();
41  stereo.charge = accumulate(amplitudesSt.begin(), amplitudesSt.end(), 0);
42  stereo.NSaturating =0;
43  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)stereo.NSaturating++;}
44 
45  mono.detId= matchedHit->monoId();
46  stereo.detId= matchedHit->stereoId();
47 
48  hits.push_back(mono);
49  hits.push_back(stereo);
50 
51  }else if(const ProjectedSiStripRecHit2D* projectedHit=dynamic_cast<const ProjectedSiStripRecHit2D*>(recHit)) {
52  if(!useStrip) continue;
53 
54  RawHits mono;
55 
56  mono.trajectoryMeasurement = &(*it);
57 
58  mono.angleCosine = cosine;
59  const std::vector<uint8_t> & amplitudes = projectedHit->cluster()->amplitudes();
60  mono.charge = accumulate(amplitudes.begin(), amplitudes.end(), 0);
61  mono.NSaturating =0;
62  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)mono.NSaturating++;}
63 
64  mono.detId= projectedHit->originalId();
65  hits.push_back(mono);
66 
67  }else if(const SiStripRecHit2D* singleHit=dynamic_cast<const SiStripRecHit2D*>(recHit)){
68  if(!useStrip) continue;
69 
70  RawHits mono;
71 
72  mono.trajectoryMeasurement = &(*it);
73 
74  mono.angleCosine = cosine;
75  const std::vector<uint8_t> & amplitudes = singleHit->cluster()->amplitudes();
76  mono.charge = accumulate(amplitudes.begin(), amplitudes.end(), 0);
77  mono.NSaturating =0;
78  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)mono.NSaturating++;}
79 
80  mono.detId= singleHit->geographicalId();
81  hits.push_back(mono);
82 
83  }else if(const SiStripRecHit1D* single1DHit=dynamic_cast<const SiStripRecHit1D*>(recHit)){
84  if(!useStrip) continue;
85 
86  RawHits mono;
87 
88  mono.trajectoryMeasurement = &(*it);
89 
90  mono.angleCosine = cosine;
91  const std::vector<uint8_t> & amplitudes = single1DHit->cluster()->amplitudes();
92  mono.charge = accumulate(amplitudes.begin(), amplitudes.end(), 0);
93  mono.NSaturating =0;
94  for(unsigned int i=0;i<amplitudes.size();i++){if(amplitudes[i]>=254)mono.NSaturating++;}
95 
96  mono.detId= single1DHit->geographicalId();
97  hits.push_back(mono);
98 
99 
100  }else if(const SiPixelRecHit* pixelHit=dynamic_cast<const SiPixelRecHit*>(recHit)){
101  if(!usePixel) continue;
102 
103  RawHits pixel;
104 
105  pixel.trajectoryMeasurement = &(*it);
106 
107  pixel.angleCosine = cosine;
108  pixel.charge = pixelHit->cluster()->charge();
109  pixel.NSaturating=-1;
110  pixel.detId= pixelHit->geographicalId();
111  hits.push_back(pixel);
112  }
113 
114  }
115  // return hits;
116  }
i
int i
Definition: DBlmapReader.cc:9
TrajectoryStateOnSurface::localDirection
LocalVector localDirection() const
Definition: TrajectoryStateOnSurface.h:122
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:67
PV3DBase::z
T z() const
Definition: PV3DBase.h:64
SiPixelRecHit
Our base class.
Definition: SiPixelRecHit.h:23
Generated for CMSSW Reference Manual by doxygen 1.8.5