/[PAMELA software]/calo/flight/CaloNuclei/src/CaloNuclei.cpp
ViewVC logotype

Diff of /calo/flight/CaloNuclei/src/CaloNuclei.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.1.1.1 by mocchiut, Tue Apr 3 15:55:16 2007 UTC revision 1.3 by mocchiut, Wed Apr 4 11:31:27 2007 UTC
# Line 20  CaloNuclei::CaloNuclei(PamLevel2 *l2p){ Line 20  CaloNuclei::CaloNuclei(PamLevel2 *l2p){
20    PKT = 0;    PKT = 0;
21    atime = 0;    atime = 0;
22    N = 5;    N = 5;
23      R = 3;
24    //    //
25  };  };
26    
27  void CaloNuclei::Clear(){  void CaloNuclei::Clear(){
28    //    //
29      tr = 0;
30    interplane = 0;    interplane = 0;
31    preq = 0.;    preq = 0.;
32    postq = 0.;    postq = 0.;
# Line 40  void CaloNuclei::Clear(){ Line 42  void CaloNuclei::Clear(){
42    
43  void CaloNuclei::Print(){  void CaloNuclei::Print(){
44    //    //
45    printf("==============================\n");    Process();
46    printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);    //
47    printf(" interplane: %i\n",interplane);    printf("========================================================================\n");
48    printf(" ethr: %f \n",ethr);    printf(" OBT: %u PKT: %u ATIME: %u Track %i \n",OBT,PKT,atime,tr);
49    printf(" dedx1: %f \n",dedx1);    printf(" interplane [number of available dE/dx before interaction]: %i\n",interplane);
50    printf(" dedx3: %f \n",dedx3);    printf(" ethr [threshold used to determine interplane]:............ %f \n",ethr);
51    printf(" multhit: %i \n",multhit);    printf(" dedx1 [dE/dx from the first calorimeter plane]:........... %f \n",dedx1);
52    printf(" gap: %i \n",gap);    printf(" dedx3 [dE/dx (average) if the first 3 Si planes]:......... %f \n",dedx3);
53    printf(" preq: %f \n",preq);    printf(" multhit [true if interplane determined by multiple hits]:. %i \n",multhit);
54    printf(" postq: %f \n",postq);    printf(" gap [true if interplane determined by a gap]:............. %i \n",gap);
55    printf(" qpremean: %f \n",qpremean);    printf(" preq [total energy in MIP before the interaction plane]:.. %f \n",preq);
56    printf(" N: %i \n",N);    printf(" postq [total energy in MIP after the interaction plane]:.. %f \n",postq);
57    printf(" qpremeanN: %f \n",qpremeanN);    printf(" qpremean [truncated mean using 3 planes and 3 strips]:.... %f \n",qpremean);
58    printf("==============================\n");    printf(" N [no of used plane]:..................................... %i \n",N);
59      printf(" R [no strip used per plane ]:............................. %i \n",R);
60      printf(" qpremeanN [truncated mean using N planes and R strips]:... %f \n",qpremeanN);
61      printf("========================================================================\n");
62    //    //
63  };  };
64    
# Line 64  void CaloNuclei::Delete(){ Line 69  void CaloNuclei::Delete(){
69    
70    
71  void CaloNuclei::Process(){  void CaloNuclei::Process(){
72    //    Process(0);
73    };
74    
75    void CaloNuclei::Process(Int_t ntr){
76      //  
77    if ( !L2 ){    if ( !L2 ){
78      printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");      printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");
79      printf(" ERROR: CaloNuclei variables not filled \n");      printf(" ERROR: CaloNuclei variables not filled \n");
# Line 86  void CaloNuclei::Process(){ Line 95  void CaloNuclei::Process(){
95    //    //
96    if ( !newentry ) return;    if ( !newentry ) return;
97    //    //
98      tr = ntr;
99      //
100    //  printf(" Processing event at OBT %u PKT %u time %u \n",OBT,PKT,atime);    //  printf(" Processing event at OBT %u PKT %u time %u \n",OBT,PKT,atime);
101    //    //
102    Clear();    Clear();
103    //    //
104    PamTrack *track = 0;    PamTrack *track = 0;
105    track = L2->GetTrack(0);    track = L2->GetTrack(ntr);
106    //    //
107    Int_t view = 0;    Int_t view = 0;
108    Int_t plane = 0;    Int_t plane = 0;
# Line 236  void CaloNuclei::Process(){ Line 247  void CaloNuclei::Process(){
247      Int_t ind = -1;      Int_t ind = -1;
248      Int_t qsplane = -1;      Int_t qsplane = -1;
249      Int_t qsview = -1;      Int_t qsview = -1;
250        Int_t ind2 = -1;
251        Int_t qsplane2 = -1;
252        Int_t qsview2 = -1;
253      Float_t qme[200];      Float_t qme[200];
254      memset(qme,0,200*sizeof(Float_t));      memset(qme,0,200*sizeof(Float_t));
255        Float_t qme2[2112];
256        memset(qme2,0,2112*sizeof(Float_t));
257      //      //
258      while ( ii<L2->GetCaloLevel1()->istrip ){      while ( ii<L2->GetCaloLevel1()->istrip ){
259        //        //
# Line 263  void CaloNuclei::Process(){ Line 279  void CaloNuclei::Process(){
279              };              };
280              qme[ind] += mip;              qme[ind] += mip;
281            };            };
282              for ( Int_t ns = 0; ns < R ; ns++){
283                Int_t ms =  track->GetCaloTrack()->tibar[plane][view] - 1 - ns + (R - 1)/2;
284                if ( strip == ms ){
285                  if ( qsplane2 != plane || qsview2 != view ){
286                    qsplane2 = plane;
287                    qsview2 = view;
288                    ind2++;
289                    if ( ind2 > 2112 ) printf(" AAAGH!! \n");
290                    qme2[ind2] = 0.;
291                  };
292                  qme2[ind2] += mip;
293                };
294              };    
295          };                };      
296          //          //
297        };        };
# Line 273  void CaloNuclei::Process(){ Line 302  void CaloNuclei::Process(){
302      //      //
303      // here we must calculate qpremean, order vector qme, select 3 lowest measurements and caculate the mean...      // here we must calculate qpremean, order vector qme, select 3 lowest measurements and caculate the mean...
304      //          //    
305      //    for (Int_t l=0; l < interplane; l++){  //     for (Int_t l=0; l < interplane; l++){
306      //       printf(" qme[%i] = %f \n",l,qme[l]);    //       printf(" qme[%i] = %f \n",l,qme[l]);  
307      //       //      if ( qme[ind] < 6.25 ) qme[ind] = 10000.;        //     //       //      if ( qme[ind] < 6.25 ) qme[ind] = 10000.;      
308      //     };  //     };
309      Long64_t work[200];      Long64_t work[200];
310      ind = 0;      ind = 0;
311      Int_t l = 0;      Int_t l = 0;
312      Float_t qm = 0.;      Float_t qm = 0.;
313        Float_t qm2 = 0.;
314        //
315        Float_t qmt = ethr*0.8; // *0.9
316        //
317      Int_t uplim = TMath::Max(3,N);      Int_t uplim = TMath::Max(3,N);
318        //
319      while ( l < uplim && ind < interplane ){      while ( l < uplim && ind < interplane ){
320        qm = TMath::KOrdStat(interplane,qme,ind,work);        qm = TMath::KOrdStat(interplane,qme,ind,work);
321        if ( qm >= mesethr*0.9 ){        if ( qm >= qmt ){
322          if ( l < 3 ) qpremean += qm;          if ( l < 3 ) qpremean += qm;
         if ( l < N ) qpremeanN += qm;  
323          l++;          l++;
324          //      printf(" value no %i qm %f \n",l,qm);          //      printf(" value no %i qm %f qmt %f \n",l,qm,qmt);
325          };
326          ind++;
327        };
328        //
329        ind = 0;
330        l = 0;
331        while ( l < uplim && ind < interplane ){
332          qm2 = TMath::KOrdStat(interplane,qme2,ind,work);
333          if ( qm2 >= qmt ){        
334            if ( l < N ) qpremeanN += qm2;
335            l++;
336            //      printf(" value no %i qm %f qmt %f \n",l,qm,qmt);
337        };        };
338        ind++;        ind++;
339      };      };
Legend:
Removed from v.1.1.1.1  
changed lines
  Added in v.1.3
  ViewVC Help
Powered by ViewVC 1.1.23