FairRoot/PandaRoot
PndGemSensor.cxx
Go to the documentation of this file.
1 // -------------------------------------------------------------------------
2 // ----- PndGemSensor source file -----
3 // ----- Created 12/02/2009 by R. Karabowicz -----
4 // -------------------------------------------------------------------------
5 
6 #include "PndGemSensor.h"
7 
8 #include "TMath.h"
9 #include "TRandom.h"
10 
11 #include <iostream>
12 #include <list>
13 #include <vector>
14 
15 using std::cout;
16 using std::endl;
17 using std::flush;
18 using std::pair;
19 using std::vector;
20 
21 
22 // ----- Default constructor -------------------------------------------
24 
25  fDetectorId = 0;
26  fType = 0;
29  fSigmaX = fSigmaY = fSigmaXY = 0.;
30 }
31 // -------------------------------------------------------------------------
32 
33 
34 
35 // ----- Enhanced constructor (by z0 and d) ------------------------------------------
36 PndGemSensor::PndGemSensor(TString tempName, Int_t detId, Int_t iType,
39  Double_t innerRad, Double_t outerRad,
40  Double_t d,
41  Double_t stripAngle0, Double_t stripAngle1,
42  Double_t pitch0, Double_t pitch1)
43 {
44 
45  fName = tempName.Data();
46  fDetectorId = detId;
47  fType = iType;
48  fPosition[0] = x0;
49  fPosition[1] = y0;
50  fPosition[2] = z0; // z position of the station
52  fInnerRadius = innerRad;
53  fOuterRadius = outerRad;
54  fD = d; // thickness of the station
55  fStripAngle[0]= stripAngle0; // strip angle
56  fStripAngle[1]= stripAngle1; // strip angle
57  fPitch[0] = pitch0; // strip pitch
58  fPitch[1] = pitch1; // strip pitch
59 
60  if ( fType == 0 ) { // r phi version
61 // fNChannelsFront = (Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[0]));
62 // fNChannelsBack = 2*(Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[1]));
63  fNChannelsFront = (Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[0] + 0.5);
64  fNChannelsBack = 2*(Int_t)((fOuterRadius-fInnerRadius)/fPitch[1]+0.5);
65  }
66  if ( fType == 1 ) { // tilted version, should not be used any more
67 // fNChannelsFront = (Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[0]));
68 // fNChannelsBack = (Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[1]));
69  fNChannelsFront = (Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[0]+0.5);
70  fNChannelsBack = (Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[1]+0.5);
71  }
72  if ( fType == 2 ) { // x y version
73 // fNChannelsFront = (Int_t)(TMath::Ceil(2.*fOuterRadius/fPitch[0]))+
74 // (Int_t)(TMath::Ceil(2.*fInnerRadius/fPitch[0]));
75 // fNChannelsBack = 2*(Int_t)(TMath::Ceil(2.*fOuterRadius/fPitch[1]));
76  fNChannelsFront = (Int_t)(2.*fOuterRadius/fPitch[0] + 2.*fInnerRadius/fPitch[0] + 0.5);
77  fNChannelsBack = 2*(Int_t)(2.*fOuterRadius/fPitch[1]+0.5);
78  }
79  if ( fType == 3 ) { // r_tree phi version
80  Int_t ifac=1;
81  if( 1 < fOuterRadius / fInnerRadius && fOuterRadius / fInnerRadius < 100){
82  //extract highest bit number to get radial tree factor
83  ifac = (Int_t) (fOuterRadius / fInnerRadius);
84  ifac |= (ifac >> 1);
85  ifac |= (ifac >> 2);
86  ifac |= (ifac >> 4);
87  ifac -= (ifac >> 1);
88 // fNChannelsFront = ifac *(Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[0]));
89  fNChannelsFront = ifac *(Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[0] + 0.5);
90  }
91  else{
92  cout << tempName.Data() << "-W- !!! type " << fType << " with OuterRadius/InnerRadius " << fOuterRadius/fInnerRadius << " is currently not supported !!!" << endl;
93  cout << " Please use type 0 instead !!!" << endl;
94 // fNChannelsFront = (Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[0]));
95  fNChannelsFront = (Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[0]+0.5);
96  }
97 // fNChannelsBack = 2*(Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[1]));
98  fNChannelsBack = 2*(Int_t)((fOuterRadius-fInnerRadius)/fPitch[1] + 0.5);
99  }
100 
101  // cout << tempName.Data() << " type " << fType << " has " << fNChannelsFront << " front and " << fNChannelsBack << " back channels" << endl;
102 
103  fSigmaX = fSigmaY = fSigmaXY = 0.;
104 }
105 // -------------------------------------------------------------------------
106 
107 // ----- Enhanced constructor (by z0 and d) ------------------------------------------
108 PndGemSensor::PndGemSensor(TString tempName, Int_t stationNr, Int_t sectorNr, Int_t iType,
111  Double_t innerRad, Double_t outerRad,
112  Double_t d,
113  Double_t stripAngle0, Double_t stripAngle1,
114  Double_t pitch0, Double_t pitch1)
115 {
116 
117  fName = tempName.Data();
118  SetDetectorId(stationNr, sectorNr);
119  fType = iType;
120  fPosition[0] = x0;
121  fPosition[1] = y0;
122  fPosition[2] = z0; // z position of the station
124  fInnerRadius = innerRad;
125  fOuterRadius = outerRad;
126  fD = d; // thickness of the station
127  fStripAngle[0]= stripAngle0; // strip angle
128  fStripAngle[1]= stripAngle1; // strip angle
129  fPitch[0] = pitch0; // strip pitch
130  fPitch[1] = pitch1; // strip pitch
131 
132  if ( fType == 0 ) { // r phi version
133 // fNChannelsFront = (Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[0]));
134 // fNChannelsBack = 2*(Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[1]));
135  fNChannelsFront = (Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[0] + 0.5);
136  fNChannelsBack = 2*(Int_t)((fOuterRadius-fInnerRadius)/fPitch[1]+0.5);
137  }
138  if ( fType == 1 ) { // tilted version, should not be used any more
139 // fNChannelsFront = (Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[0]));
140 // fNChannelsBack = (Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[1]));
141  fNChannelsFront = (Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[0]+0.5);
142  fNChannelsBack = (Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[1]+0.5);
143  }
144  if ( fType == 2 ) { // x y version
145 // fNChannelsFront = 2*(Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]))+
146 // 4*(Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
147 // fNChannelsBack = 2*(Int_t)(TMath::Ceil(2.*fOuterRadius/fPitch[1]));
148  fNChannelsFront = (Int_t)( 2*(fOuterRadius-fInnerRadius)/fPitch[0] + 4*fInnerRadius/fPitch[0] + 0.5);
149  fNChannelsBack = 2*(Int_t)(2.*fOuterRadius/fPitch[1]+0.5);
150  }
151  if ( fType == 3 ) { // r_tree phi version
152  Int_t ifac=1;
153  if( 1 < fOuterRadius / fInnerRadius && fOuterRadius / fInnerRadius < 100){
154  //extract highest bit number
155  ifac = (Int_t) (fOuterRadius / fInnerRadius);
156  ifac |= (ifac >> 1);
157  ifac |= (ifac >> 2);
158  ifac |= (ifac >> 4);
159  ifac -= (ifac >> 1);
160 // fNChannelsFront = ifac *(Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[0]));
161  fNChannelsFront = ifac *(Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[0] + 0.5);
162  }
163  else{
164  cout << tempName.Data() << "-W- !!! type " << fType << " with OuterRadius/InnerRadius " << fOuterRadius/fInnerRadius << " is currently not supported !!!" << endl;
165  cout << " Please use type 0 instead !!!" << endl;
166 // fNChannelsFront = (Int_t)(TMath::Ceil(2.*TMath::Pi()*fInnerRadius/fPitch[0]));
167  fNChannelsFront = (Int_t)(2.*TMath::Pi()*fInnerRadius/fPitch[0]+0.5);
168  }
169 // fNChannelsBack = 2*(Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[1]));
170  fNChannelsBack = 2*(Int_t)((fOuterRadius-fInnerRadius)/fPitch[1] + 0.5);
171  }
172 
173  // cout << tempName.Data() << " type " << fType << " has " << fNChannelsFront << " front and " << fNChannelsBack << " back channels" << endl;
174 
175  fSigmaX = fSigmaY = fSigmaXY = 0.;
176 }
177 // -------------------------------------------------------------------------
178 
179 // ----- Copy constructor -----------------------------------------------
181  TNamed(tempSensor),
182  fDetectorId (tempSensor.fDetectorId),
183  fType (tempSensor.fType),
184  // fPosition {tempSensor.fPosition[0],
185  // tempSensor.fPosition[1],
186  // tempSensor.fPosition[2]},
187  fRotation (tempSensor.fRotation),
188  fInnerRadius (tempSensor.fInnerRadius),
189  fOuterRadius (tempSensor.fOuterRadius),
190  fD (tempSensor.fD),
191  // fStripAngle[0] (tempSensor.fStripAngle[0]),
192  // fStripAngle[1] (tempSensor.fStripAngle[1]),
193  // fPitch[0] (tempSensor.fPitch[0]),
194  // fPitch[1] (tempSensor.fPitch[1]),
195  fNChannelsFront(tempSensor.fNChannelsFront),
196  fNChannelsBack (tempSensor.fNChannelsBack),
197  fSigmaX (tempSensor.fSigmaX),
198  fSigmaY (tempSensor.fSigmaY),
199  fSigmaXY (tempSensor.fSigmaXY)
200 {
201 }
202 // -------------------------------------------------------------------------
203 
204 
205 // ----- Destructor ----------------------------------------------------
207 // -------------------------------------------------------------------------
208 
209 
210 
211 
212 
213 // ----- Public method Reset -------------------------------------------
215 }
216 // -------------------------------------------------------------------------
217 
218 
219 
220 // ----- Public method Print -------------------------------------------
222  cout << fName.Data() << ": Sensor #";
223  cout.width(3);
224  cout << GetSensorNr() << ", Type ";
225  cout.width(1);
226  cout << fType << ", centre (";
227  cout.width(6);
228  cout << fPosition[0] << ", ";
229  cout.width(6);
230  cout << fPosition[1] << ", ";
231  cout.width(6);
232  cout << fPosition[2] << ") cm, rotation " ;
233  cout.width(6);
234  cout << fRotation*180./TMath::Pi() << " deg, radius from ";
235  cout.width(6);
236  cout << fInnerRadius << " to ";
237  cout.width(6);
238  cout << fOuterRadius << ", thickness ";
239  cout.width(6);
240  cout << fD << " cm, strip angle (";
241  cout.width(6);
242  cout << fStripAngle[0] << ", ";
243  cout.width(6);
244  cout << fStripAngle[1] << ") cm, pitch (";
245  cout.width(6);
246  cout << fPitch[0] << ", ";
247  cout.width(6);
248  cout << fPitch[1] << ") cm. ";
249  cout << endl;
250 
251  Int_t tempSId = fDetectorId;
252  Int_t bc = 0;
253  cout << "SENSOR : " << flush;
254  while ( tempSId > 0 ) {
255  bc++;
256  cout << "\b" << tempSId%2 << "\b" << flush;
257  if ( bc == 27 || bc == 21 || bc == 8 || bc == 6 || bc == 5 )
258  cout << "\b|\b" << flush;
259  tempSId = tempSId/2;
260  }
261  cout << endl;
262 
263 }
264 // -------------------------------------------------------------------------
265 
266 // ----- Public method GetChannel --------------------------------------
268 
269  if (iSide !=0 && iSide != 1) {
270  cout << "-W- PndGemSensor::GetChannel: Illegal side number "
271  << iSide << endl;
272  return -1;
273  }
274 
275  if ( !Inside(x,y) ) return -1;
276  Double_t radius = TMath::Sqrt(x*x+y*y);
277 
278  if ( fType == 1 ) {
279  cout << "do not use this type anymore" << endl;
280  return -1;
281  }
282  if ( fType == 0 ) { // angle info for iSide 0, radius info for iSide 1
283  if ( iSide == 0 ) {
284  Double_t hitPhi = TMath::ACos(y/radius);
285  if ( x < 0. )
286  hitPhi = 2.*TMath::Pi()-hitPhi;
287  hitPhi = 2.*TMath::Pi()-hitPhi;
288  return (Int_t)(TMath::Ceil((hitPhi*fInnerRadius)/fPitch[iSide]));
289  }
290  if ( iSide == 1 ) {
291  if ( x < 0. )
292  return (Int_t)(TMath::Ceil((radius-fInnerRadius)/fPitch[1]));
293  else
294  return (Int_t)(TMath::Ceil((radius-fInnerRadius)/fPitch[1])) + fNChannelsBack/2;
295  }
296  }
297  if ( fType == 2 ) { // x y strips
298  if ( iSide == 0 ) { // x information encoded
299 // Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
300 // Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
301  Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5);
302  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] + 0.5);
303  if ( x <= -fInnerRadius )
304  return (Int_t)((x+fOuterRadius)/fPitch[0]);
305  if ( x < 0. && y >= 0. )
306  return nlStrips +(Int_t)((x+fInnerRadius)/fPitch[0]);
307  if ( x < 0. && y < 0. )
308  return nlStrips+ nsStrips+(Int_t)((x+fInnerRadius)/fPitch[0]);
309  // now x can't be smaller than 0.
310  if ( x >= fInnerRadius )
311  return nlStrips+4*nsStrips+(Int_t)((x-fInnerRadius)/fPitch[0]);
312  if ( y >= 0. )
313  return nlStrips+2*nsStrips+(Int_t)((x)/fPitch[0]);
314  if ( y < 0. )
315  return nlStrips+3*nsStrips+(Int_t)((x)/fPitch[0]);
316  }
317  if ( iSide == 1 ) { // y information encoded
318  if ( x < 0. )
319  return (Int_t)((y+fOuterRadius)/fPitch[1]);
320  if ( x >= 0. )
321  return (Int_t)((y+fOuterRadius)/fPitch[1])+fNChannelsBack/2;
322  }
323  }
324  if ( fType == 3 ) { // angle info for iSide 0, radius info for iSide 1
325  if ( iSide == 0 ) {
326  Double_t hitPhi = TMath::ACos(y/radius);
327  if ( x < 0. )
328  hitPhi = 2.*TMath::Pi()-hitPhi;
329  hitPhi = 2.*TMath::Pi()-hitPhi;
330  //extract radial tree factor at the OuterRadius
331  Int_t ifac=1;
332  //extract highest bit number
333  ifac = (Int_t) (fOuterRadius / fInnerRadius);
334  ifac |= (ifac >> 1);
335  ifac |= (ifac >> 2);
336  ifac |= (ifac >> 4);
337  ifac -= (ifac >> 1);
338  //extract radial tree factor at the radius
339  Int_t ifac_r=1;
340  ifac_r = (Int_t) ( radius / fInnerRadius);
341  ifac_r |= (ifac_r >> 1);
342  ifac_r |= (ifac_r >> 2);
343  ifac_r |= (ifac_r >> 4);
344  ifac_r -= (ifac_r >> 1);
345  Int_t ich_step;
346  if ( ifac_r == 0 ) return -1;
347  else ich_step = (Int_t) (ifac / ifac_r);
348  return (Int_t)( ich_step * TMath::Ceil(hitPhi*fNChannelsFront/ich_step/2./TMath::Pi()));
349  }
350  if ( iSide == 1 ) {
351  if ( x < 0. )
352  return (Int_t)(TMath::Ceil((radius-fInnerRadius)/fPitch[1]));
353  else
354  return (Int_t)(TMath::Ceil((radius-fInnerRadius)/fPitch[1])) + fNChannelsBack/2;
355  }
356  }
357  return -1;
358 }
359 // -------------------------------------------------------------------------
360 
361 // ----- Public method GetNeighbours -----------------------------------
362 Int_t PndGemSensor::GetNeighbours(Int_t iSide, Int_t iChan, Int_t& nChan1, Int_t& nChan2, Int_t& nChan3) {
363  nChan1=-1;
364  nChan2=-1;
365  nChan3=-1;
366  if ( fType == 0 ) {
367  nChan1 = iChan-1;
368  nChan2 = iChan+1;
369  return 2;
370  }
371  if ( fType == 1 ) {
372  nChan1 = iChan-1;
373  nChan2 = iChan+1;
374  return 2;
375  }
376  if ( fType == 2 ) {
377  if ( iSide == 1 ) {
378  nChan1 = iChan-1;
379  nChan2 = iChan+1;
380  return 2;
381  }
382  if ( iSide == 0 ) {
383  nChan1 = iChan-1;
384  nChan2 = iChan+1;
385  return 2;
386  }
387  }
388  if ( fType == 3 ) {
389  if ( iSide == 0 ) {// radial tree
390  // need to decide what to return...(?)
391  nChan1 = iChan-1;
392  nChan2 = iChan+1;
393  return 2;
394  }
395  if ( iSide == 1 ) {
396  nChan1 = iChan-1;
397  nChan2 = iChan+1;
398  }
399  }
400  return 0;
401 }
402 // -------------------------------------------------------------------------
403 
404 // ----- Public method GetSensorPart --------------------------------------
405 Int_t PndGemSensor::GetSensorPart(Int_t iSide, Int_t chan) {
406  if ( fType == 0 ) {
407  return -1;
408  }
409  if ( fType == 1 ) {
410  return -1;
411  }
412  if ( fType == 2 ) {
413  if ( iSide == 1 ) {
414  return -1;
415  }
416  if ( iSide == 0 ) {
417 // Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
418  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0]+0.5);
419  //if ( TMath::Abs(chan-fNChannelsFront/2) > nsStrips*2 )
420  if ( TMath::Abs(chan-fNChannelsFront/2 +0.5 ) > nsStrips*2 )
421  return 0;
422  //else if ( TMath::Abs(chan-fNChannelsFront/2) < nsStrips )
423  else if ( TMath::Abs(chan-fNChannelsFront/2 +0.5) < nsStrips )
424  return 2;
425  else
426  return 1;
427  }
428  }
429  if ( fType == 3 ) {
430  //should it be irregular??
431  return -1;
432  }
433  return -1;
434 }
435 // -------------------------------------------------------------------------
436 
437 // ----- Public method GetStripOrientation --------------------------------------
439  if (iSide !=0 && iSide != 1) return -1.;
440  if ( !Inside(x,y) ) return -1.;
441 
442  if ( fType == 1 ) {
443  return -1.;
444  }
445  if ( fType == 0 ) { // angle info for iSide 0, radius info for iSide 1
446  if ( iSide == 0 ) {
447  Double_t radius = TMath::Sqrt(x*x+y*y);
448  Double_t hitPhi = TMath::ACos(y/radius);
449  if ( x < 0. )
450  hitPhi = 2.*TMath::Pi()-hitPhi;
451  hitPhi = 2.*TMath::Pi()-hitPhi;
452  return hitPhi;
453  }
454  if ( iSide == 1 ) {
455  Double_t radius = TMath::Sqrt(x*x+y*y);
456  Double_t hitPhi = TMath::ACos(y/radius);
457  if ( x < 0. )
458  hitPhi = 2.*TMath::Pi()-hitPhi;
459  hitPhi = 2.*TMath::Pi()-hitPhi;
460  return hitPhi+TMath::Pi()/2.;
461  }
462  }
463  if ( fType == 2 ) { // x y strips
464  if ( iSide == 0 ) { // x information encoded
465  return 0.;
466  }
467  if ( iSide == 1 ) { // y information encoded
468  return TMath::Pi()/2.;
469  }
470  }
471  if ( fType == 3 ) { // angle info for iSide 0, radius info for iSide 1
472  if ( iSide == 0 ) {
473  Double_t radius = TMath::Sqrt(x*x+y*y);
474  Double_t hitPhi = TMath::ACos(y/radius);
475  if ( x < 0. )
476  hitPhi = 2.*TMath::Pi()-hitPhi;
477  hitPhi = 2.*TMath::Pi()-hitPhi;
478  return hitPhi;
479  }
480  if ( iSide == 1 ) {
481  Double_t radius = TMath::Sqrt(x*x+y*y);
482  Double_t hitPhi = TMath::ACos(y/radius);
483  if ( x < 0. )
484  hitPhi = 2.*TMath::Pi()-hitPhi;
485  hitPhi = 2.*TMath::Pi()-hitPhi;
486  return hitPhi+TMath::Pi()/2.;
487  }
488  }
489  return -1.;
490 }
491 // -------------------------------------------------------------------------
492 
493 // ----- Public method GetDistance --------------------------------------
495  if ( chan1 == -1 || chan2 == -1 ){
496  //cout<<"-W- !!! GetDistance(): channel number is -1 !!"<<endl;
497  return -1.;
498  }
499  if ( iSide == 0 && ( chan1 - fNChannelsFront/2. - 1 ) * ( chan2 - fNChannelsFront/2. ) < 0. ) return -1.;
500  if ( iSide == 1 && ( chan1 - fNChannelsBack /2. - 1 ) * ( chan2 - fNChannelsBack /2. ) < 0. ) return -1.;
501  if ( fType == 0 ) {
502  return TMath::Abs(chan1-chan2);
503  }
504  if ( fType == 1 ) {
505  return TMath::Abs(chan1-chan2);
506  }
507  if ( fType == 2 ) {
508  if ( iSide == 0 ) {
509  //Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
510 // Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
511  //Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5); //[R.K. 01/2017] unused variable?
512  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] +0.5);
513  Int_t part1 = 666;
514  Int_t part2 = 666;
515  //if ( TMath::Abs(chan1-fNChannelsFront/2) >= nsStrips*2 ) part1 = 0;
516  if ( TMath::Abs(chan1-fNChannelsFront/2 +0.5 ) > nsStrips*2 ) part1 = 0;
517  //else if ( TMath::Abs(chan1-fNChannelsFront/2) < nsStrips ) {
518  else if ( TMath::Abs(chan1-fNChannelsFront/2 +0.5 ) < nsStrips ) {
519  part1 = -1;
520  if ( chan1 < fNChannelsFront/2 ) chan1 = chan1-nsStrips;
521  else chan1 = chan1+nsStrips;
522  }
523  else part1 = 1;
524  //if ( TMath::Abs(chan2-fNChannelsFront/2) >= nsStrips*2 ) part2 = 0;
525  if ( TMath::Abs(chan2-fNChannelsFront/2 +0.5 ) > nsStrips*2 ) part2 = 0;
526  //else if ( TMath::Abs(chan2-fNChannelsFront/2) < nsStrips ) {
527  else if ( TMath::Abs(chan2-fNChannelsFront/2 +0.5 ) < nsStrips ) {
528  part2 = -1;
529  if ( chan2 < fNChannelsFront/2 ) chan2 = chan2-nsStrips;
530  else chan2 = chan2+nsStrips;
531  }
532  else part2 = 1;
533 
534  if ( part1*part2 == -1 ) return -1;
535  return TMath::Abs(chan1-chan2);
536  }
537  if ( iSide == 1 ) {
538  return TMath::Abs(chan1-chan2);
539  }
540  }
541  if ( fType == 3 ) {
542  if ( iSide == 0 ) {
543  // radial tree
544  // it returns distance assuming the inner ring as possible minimum situation.
545  Int_t ifac=1;
546  Double_t retval;
547  //extract channel step at the inner radius as radial tree factor
548  ifac = (Int_t) (fOuterRadius / fInnerRadius);
549  ifac |= (ifac >> 1);
550  ifac |= (ifac >> 2);
551  ifac |= (ifac >> 4);
552  ifac -= (ifac >> 1);
553  retval=(Double_t) TMath::Abs( (Int_t)(chan1/ifac) - (Int_t)(chan2/ifac) );
554 
555  //cout << "-I- PndGemSensor::GetDistance("<<iSide<<","<<chan1<<","<<chan2<<") returns " << retval <<endl;
556 
557  return retval;
558  }
559  if ( iSide == 1 ) {
560  return TMath::Abs(chan1-chan2);
561  }
562  }
563  return -1; // ?ADDED (Stefano)
564 }
565 // -------------------------------------------------------------------------
566 
567 // ----- Public method GetDistance -------------------------------------
568 Int_t PndGemSensor::GetDistance(Int_t iSide, Int_t chanMin, Int_t chanMax, Int_t chanTest) {
569  if ( chanMin == -1 || chanMax == -1 || chanTest ==-1 ){
570  //cout<<"-W- !!! GetDistance(): channel number is -1 !!"<<endl;
571  return -1.;
572  }
573  if ( iSide == 0 && ( chanMin - fNChannelsFront/2 ) * ( chanTest - fNChannelsFront/2 ) < 0. ) return -1;
574  if ( iSide == 1 && ( chanMin - fNChannelsBack /2 ) * ( chanTest - fNChannelsBack /2 ) < 0. ) return -1;
575  if ( fType == 0 ) {
576  if ( chanTest > chanMin && chanTest < chanMax ) return 0;
577  return TMath::Min(TMath::Abs(chanTest-chanMax),TMath::Abs(chanTest-chanMin));
578  }
579  if ( fType == 1 ) {
580  if ( chanTest > chanMin && chanTest < chanMax ) return 0;
581  return TMath::Min(TMath::Abs(chanTest-chanMax),TMath::Abs(chanTest-chanMin));
582  }
583  if ( fType == 2 ) {
584  if ( iSide == 0 ) {
585  //Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
586 // Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
587  //Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5); //[R.K. 01/2017] unused variable?
588  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] + 0.5);
589  // THIS BROKEN VERTICAL STRIPS NEED SPECIAL ATTENTION
590  // DO NOT YET CONSIDER ALL THE PROBLEMS...
591  //if ( TMath::Abs(chanMin-fNChannelsFront/2) < nsStrips ) {
592  if ( TMath::Abs(chanMin-fNChannelsFront/2 +0.5) < nsStrips ) {
593  // cout << "chanMin " << chanMin << " is in bad position, consider " << fNChannelsFront << "/" << nlStrips << "/" << nsStrips << flush;
594  if ( chanMin < fNChannelsFront/2 ) chanMin -= nsStrips;
595  else chanMin += nsStrips;
596  // cout << " -----> " << chanMin << " ( " << chanMax << " , " << chanTest << " )" << endl;
597  }
598  //if ( TMath::Abs(chanMax-fNChannelsFront/2) < nsStrips ) {
599  if ( TMath::Abs(chanMax-fNChannelsFront/2 +0.5) < nsStrips ) {
600  // cout << "chanMax " << chanMax << " is in bad position, consider " << fNChannelsFront << "/" << nlStrips << "/" << nsStrips << flush;
601  if ( chanMax < fNChannelsFront/2 ) chanMax -= nsStrips;
602  else chanMax += nsStrips;
603  // cout << " -----> " << chanMax << " ( " << chanMin << " , " << chanTest << " )" << endl;
604  }
605  //if ( TMath::Abs(chanTest-fNChannelsFront/2) < nsStrips ) {
606  if ( TMath::Abs(chanTest-fNChannelsFront/2 +0.5) < nsStrips ) {
607  // cout << "chanTest " << chanTest << " is in bad position, consider " << fNChannelsFront << "/" << nlStrips << "/" << nsStrips << flush;
608  if ( chanTest < fNChannelsFront/2 ) chanTest -= nsStrips;
609  else chanTest += nsStrips;
610  // cout << " -----> " << chanTest << " ( " << chanMin << " , " << chanMax << " )" << endl;
611  }
612  if ( chanTest > chanMin && chanTest < chanMax ) return 0;
613  return TMath::Min(TMath::Abs(chanTest-chanMax),TMath::Abs(chanTest-chanMin));
614  }
615  if ( iSide == 1 ) {
616  if ( chanTest > chanMin && chanTest < chanMax ) return 0;
617  return TMath::Min(TMath::Abs(chanTest-chanMax),TMath::Abs(chanTest-chanMin));
618  }
619  }
620  if ( fType == 3 ) {
621  if ( iSide == 0 ) {
622  // radial tree
623  if ( chanTest > chanMin && chanTest < chanMax ) return 0;
624  // it returns distance assuming the inner ring as possible minimum situation.
625  Int_t ifac=1;
626  Double_t disMin, disMax, retval;
627  //extract channel step at the inner radius as radial tree factor
628  ifac = (Int_t) (fOuterRadius / fInnerRadius);
629  ifac |= (ifac >> 1);
630  ifac |= (ifac >> 2);
631  ifac |= (ifac >> 4);
632  ifac -= (ifac >> 1);
633  disMin=(Double_t) TMath::Abs( (Int_t)(chanTest/ifac) - (Int_t)(chanMin/ifac) );
634  disMax=(Double_t) TMath::Abs( (Int_t)(chanTest/ifac) - (Int_t)(chanMax/ifac) );
635  retval = TMath::Min(disMin,disMax);
636  //cout << "-I- PndGemSensor::GetDistance("<<iSide<<","<<chanMin<<","<<chanMax<<","<<chanTest<<") returns " << retval <<endl;
637 
638  return retval;
639  }
640  if ( iSide == 1 ) {
641  if ( chanTest > chanMin && chanTest < chanMax ) return 0;
642  return TMath::Min(TMath::Abs(chanTest-chanMax),TMath::Abs(chanTest-chanMin));
643  }
644  }
645  return -1; // ? ADDED (Stefano)
646 }
647 // -------------------------------------------------------------------------
648 
649 
650 // ----- Public method -------------------------------------------------
651 // Similar to GetChannel but for comparison of Digi and Digi.
652 // Type 3 considers r-info of r-strips in this method.
654  if ( chan1 == -1 || chan2 == -1 ) return -1.;
655  if ( iSide == 0 && ( chan1 - fNChannelsFront/2. ) * ( chan2 - fNChannelsFront/2. ) < 0. ) return -1.;
656  if ( iSide == 1 && ( chan1 - fNChannelsBack /2. ) * ( chan2 - fNChannelsBack /2. ) < 0. ) return -1.;
657  if ( fType == 0 ) {
658  return TMath::Abs(chan1-chan2);
659  }
660  if ( fType == 1 ) {
661  return TMath::Abs(chan1-chan2);
662  }
663  if ( fType == 2 ) {
664  if ( iSide == 0 ) {
665 // Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
666 // Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
667  //Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5); //[R.K. 01/2017] unused variable?
668  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] + 0.5);
669  Int_t part1 = 666;
670  Int_t part2 = 666;
671  //if ( TMath::Abs(chan1-fNChannelsFront/2) > nsStrips*2 ) part1 = 0;
672  if ( TMath::Abs(chan1-fNChannelsFront/2 +0.5 ) > nsStrips*2 ) part1 = 0;
673  //else if ( TMath::Abs(chan1-fNChannelsFront/2) < nsStrips ) {
674  else if ( TMath::Abs(chan1-fNChannelsFront/2 +0.5 ) < nsStrips ) {
675  part1 = -1;
676  if ( chan1 < fNChannelsFront/2 ) chan1 = chan1-nsStrips;
677  else chan1 = chan1+nsStrips;
678  }
679  else part1 = 1;
680  //if ( TMath::Abs(chan2-fNChannelsFront/2) > nsStrips*2 ) part2 = 0;
681  if ( TMath::Abs(chan2-fNChannelsFront/2 +0.5 ) > nsStrips*2 ) part2 = 0;
682  //else if ( TMath::Abs(chan2-fNChannelsFront/2) < nsStrips ) {
683  else if ( TMath::Abs(chan2-fNChannelsFront/2 +0.5) < nsStrips ) {
684  part2 = -1;
685  if ( chan2 < fNChannelsFront/2 ) chan2 = chan2-nsStrips;
686  else chan2 = chan2+nsStrips;
687  }
688  else part2 = 1;
689 
690  if ( part1*part2 == -1 ) return -1;
691  return TMath::Abs(chan1-chan2);
692  }
693  if ( iSide == 1 ) {
694  return TMath::Abs(chan1-chan2);
695  }
696  }
697  // radial tree
698  if ( fType == 3 ) {
699  if ( iSide == 0 ) {
700  Int_t ifac=1,ifac_test=1, ifac1=1, ifac2=1,icom_ifac=0;
701  Int_t ideno,iring_test=0,iring1=0,iring2=0; //,icom_ring=0; //[R.K. 01/2017] unused variable?
702  Int_t itest1,itest2,i1ok=0,i2ok=0;
703  Int_t icom_chstep=0;
704  Double_t retval;
705  //extract the radial tree factor of outer radius
706  ifac = (Int_t) (fOuterRadius / fInnerRadius);
707  ifac |= (ifac >> 1);
708  ifac |= (ifac >> 2);
709  ifac |= (ifac >> 4);
710  ifac -= (ifac >> 1);
711  //find minimum ring for ch1 and ch2
712  itest1=(Int_t)chan1;
713  itest2=(Int_t)chan2;
714  for( ideno=ifac; ideno>=1; ideno/=2){
715  iring_test++;
716  //cout << "ideno="<< ideno <<" iring_test="<<iring_test<<" ifac_test="<<ifac_test<<endl;
717  itest1%=ideno;
718  if( itest1 == 0 && i1ok == 0){
719  //cout <<"Minimum ring of ch1 found at #"<<iring_test<<endl;
720  iring1=iring_test;
721  ifac1=ifac_test;
722  i1ok=1;
723  }
724  itest2%=ideno;
725  if( itest2 == 0 && i2ok == 0){
726  //cout <<"Minimum ring of ch2 found at #"<<iring_test<<endl;
727  iring2=iring_test;
728  ifac2=ifac_test;
729  i2ok=1;
730  }
731  ifac_test*=2;
732  }
733  //take minimum common ring
734  if( iring1 >= iring2 ){
735  //icom_ring=iring1; //[R.K. 01/2017] unused variable?
736  icom_ifac=ifac1;
737  }
738  else{
739  //icom_ring=iring2; //[R.K. 01/2017] unused variable?
740  icom_ifac=ifac2;
741  }
742  //channel step at the common ring
743  if(ifac%icom_ifac == 0){
744  icom_chstep=ifac/icom_ifac;
745  }
746  else{
747  cout<<"-W- !!! something wrong on ifac..in PndGemSensor::GetDistance2() !!!"<<endl;
748  }
749  retval = TMath::Abs(chan1-chan2)/(Double_t)icom_chstep;
750 
751  //cout << "-I- PndGemSensor::GetDistance2("<<iSide<<","<<chan1<<","<<chan2<<") returns "<< retval <<" for Common ring #"<<icom_ring<<" ifac="<<icom_ifac<<" ch_step="<<icom_chstep<<endl;
752 
753  return retval;
754  }
755  if ( iSide == 1 ) {
756  return TMath::Abs(chan1-chan2);
757  }
758  }
759  return -1; // ? ADDED (Stefano)
760 }
761 // -------------------------------------------------------------------------
762 
763 // ----- Public method GetMeanChannel --------------------------------------
764 Double_t PndGemSensor::GetMeanChannel(Int_t iSide, Double_t chan1, Double_t weight1, Double_t chan2, Double_t weight2) {
765  if ( chan1 == -1 || chan2 == -1 ) return -1.;
766  if ( iSide == 0 && ( chan1 - fNChannelsFront/2. ) * ( chan2 - fNChannelsFront/2. ) < 0. ) return -1.;
767  if ( iSide == 1 && ( chan1 - fNChannelsBack /2. ) * ( chan2 - fNChannelsBack /2. ) < 0. ) return -1.;
768  if ( fType == 0 ) {
769  return (chan1*weight1+chan2*weight2)/(weight1+weight2);
770  }
771  if ( fType == 1 ) {
772  return (chan1*weight1+chan2*weight2)/(weight1+weight2);
773  }
774  if ( fType == 2 ) {
775  if ( iSide == 0 ) {
776  //Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
777 // Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
778  Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5);
779  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] + 0.5);
780  Int_t part1 = 666;
781  Int_t part2 = 666;
782  //if ( TMath::Abs(chan1-fNChannelsFront/2) >= nsStrips*2 ) part1 = 0;
783  if ( TMath::Abs(chan1-fNChannelsFront/2 +0.5 ) > nsStrips*2 ) part1 = 0;
784  //else if ( TMath::Abs(chan1-fNChannelsFront/2) < nsStrips ) {
785  else if ( TMath::Abs(chan1-fNChannelsFront/2 +0.5 ) < nsStrips ) {
786  part1 = -1;
787  if ( chan1 < fNChannelsFront/2 ) chan1 = chan1-nsStrips;
788  else chan1 = chan1+nsStrips;
789  }
790  else part1 = 1;
791  //if ( TMath::Abs(chan2-fNChannelsFront/2) >= nsStrips*2 ) part2 = 0;
792  if ( TMath::Abs(chan2-fNChannelsFront/2 +0.5) > nsStrips*2 ) part2 = 0;
793  //else if ( TMath::Abs(chan2-fNChannelsFront/2) < nsStrips ) {
794  else if ( TMath::Abs(chan2-fNChannelsFront/2 +0.5) < nsStrips ) {
795  part2 = -1;
796  if ( chan2 < fNChannelsFront/2 ) chan2 = chan2-nsStrips;
797  else chan2 = chan2+nsStrips;
798  }
799  else part2 = 1;
800 
801  if ( part1*part2 == -1 ) return -1;
802 
803  Double_t meanCh = (chan1*weight1+chan2*weight2)/(weight1+weight2);
804 
805  if ( part1+part2 < 0 ) {
806  if ( meanCh > nlStrips && meanCh < nlStrips+nsStrips ) meanCh += nsStrips;
807  if ( meanCh > fNChannelsFront/2+nsStrips && meanCh < fNChannelsFront-nlStrips ) meanCh -= nsStrips;
808  }
809 
810  return meanCh;
811  }
812  if ( iSide == 1 ) {
813  return (chan1*weight1+chan2*weight2)/(weight1+weight2);
814  }
815  }
816  if ( fType == 3 ) {
817  return (chan1*weight1+chan2*weight2)/(weight1+weight2);
818  }
819  return -1; // ? ADDED (Stefano)
820 }
821 // -------------------------------------------------------------------------
822 
823 // ----- Public method GetChannel2 --------------------------------------
824 Int_t PndGemSensor::GetChannel2(Double_t x, Double_t y, Int_t iSide, Double_t& feeDist) {
825  feeDist = -1.;
826 
827  if (iSide !=0 && iSide != 1) {
828  cout << "-W- PndGemSensor::GetChannel: Illegal side number "
829  << iSide << endl;
830  return -1;
831  }
832 
833  if ( !Inside(x,y) ) return -1;
834  Double_t radius = TMath::Sqrt(x*x+y*y);
835 
836  if ( fType == 1 ) {
837  cout << "-E- !!! do not use this type anymore" << endl;
838  return -1;
839  }
840  if ( fType == 0 ) { // angle info for iSide 0, radius info for iSide 1
841  if ( iSide == 0 ) {
842  feeDist = fOuterRadius - radius;
843 
844  Double_t hitPhi = TMath::ACos(y/radius);
845  if ( x < 0. )
846  hitPhi = 2.*TMath::Pi()-hitPhi;
847  hitPhi = 2.*TMath::Pi()-hitPhi;
848  return (Int_t)(TMath::Ceil((hitPhi*fInnerRadius)/fPitch[iSide]));
849  }
850  if ( iSide == 1 ) {
851  Double_t hitPhi = TMath::ACos(y/radius);
852  if ( x < 0. )
853  hitPhi = 2.*TMath::Pi()-hitPhi;
854  hitPhi = 2.*TMath::Pi()-hitPhi;
855  // feeDist = hitPhi;
856  if ( hitPhi < TMath::Pi()/2. )
857  feeDist = hitPhi*radius;
858  else if ( hitPhi < TMath::Pi() )
859  feeDist = (TMath::Pi()-hitPhi)*radius;
860  else if ( hitPhi < 3.*TMath::Pi()/2. )
861  feeDist = (hitPhi-TMath::Pi())*radius;
862  else
863  feeDist = (2.*TMath::Pi()-hitPhi)*radius;
864  feeDist -= fStripAngle[1]/2.; // that's the width of middle bar
865 
866  if ( x < 0. )
867  return (Int_t)(TMath::Ceil((radius-fInnerRadius)/fPitch[1]));
868  else
869  return (Int_t)(TMath::Ceil((radius-fInnerRadius)/fPitch[1])) + fNChannelsBack/2;
870  }
871  }
872  if ( fType == 2 ) { // x y strips
873  if ( iSide == 0 ) { // x information encoded
874  Double_t hitPhi = TMath::ACos(y/radius);
875  if ( x < 0. )
876  hitPhi = 2.*TMath::Pi()-hitPhi;
877  hitPhi = 2.*TMath::Pi()-hitPhi;
878  Double_t feeY = TMath::Cos(hitPhi)*fOuterRadius;
879  feeDist = TMath::Abs(feeY)-TMath::Abs(y);
880 
881  //Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
882  //Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
883  Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5);
884  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] + 0.5);
885  if ( x <= -fInnerRadius )
886  return (Int_t)((x+fOuterRadius)/fPitch[0]);
887  if ( x < 0. && y >= 0. )
888  return (Int_t)(nlStrips +(x+fInnerRadius)/fPitch[0]);
889  if ( x < 0. && y < 0. )
890  return (Int_t)(nlStrips+ nsStrips+(x+fInnerRadius)/fPitch[0]);
891  // now x can't be smaller than 0.
892  if ( x >= fInnerRadius )
893  return (Int_t)(nlStrips+4*nsStrips+(x-fInnerRadius)/fPitch[0]);
894  if ( y >= 0. )
895  return (Int_t)(nlStrips+2*nsStrips+(x)/fPitch[0]);
896  if ( y < 0. )
897  return (Int_t)(nlStrips+3*nsStrips+(x)/fPitch[0]);
898  }
899  if ( iSide == 1 ) { // y information encoded
900  feeDist = TMath::Abs(x) - fStripAngle[1]/2.;
901 
902  if ( x < 0. )
903  return (Int_t)((y+fOuterRadius)/fPitch[1]);
904  if ( x >= 0. )
905  return (Int_t)((y+fOuterRadius)/fPitch[1]+fNChannelsBack/2);
906  }
907  }
908  if ( fType == 3 ) { // angle info for iSide 0, radius info for iSide 1
909  if ( iSide == 0 ) {
910  feeDist = fOuterRadius - radius;
911 
912  Double_t hitPhi = TMath::ACos(y/radius);
913  if ( x < 0. )
914  hitPhi = 2.*TMath::Pi()-hitPhi;
915  hitPhi = 2.*TMath::Pi()-hitPhi;
916  //extract channel factor at the Outer radius
917  Int_t ifac=1;
918  ifac = (Int_t) (fOuterRadius / fInnerRadius);
919  ifac |= (ifac >> 1);
920  ifac |= (ifac >> 2);
921  ifac |= (ifac >> 4);
922  ifac -= (ifac >> 1);
923  //extract channel factor at the hit radius
924  Int_t ifac_r=1;
925  ifac_r = (Int_t) ( radius / fInnerRadius);
926  ifac_r |= (ifac_r >> 1);
927  ifac_r |= (ifac_r >> 2);
928  ifac_r |= (ifac_r >> 4);
929  ifac_r -= (ifac_r >> 1);
930  //calculate channel step at the radius
931  Int_t ich_step;
932  if ( ifac_r == 0 ) return -1;
933  else ich_step = (Int_t) (ifac / ifac_r);
934  //return channel number exists at the radius (use ceil as type0)
935  return (Int_t)( ich_step * TMath::Ceil(hitPhi*fNChannelsFront/ich_step/2./TMath::Pi()));
936  }
937  if ( iSide == 1 ) {
938  Double_t hitPhi = TMath::ACos(y/radius);
939  if ( x < 0. )
940  hitPhi = 2.*TMath::Pi()-hitPhi;
941  hitPhi = 2.*TMath::Pi()-hitPhi;
942  if ( hitPhi < TMath::Pi()/2. )
943  feeDist = hitPhi*radius;
944  else if ( hitPhi < TMath::Pi() )
945  feeDist = (TMath::Pi()-hitPhi)*radius;
946  else if ( hitPhi < 3.*TMath::Pi()/2. )
947  feeDist = (hitPhi-TMath::Pi())*radius;
948  else
949  feeDist = (2.*TMath::Pi()-hitPhi)*radius;
950  feeDist -= fStripAngle[1]/2.; // that's the width of middle bar
951  if ( x < 0. )
952  return (Int_t)(TMath::Ceil((radius-fInnerRadius)/fPitch[1]));
953  else
954  return (Int_t)(TMath::Ceil((radius-fInnerRadius)/fPitch[1])) + fNChannelsBack/2;
955  }
956  }
957  return -1.;
958 }
959 
960 // -------------------------------------------------------------------------
961 // COMMENTED to avoid confusion with Int_t GetChannel()
962 //
963 // -------------------------------------------------------------------------
964 // ----- Public method GetChannel --------------------------------------
965 // Double_t PndGemSensor::GetChannel(Double_t x, Double_t y, Int_t iSide, Double_t& stripWidth) {
966 // stripWidth = fPitch[iSide];
967 
968 // if (iSide !=0 && iSide != 1) {
969 // cout << "-W- PndGemSensor::GetChannel: Illegal side number "
970 // << iSide << endl;
971 // return -1;
972 // }
973 
974 // if ( !Inside(x,y) ) return -1;
975 // Double_t radius = TMath::Sqrt(x*x+y*y);
976 
977 // if ( fType == 1 ) {
978 // cout << "do not use this type anymore" << endl;
979 // return -1;
980 // }
981 // if ( fType == 0 ) { // angle info for iSide 0, radius info for iSide 1
982 // if ( iSide == 0 ) {
983 // stripWidth = stripWidth*radius/fInnerRadius;
984 // Double_t hitPhi = TMath::ACos(y/radius);
985 // if ( x < 0. )
986 // hitPhi = 2.*TMath::Pi()-hitPhi;
987 // hitPhi = 2.*TMath::Pi()-hitPhi;
988 // return (hitPhi*fInnerRadius)/fPitch[iSide];
989 // }
990 // if ( iSide == 1 ) {
991 // if ( x < 0. )
992 // return (radius-fInnerRadius)/fPitch[1];
993 // else
994 // return (radius-fInnerRadius)/fPitch[1] + fNChannelsBack/2;
995 // }
996 // }
997 // if ( fType == 2 ) { // x y strips
998 // if ( iSide == 0 ) { // x information encoded
999 // Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
1000 // Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
1001 // if ( x <= -fInnerRadius )
1002 // return (x+fOuterRadius)/fPitch[0];
1003 // if ( x < 0. && y >= 0. )
1004 // return nlStrips +(x+fInnerRadius)/fPitch[0];
1005 // if ( x < 0. && y < 0. )
1006 // return nlStrips+ nsStrips+(x+fInnerRadius)/fPitch[0];
1007 // // now x can't be smaller than 0.
1008 // if ( x >= fInnerRadius )
1009 // return nlStrips+4*nsStrips+(x-fInnerRadius)/fPitch[0];
1010 // if ( y >= 0. )
1011 // return nlStrips+2*nsStrips+(x)/fPitch[0];
1012 // if ( y < 0. )
1013 // return nlStrips+3*nsStrips+(x)/fPitch[0];
1014 // }
1015 // if ( iSide == 1 ) { // y information encoded
1016 // if ( x < 0. )
1017 // return (y+fOuterRadius)/fPitch[1];
1018 // if ( x >= 0. )
1019 // return (y+fOuterRadius)/fPitch[1]+fNChannelsBack/2;
1020 // }
1021 // }
1022 // return -1.;
1023 // }
1024 // -------------------------------------------------------------------------
1025 
1026 // ----- Public method Inside ------------------------------------------
1028  if ( TMath::Abs(x) < fStripAngle[1]/2. ) return kFALSE;
1029  Double_t radSq = x*x+y*y;
1030  if ( radSq < fInnerRadius*fInnerRadius ) return kFALSE;
1031  if ( radSq > fOuterRadius*fOuterRadius ) return kFALSE;
1032  return kTRUE;
1033 }
1034 // -------------------------------------------------------------------------
1035 
1036 // ----- Public method Inside ------------------------------------------
1038  if ( radius < fInnerRadius ) return kFALSE;
1039  if ( radius > fOuterRadius ) return kFALSE;
1040  return kTRUE;
1041 }
1042 // -------------------------------------------------------------------------
1043 
1044 // ----- Public method Intersect ---------------------------------------
1045 Int_t PndGemSensor::Intersect(Double_t iFStrip, Double_t iBStrip, Double_t& xCross, Double_t& yCross, Double_t& zCross) {
1046  // cout << "trying to find intersection of strip " << iFStrip << " and " << iBStrip << endl;
1047 
1048  if ( fType == -1 ) {
1049  cout << "-E- !!! not supported anymore" << endl;
1050  return -1;
1051  }
1052  // the hits are on different sides
1053  if ( iFStrip < fNChannelsFront/2 && iBStrip >= fNChannelsBack/2 ) return -1;
1054  if ( iFStrip >= fNChannelsFront/2 && iBStrip < fNChannelsBack/2 ) return -1;
1055 
1056  Double_t bs = iBStrip;
1057  if ( bs >= fNChannelsBack/2 ) bs -= fNChannelsBack/2;
1058  if ( fType == 0 ) { // r phi strips
1059  Double_t phi = fPitch[0]*((Double_t)iFStrip-0.5) / fInnerRadius;
1060  Double_t radius = fPitch[1]*((Double_t) bs -0.5) + fInnerRadius;
1061  yCross = radius*TMath::Cos(phi);
1062  xCross = -radius*TMath::Sin(phi);
1063  zCross = fPosition[2];
1064  if ( Inside(xCross,yCross) )
1065  return fDetectorId;
1066  else
1067  return -1;
1068  }
1069  if ( fType == 2 ) { // x y strips
1070  //Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
1071  //Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
1072  Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5);
1073  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] + 0.5);
1074  //Double_t x = -666.; //[R.K. 01/2017] unused variable?
1075  yCross = -fOuterRadius+(Double_t(bs)+0.5)*fPitch[1];
1076  zCross = fPosition[2];
1077  if ( iFStrip < nlStrips ) {
1078  xCross = -fOuterRadius+(Double_t(iFStrip)+0.5)*fPitch[0];
1079  if ( !Inside(xCross,yCross) ) return -1;
1080  return fDetectorId;
1081  }
1082  if ( iFStrip < nlStrips+ nsStrips ) {
1083  if ( bs < fNChannelsBack/4 ) return -1;
1084  xCross = -fOuterRadius+(Double_t(iFStrip)+0.5)*fPitch[0];
1085  if ( !Inside(xCross,yCross) ) return -1;
1086  return fDetectorId;
1087  }
1088  if ( iFStrip < nlStrips+2*nsStrips ) {
1089  if ( bs >= fNChannelsBack/4 ) return -1;
1090  xCross = -fOuterRadius+(Double_t(iFStrip-nsStrips)+0.5)*fPitch[0];
1091  if ( !Inside(xCross,yCross) ) return -1;
1092  return fDetectorId;
1093  }
1094  if ( iFStrip < nlStrips+3*nsStrips ) {
1095  if ( bs < fNChannelsBack/4 ) return -1;
1096  xCross = -fOuterRadius+(Double_t(iFStrip-nsStrips)+0.5)*fPitch[0];
1097  if ( !Inside(xCross,yCross) ) return -1;
1098  return fDetectorId;
1099  }
1100  if ( iFStrip < nlStrips+4*nsStrips ) {
1101  if ( bs >= fNChannelsBack/4 ) return -1;
1102  xCross = -fOuterRadius+(Double_t(iFStrip-3*nsStrips)+0.5)*fPitch[0];
1103  if ( !Inside(xCross,yCross) ) return -1;
1104  return fDetectorId;
1105  }
1106  }
1107  if ( fType == 3 ) {
1108  Double_t phi = 2. * TMath::Pi() * ((Double_t)iFStrip-0.5) / (Double_t)fNChannelsFront;
1109  Double_t radius = fPitch[1]*((Double_t) bs -0.5) + fInnerRadius;
1110  yCross = radius*TMath::Cos(phi);
1111  xCross = -radius*TMath::Sin(phi);
1112  zCross = fPosition[2];
1113  if ( Inside(xCross,yCross) )
1114  return fDetectorId;
1115  else
1116  return -1;
1117  }
1118  return -1;
1119 }
1120 // -------------------------------------------------------------------------
1121 
1122 
1123 // ----- Public method Intersect ---------------------------------------
1124 Int_t PndGemSensor::Intersect(Double_t iFStrip, Double_t iBStrip, Double_t& xCross, Double_t& yCross, Double_t& zCross,
1125  Double_t& dr, Double_t& dp) {
1126  // cout << "trying to find intersection of strip " << iFStrip << " and " << iBStrip << endl;
1127 
1128  if ( fType == -1 ) {
1129  cout << "-E- !!! not supported anymore" << endl;
1130  return -1;
1131  }
1132  // the hits are on different sides
1133  //cout << iFStrip << " of " << fNChannelsFront << " and " << iBStrip << " of " << fNChannelsBack << endl;
1134  if ( iFStrip < fNChannelsFront/2 && iBStrip >= fNChannelsBack/2 ) return -1;
1135  if ( iFStrip >= fNChannelsFront/2 && iBStrip < fNChannelsBack/2 ) return -1;
1136  Double_t bs = iBStrip;
1137  if ( bs >= fNChannelsBack/2 ) bs -= fNChannelsBack/2;
1138  if ( fType == 0 ) { // r phi strips
1139  Double_t phi = fPitch[0]*((Double_t)iFStrip-0.5) / fInnerRadius;
1140  Double_t radius = fPitch[1]*((Double_t) bs -0.5) + fInnerRadius;
1141  yCross = radius*TMath::Cos(phi);
1142  xCross = -radius*TMath::Sin(phi);
1143  zCross = fPosition[2];
1144 
1145  dp = fPitch[0]*radius/fInnerRadius/TMath::Sqrt(12.);
1146  dr = fPitch[1]/TMath::Sqrt(12.);
1147 
1148  if ( Inside(xCross,yCross) )
1149  return fDetectorId;
1150  else
1151  return -1;
1152  }
1153  if ( fType == 2 ) { // x y strips
1154  //Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
1155 // Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
1156  Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5);
1157  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] + 0.5);
1158  //Double_t x = -666.; //[R.K. 01/2017] unused variable?
1159  yCross = -fOuterRadius+(Double_t(bs)+0.5)*fPitch[1];
1160  zCross = fPosition[2];
1161 
1162  if ( iFStrip < nlStrips )
1163  xCross = -fOuterRadius+(Double_t(iFStrip)+0.5)*fPitch[0];
1164  else if ( iFStrip < nlStrips+ nsStrips ) {
1165  if ( bs < fNChannelsBack/4 ) return -1;
1166  xCross = (Double_t(iFStrip-nlStrips-1*nsStrips)+0.5)*fPitch[0];
1167  }
1168  else if ( iFStrip < nlStrips+2*nsStrips ) {
1169  if ( bs >= fNChannelsBack/4 ) return -1;
1170  xCross = (Double_t(iFStrip-nlStrips-2*nsStrips)+0.5)*fPitch[0];
1171  }
1172  else if ( iFStrip < nlStrips+3*nsStrips ) {
1173  if ( bs < fNChannelsBack/4 ) return -1;
1174  xCross = (Double_t(iFStrip-nlStrips-2*nsStrips)+0.5)*fPitch[0];
1175  }
1176  else if ( iFStrip < nlStrips+4*nsStrips ) {
1177  if ( bs >= fNChannelsBack/4 ) return -1;
1178  xCross = (Double_t(iFStrip-nlStrips-3*nsStrips)+0.5)*fPitch[0];
1179  }
1180  else
1181  xCross = (Double_t(iFStrip-nlStrips-3*nsStrips)+0.5)*fPitch[0];
1182 
1183  if ( !Inside(xCross,yCross) ) return -1;
1184 
1185  /*
1186  Double_t rMax = TMath::Sqrt( (TMath::Abs(xCross)+fPitch[0])*(TMath::Abs(xCross)+fPitch[0])+
1187  (TMath::Abs(yCross)+fPitch[1])*(TMath::Abs(yCross)+fPitch[1]) );
1188  Double_t rMin = TMath::Sqrt( (TMath::Abs(xCross)-fPitch[0])*(TMath::Abs(xCross)-fPitch[0])+
1189  (TMath::Abs(yCross)-fPitch[1])*(TMath::Abs(yCross)-fPitch[1]) );
1190 ` Double_t phi1 = TMath::ATan( (TMath::Abs(yCross)-fPitch[1])/(TMath::Abs(xCross)+fPitch[0]) );
1191  Double_t phi2 = TMath::ATan( (TMath::Abs(yCross)+fPitch[1])/(TMath::Abs(xCross)-fPitch[0]) );
1192  dp = TMath::Abs(phi1-phi2);
1193  dp = TMath::Tan(dp)*(rMax+rMin)/2./TMath::Sqrt(12.);
1194  dr = (rMax-rMin)/TMath::Sqrt(12.);
1195  // if ( dp > 350 )
1196  // dp = TMath::Abs(360.-dp);
1197  */
1198  dr = fPitch[0]/TMath::Sqrt(12.);
1199  dp = fPitch[1]/TMath::Sqrt(12.);
1200 
1201  return fDetectorId;
1202  }
1203  if ( fType == 3 ) { // r_a phi strips
1204  Double_t phi = 2. * TMath::Pi() * ((Double_t)iFStrip-0.5) / (Double_t)fNChannelsFront;
1205  Double_t radius = fPitch[1]*((Double_t) bs -0.5) + fInnerRadius;
1206  yCross = radius*TMath::Cos(phi);
1207  xCross = -radius*TMath::Sin(phi);
1208  zCross = fPosition[2];
1209  //extract channel factor at the hit radius
1210  Int_t ifac_r=1;
1211  ifac_r = (Int_t) ( radius / fInnerRadius);
1212  ifac_r |= (ifac_r >> 1);
1213  ifac_r |= (ifac_r >> 2);
1214  ifac_r |= (ifac_r >> 4);
1215  ifac_r -= (ifac_r >> 1);
1216  dp = fPitch[0]*radius/fInnerRadius/ifac_r/TMath::Sqrt(12.);
1217  dr = fPitch[1]/TMath::Sqrt(12.);
1218 
1219  if ( Inside(xCross,yCross) )
1220  return fDetectorId;
1221  else
1222  return -1;
1223  }
1224 
1225  return -1;
1226 }
1227 // -------------------------------------------------------------------------
1228 
1229 // ----- Public method Intersect ---------------------------------------
1230 Int_t PndGemSensor::Intersect(Double_t iFStrip, Double_t iBStrip, Double_t& xCross, Double_t& yCross, Double_t& zCross,
1231  Double_t& dx, Double_t& dy, Double_t& dr, Double_t& dp) {
1232  // cout << "trying to find intersection of strip " << iFStrip << " and " << iBStrip << endl;
1233  if ( fType == -1 ) {
1234  cout << "-E- !!! not supported anymore" << endl;
1235  return -1;
1236  }
1237  // the hits are on different sides
1238  //cout << iFStrip << " of " << fNChannelsFront << " and " << iBStrip << " of " << fNChannelsBack << endl;
1239  if ( iFStrip < fNChannelsFront/2 && iBStrip >= fNChannelsBack/2 ) return -1;
1240  if ( iFStrip >= fNChannelsFront/2 && iBStrip < fNChannelsBack/2 ) return -1;
1241  Double_t bs = iBStrip;
1242  if ( bs >= fNChannelsBack/2 ) bs -= fNChannelsBack/2;
1243  if ( fType == 0 ) { // r phi strips
1244  Double_t phi = fPitch[0]*((Double_t)iFStrip-0.5) / fInnerRadius; // the angle is counted from Y axis
1245  Double_t radius = fPitch[1]*((Double_t) bs -0.5) + fInnerRadius;
1246  yCross = radius*TMath::Cos(phi);
1247  xCross = -radius*TMath::Sin(phi);
1248  zCross = fPosition[2];
1249 
1250  dp = fPitch[0]*radius/fInnerRadius/TMath::Sqrt(12.);
1251  dr = fPitch[1]/TMath::Sqrt(12.);
1252  dx = dr*TMath::Abs(TMath::Sin(phi))+dp*TMath::Abs(TMath::Cos(phi));
1253  dy = dr*TMath::Abs(TMath::Cos(phi))+dp*TMath::Abs(TMath::Sin(phi));
1254 
1255  if ( Inside(xCross,yCross) )
1256  return fDetectorId;
1257  else
1258  return -1;
1259  }
1260  if ( fType == 2 ) { // x y strips
1261  //Int_t nlStrips = (Int_t)(TMath::Ceil((fOuterRadius-fInnerRadius)/fPitch[0]));
1262  //Int_t nsStrips = (Int_t)(TMath::Ceil(fInnerRadius/fPitch[0]));
1263  Int_t nlStrips = (Int_t)((fOuterRadius-fInnerRadius)/fPitch[0] + 0.5);
1264  Int_t nsStrips = (Int_t)(fInnerRadius/fPitch[0] + 0.5);
1265  //cout << " nlStrips="<<nlStrips<<" nsStrips="<<nsStrips<<endl;
1266  //Double_t x = -666.; //[R.K. 01/2017] unused variable?
1267  yCross = -fOuterRadius+(Double_t(bs)+0.5)*fPitch[1];
1268  zCross = fPosition[2];
1269 
1270  if ( iFStrip < nlStrips )
1271  xCross = -fOuterRadius+(Double_t(iFStrip)+0.5)*fPitch[0];
1272  else if ( iFStrip < nlStrips+ nsStrips ) {
1273  if ( bs < fNChannelsBack/4 ) return -1;
1274  xCross = (Double_t(iFStrip-nlStrips-1*nsStrips)+0.5)*fPitch[0];
1275  }
1276  else if ( iFStrip < nlStrips+2*nsStrips ) {
1277  if ( bs >= fNChannelsBack/4 ) return -1;
1278  xCross = (Double_t(iFStrip-nlStrips-2*nsStrips)+0.5)*fPitch[0];
1279  }
1280  else if ( iFStrip < nlStrips+3*nsStrips ) {
1281  if ( bs < fNChannelsBack/4 ) return -1;
1282  xCross = (Double_t(iFStrip-nlStrips-2*nsStrips)+0.5)*fPitch[0];
1283  }
1284  else if ( iFStrip < nlStrips+4*nsStrips ) {
1285  if ( bs >= fNChannelsBack/4 ) return -1;
1286  xCross = (Double_t(iFStrip-nlStrips-3*nsStrips)+0.5)*fPitch[0];
1287  }
1288  else
1289  xCross = (Double_t(iFStrip-nlStrips-3*nsStrips)+0.5)*fPitch[0];
1290 
1291 
1292  if ( !Inside(xCross,yCross) ) return -1;
1293 
1294  /*
1295  Double_t rMax = TMath::Sqrt( (TMath::Abs(xCross)+fPitch[0])*(TMath::Abs(xCross)+fPitch[0])+
1296  (TMath::Abs(yCross)+fPitch[1])*(TMath::Abs(yCross)+fPitch[1]) );
1297  Double_t rMin = TMath::Sqrt( (TMath::Abs(xCross)-fPitch[0])*(TMath::Abs(xCross)-fPitch[0])+
1298  (TMath::Abs(yCross)-fPitch[1])*(TMath::Abs(yCross)-fPitch[1]) );
1299 ` Double_t phi1 = TMath::ATan( (TMath::Abs(yCross)-fPitch[1])/(TMath::Abs(xCross)+fPitch[0]) );
1300  Double_t phi2 = TMath::ATan( (TMath::Abs(yCross)+fPitch[1])/(TMath::Abs(xCross)-fPitch[0]) );
1301  dp = TMath::Abs(phi1-phi2);
1302  dp = TMath::Tan(dp)*(rMax+rMin)/2./TMath::Sqrt(12.);
1303  dr = (rMax-rMin)/TMath::Sqrt(12.);
1304  // if ( dp > 350 )
1305  // dp = TMath::Abs(360.-dp);
1306  */
1307  dx = fPitch[0]/TMath::Sqrt(12.);
1308  dy = fPitch[1]/TMath::Sqrt(12.);
1309  dr = 0.;
1310  dp = 0.;
1311 
1312  return fDetectorId;
1313  }
1314  if ( fType == 3 ) { // r_a phi strips
1315  Double_t radius = fPitch[1]*((Double_t) bs -0.5) + fInnerRadius;
1316  //extract channel factor at the Outer radius
1317  Int_t ifac=1;
1318  ifac = (Int_t) (fOuterRadius / fInnerRadius);
1319  ifac |= (ifac >> 1);
1320  ifac |= (ifac >> 2);
1321  ifac |= (ifac >> 4);
1322  ifac -= (ifac >> 1);
1323  //extract channel factor at the hit radius
1324  Int_t ifac_r=1;
1325  ifac_r = (Int_t) ( radius / fInnerRadius);
1326  if ( ifac_r == 0 ) {
1327  cout<<"-W- PndGemSensor::Intersect() radius (" << radius << ") in smaller than InnerRadius (" << fInnerRadius << ").. return -1! (strip " << iFStrip << " / " << iBStrip << " )" <<endl;
1328  return -1;
1329  }
1330  ifac_r |= (ifac_r >> 1);
1331  ifac_r |= (ifac_r >> 2);
1332  ifac_r |= (ifac_r >> 4);
1333  ifac_r -= (ifac_r >> 1);
1334  //calculate channel step at the radius
1335  Int_t ich_step;
1336  if ( ifac_r == 0 ) return -1;
1337  else ich_step = (Int_t) (ifac / ifac_r);
1338  Double_t phi = 2. * TMath::Pi() * ((Double_t)iFStrip-0.5*(Double_t)ich_step) / (Double_t)fNChannelsFront;
1339  yCross = radius*TMath::Cos(phi);
1340  xCross = -radius*TMath::Sin(phi);
1341  zCross = fPosition[2];
1342  dp = fPitch[0]*radius/fInnerRadius/ifac_r/TMath::Sqrt(12.);
1343  dr = fPitch[1]/TMath::Sqrt(12.);
1344  dx = dr*TMath::Abs(TMath::Sin(phi))+dp*TMath::Abs(TMath::Cos(phi));
1345  dy = dr*TMath::Abs(TMath::Cos(phi))+dp*TMath::Abs(TMath::Sin(phi));
1346 
1347  if ( Inside(xCross,yCross) )
1348  return fDetectorId;
1349  else
1350  return -1;
1351  }
1352 
1353  return -1;
1354 }
1355 // -------------------------------------------------------------------------
1356 
1357 
Double_t z0
Definition: checkhelixhit.C:62
Double_t x0
Definition: checkhelixhit.C:70
double dy
TObjArray * d
Int_t Intersect(Double_t iFStrip, Double_t iBStrip, Double_t &xCross, Double_t &yCross, Double_t &zCross)
static T Sqrt(const T &x)
Definition: PndCAMath.h:37
static T Sin(const T &x)
Definition: PndCAMath.h:42
TGeoRotation rotation
virtual ~PndGemSensor()
Double_t fOuterRadius
Definition: PndGemSensor.h:204
static T Cos(const T &x)
Definition: PndCAMath.h:43
Double_t fPosition[3]
Definition: PndGemSensor.h:201
Double_t GetStripOrientation(Double_t x, Double_t y, Int_t iSide)
Double_t fSigmaXY
Definition: PndGemSensor.h:220
static T Abs(const T &x)
Definition: PndCAMath.h:39
void SetDetectorId(Int_t stationNr, Int_t sensorNr)
Definition: PndGemSensor.h:87
Int_t GetChannel(Double_t x, Double_t y, Int_t iSide)
static T Min(const T &x, const T &y)
Definition: PndCAMath.h:35
Double_t fPitch[2]
Definition: PndGemSensor.h:207
Double_t
Double_t y0
Definition: checkhelixhit.C:71
Double_t GetDistance(Int_t iSide, Double_t chan1, Double_t chan2)
Double_t fInnerRadius
Definition: PndGemSensor.h:203
Int_t GetChannel2(Double_t x, Double_t y, Int_t iSide, Double_t &feeDist)
Double_t fD
Definition: PndGemSensor.h:205
Double_t fStripAngle[2]
Definition: PndGemSensor.h:206
double dx
Bool_t Inside(Double_t x, Double_t y)
Int_t GetSensorPart(Int_t iSide, Int_t chan)
Double_t fSigmaY
Definition: PndGemSensor.h:219
Double_t fSigmaX
Definition: PndGemSensor.h:218
Int_t fDetectorId
Definition: PndGemSensor.h:199
Double_t x
Int_t GetSensorNr() const
Definition: PndGemSensor.h:95
Int_t fNChannelsBack
Definition: PndGemSensor.h:211
ClassImp(PndAnaContFact)
Double_t fRotation
Definition: PndGemSensor.h:202
Double_t GetDistance2(Int_t iSide, Double_t chan1, Double_t chan2)
Int_t GetNeighbours(Int_t iSide, Int_t iChan, Int_t &nChan1, Int_t &nChan2, Int_t &nChan3)
Int_t fNChannelsFront
Definition: PndGemSensor.h:210
Double_t y
Double_t GetMeanChannel(Int_t iSide, Double_t chan1, Double_t weight1, Double_t chan2, Double_t weight2)
Double_t Pi