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

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

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

revision 1.1.1.1 by mocchiut, Fri Nov 9 09:11:24 2007 UTC revision 1.19 by mocchiut, Tue May 18 04:03:29 2010 UTC
# Line 1  Line 1 
1  #include <CaloProfile.h>  #include <CaloProfile.h>
2    //
3    ClassImp(CaloLat);
4    ClassImp(CaloLong);
5    ClassImp(Calo2D);
6    
7    
8    int isempty(struct stack *s){
9      return (s->top == EMPTY) ? 1 : 0;
10    }
11    
12    void emptystack(struct stack* s){
13      s->top=EMPTY;
14      strcpy(s->data,"                                                 ");
15    }
16    
17    void push(struct stack* s,int item){
18      if(s->top == (MAX-1)){
19        printf("\n ERROR! STACK FULL (too many digits in SetLower/UpperLimit)\n");
20      } else {
21        s->data[++s->top]=(char)item;
22      };
23    }
24    
25    char pop(struct stack* s){
26      char ret=(char)EMPTY;
27      if(!isempty(s)){
28        ret= s->data[s->top--];
29      };
30      return ret;
31    }
32    
33    
34    int ipop(struct stack* s){
35      int ret=EMPTY;
36      if(s->top == EMPTY)
37        printf("\n ERROR! STACK EMPTY (too few digits in SetLower/UpperLimit)\n");
38      else {
39        ret= s->data[s->top--];
40      };
41      return ret;
42    }
43    
44    
45    void display(struct stack s){
46      int ss = s.top;
47      while(s.top != EMPTY){
48        printf("\n%d",s.data[s.top--]);
49      };
50      printf(" s.top %i \n",ss);
51    }
52    
53    int isoperator(char e){
54      if(e == '+' || e == '-' || e == '*' || e == '/' || e == '%')
55        return 1;
56      else
57        return 0;
58    }
59    
60    
61    int priority(char e){
62      int pri = 0;
63      if(e == '*' || e == '/' || e =='%')
64        pri = 2;
65      else {
66        if(e == '+' || e == '-') pri = 1;
67      };
68      return pri;
69    }
70    
71    void infix2postfix(char* infix, char * postfix, int insertspace){
72      char *i,*p;
73      struct stack X;
74      char n1;
75      emptystack(&X);
76      i = &infix[0];
77      p = &postfix[0];  
78      while(*i){
79        while(*i == ' ' || *i == '\t' ){
80          i++;
81        };    
82        TString c=i;
83        if( isdigit(*i) || isalpha(*i) || c.BeginsWith(".") ){
84          while( isdigit(*i) || isalpha(*i) ||  c.BeginsWith(".") ){
85            *p = *i;
86            p++;
87            i++;
88            c=i;
89          };
90          /*SPACE CODE*/
91          if(insertspace){
92            *p = ' ';
93            p++;
94          };
95          /*END SPACE CODE*/
96        };    
97        if( *i == '(' ){
98          push(&X,*i);
99          i++;
100        };
101        if( *i == ')'){
102          n1 = pop(&X);
103          while( n1 != '(' ){
104            *p = n1;
105            p++;
106            /*SPACE CODE*/
107            if(insertspace){
108              *p = ' ';
109              p++;
110            };
111            /*END SPACE CODE*/
112            n1 = pop(&X);
113          };
114          i++;
115        };    
116        if( isoperator(*i) ){
117          if(isempty(&X)){
118            push(&X,*i);
119          } else {
120            n1 = pop(&X);
121            while( priority(n1) >= priority(*i) ){
122              *p = n1;
123              p++;
124              /*SPACE CODE*/
125              if(insertspace){
126                *p = ' ';
127                p++;
128              };
129              /*END SPACE CODE*/
130              n1 = pop(&X);
131            };
132            push(&X,n1);
133            push(&X,*i);
134          };
135          i++;
136        };
137      };
138      while(!isempty(&X)){
139        n1 = pop(&X);
140        *p = n1;
141        p++;
142        /*SPACE CODE*/
143        if(insertspace){
144          *p = ' ';
145          p++;
146        };
147        /*END SPACE CODE*/
148      };
149      *p = '\0';
150    }
151    
152    
153    Float_t evaluate(char *postfix){
154      //
155      Float_t op1 = 0.;
156      Float_t op2 = 0.;
157      Float_t result = 0.;
158      //
159      TString e = postfix;
160      Float_t st[50];
161      memset(st,0,50*sizeof(Float_t));
162      TObjArray *ae = e.Tokenize(" ");
163      Int_t o = 0;
164      Int_t a = 0;    
165      Int_t ap = 0;    
166      //
167      while ( o < ae->GetEntries() ){
168        if (  ((TString)(((TObjString*)ae->At(o))->GetString())).IsFloat() ){
169          st[a] = (Float_t)((TString)(((TObjString*)ae->At(o))->GetString())).Atof();;
170          a++;
171        } else {      
172          ap = a-1;
173          op1 = st[ap--];
174          op2 = st[ap];
175          const char *p=((TString)(((TObjString*)ae->At(o))->GetString())).Data();
176          switch(*p){
177           case '+':
178            result = op2 + op1;
179            break;
180          case '-':
181            result = op2 - op1;
182            break;
183          case '/':
184            result = op2 / op1;
185            break;
186          case '*':
187            result = op2 * op1;
188            break;
189          case '%':
190            result = (Int_t)round(op2) % (Int_t)round(op1);
191            break;
192          default:
193            printf("\nInvalid Operator: %s \n",((TString)(((TObjString*)ae->At(o))->GetString())).Data());
194            return 0;
195          };
196          st[ap] = result;
197          a = ap+1;
198          //
199        };
200        o++;
201      };
202      return result;
203      //
204    }
205    
206    
207    
208  ////////////////////////////////////////////////////////////////////////              ////////////////////////////////////////////////////////////////////////            
209  /**  /**
# Line 60  CaloLat::CaloLat(){ Line 265  CaloLat::CaloLat(){
265  /**  /**
266   * Default constructor   * Default constructor
267   */   */
268  CaloLong::CaloLong(){  Calo2D::Calo2D(){
269    Clear();    Clear();
270  };  };
271    
# Line 76  CaloLat::CaloLat(PamLevel2 *l2p){   Line 281  CaloLat::CaloLat(PamLevel2 *l2p){  
281    PKT = 0;    PKT = 0;
282    atime = 0;    atime = 0;
283    //    //
284      suf = "";
285    debug = false;    debug = false;
286      usepl18x = false;
287    //    //
288  };  };
289    
290  CaloLong::CaloLong(PamLevel2 *l2p){    Calo2D::Calo2D(PamLevel2 *l2p){  
291    //    //
292    Clear();    Clear();
293    //    //
# Line 92  CaloLong::CaloLong(PamLevel2 *l2p){   Line 299  CaloLong::CaloLong(PamLevel2 *l2p){  
299    PKT = 0;    PKT = 0;
300    atime = 0;    atime = 0;
301    //    //
302      suf = "";
303    debug = false;    debug = false;
304      usepl18x = false;
305    //    //
306  };  };
307    
# Line 100  void CaloLat::Clear(){ Line 309  void CaloLat::Clear(){
309    //    //
310  };  };
311    
312  void CaloLong::Clear(){  void Calo2D::Clear(){
313    //    //
314  };  };
315    
# Line 116  void CaloLat::Print(){ Line 325  void CaloLat::Print(){
325    //    //
326  };  };
327    
328  void CaloLong::Print(){  void Calo2D::Print(){
329    //    //
330    Process();    Process();
331    //    //
332    printf("==================== Calorimeter Lateral Profile =======================\n");    printf("==================== Calorimeter 2D Profile =======================\n");
333    printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);    printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);
334  //  printf(" nx [number of X combination]:.. %i\n",nx);  //  printf(" nx [number of X combination]:.. %i\n",nx);
335  //  printf(" ny [number of Y combination]:.. %i\n",ny);  //  printf(" ny [number of Y combination]:.. %i\n",ny);
# Line 134  void CaloLat::Draw(){ Line 343  void CaloLat::Draw(){
343    Draw(-1,-1);    Draw(-1,-1);
344  };  };
345    
346  void CaloLong::Draw(){  void Calo2D::Draw(){
347    //    //
348    Process();    Process();
349    Draw(-1);    Draw(-1);
350  };  };
351    
   
352  void CaloLat::Draw(Int_t view,Int_t plane){  void CaloLat::Draw(Int_t view,Int_t plane){
353    //    //
354    Int_t minv = 0;    Int_t minv = 0;
# Line 171  void CaloLat::Draw(Int_t view,Int_t plan Line 379  void CaloLat::Draw(Int_t view,Int_t plan
379    //    //
380    for (Int_t v=minv; v<maxv;v++){    for (Int_t v=minv; v<maxv;v++){
381       for (Int_t p=minp; p<maxp;p++){       for (Int_t p=minp; p<maxp;p++){
382          TString hid = Form("clatv%ip%i",v,p);           TString hid = Form("clatv%ip%i%s",v,p,suf.Data());      
383          TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));          TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
384          if ( tc ){          if ( tc ){
385  //       tc->Clear();  //       tc->Clear();
# Line 179  void CaloLat::Draw(Int_t view,Int_t plan Line 387  void CaloLat::Draw(Int_t view,Int_t plan
387           tc = new TCanvas(hid,hid);           tc = new TCanvas(hid,hid);
388          };          };
389          //          //
390          TString thid = Form("hlatv%ip%i",v,p);            TString thid = Form("hlatv%ip%i%s",v,p,suf.Data());    
391          TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));          TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));
392          if ( th ) th->Delete();          if ( th ) th->Delete();
393  //       th->Clear();  //       th->Clear();
# Line 203  void CaloLat::Draw(Int_t view,Int_t plan Line 411  void CaloLat::Draw(Int_t view,Int_t plan
411    //    //
412  };  };
413    
414  void CaloLong::Draw(Int_t view){  void Calo2D::Draw(Int_t plane){
415    //    //
416    Int_t minv = 0;    Int_t minp = 0;
417    Int_t maxv = 0;    Int_t maxp = 0;
418    //    //
419    if ( view == -1 ){    if ( plane == -1 ){
420          maxv = -1;          minp = 0;
421            maxp = 23;
422    } else {    } else {
423          minv = view;          minp = plane;
424          maxv = view+1;          maxp = plane+1;
425    };    };
426    //    //
427    Process();    Process();
# Line 220  void CaloLong::Draw(Int_t view){ Line 429  void CaloLong::Draw(Int_t view){
429    gStyle->SetLabelSize(0.04);    gStyle->SetLabelSize(0.04);
430    gStyle->SetNdivisions(510,"XY");    gStyle->SetNdivisions(510,"XY");
431    //    //
432    if ( maxv != -1 ){    for (Int_t p=minp; p<maxp;p++){
433    for (Int_t v=minv; v<maxv;v++){      TString hid = Form("c2dp%i%s",p,suf.Data());        
434          TString hid = Form("clongv%i",v);            TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
435          TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));      if ( tc ){
436          if ( tc ){        //         tc->Clear();
437  //       tc->Clear();      } else {
438          } else {        tc = new TCanvas(hid,hid);
439           tc = new TCanvas(hid,hid);      };
440          };      //
441          //      TString thid = Form("h2dp%i%s",p,suf.Data());      
442          TString thid = Form("hlongv%i",v);            TH2F *th  = dynamic_cast<TH2F*>(gDirectory->FindObject(thid));
443          TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));      if ( th ) th->Delete();
444          if ( th ) th->Delete();      //   th->Clear();
445  //       th->Clear();      //   th->Reset();
446  //       th->Reset();      //  } else {
447  //      } else {      Int_t minx = smax[p] - 10;
448           th = new TH1F(thid,thid,22,-0.5,21.5);      if ( minx < 0 ) minx = 0;
449  //      };      Int_t maxx = minx + 20;
450          tc->cd();      if ( maxx > 95 ){
451          //        maxx = 95;
452          for (Int_t st=0;st<22;st++){        minx = 75;
453           th->Fill(st,eplane[v][st]);      };
454          };      Int_t miny = smay[p] - 10;
455          th->Draw();      if ( miny < 0 ) miny = 0;
456          tc->Modified();      Int_t maxy = miny + 20;
457          tc->Update();      if ( maxy > 95 ){
458    };        maxy = 95;
459    } else {        miny = 75;
460          //      };
461          TString hid = Form("clongvyvx");              th = new TH2F(thid,thid,20,(Float_t)minx-0.5,(Float_t)maxx-0.5,20,(Float_t)miny-0.5,(Float_t)maxy-0.5);
462          TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));      //    th = new TH2F(thid,thid,96,-0.5,95.5,96,-0.5,95.5);
463          if ( tc ){      //  };
464          } else {      tc->cd();
465           tc = new TCanvas(hid,hid);      //
466          };      for (Int_t stx=minx;stx<maxx+1;stx++){
467          //        for (Int_t sty=miny;sty<maxy+1;sty++){
468          TString thid = Form("hlongvyvx");                th->Fill(stx,sty,estrip[p][stx][sty]);
469          TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));        };
470          if ( th ) th->Delete();      };
471          th = new TH1F(thid,thid,44,-0.5,43.5);      gStyle->SetPalette(1);
472          tc->cd();      //    tc->SetLogz();    
473          Int_t pp=0;      //    th->Draw("colbox");
474          for (Int_t st=0;st<22;st++){      th->Draw("cont4");
475            for (Int_t v=1; v>=0;v--){      tc->Modified();
476                  //      tc->Update();
                 th->Fill(pp,eplane[v][st]);  
                 //  
                 pp++;  
           };  
         };  
         th->Draw();  
         tc->Modified();  
         tc->Update();  
477    };    };
478    //    //
479    gStyle->SetLabelSize(0);    gStyle->SetLabelSize(0);
# Line 285  void CaloLat::Delete(){ Line 486  void CaloLat::Delete(){
486    //delete this;    //delete this;
487  };  };
488    
489  void CaloLong::Delete(){  void Calo2D::Delete(){
490    Clear();    Clear();
491    //delete this;    //delete this;
492  };  };
493    
494    
495  void CaloLat::Process(){  void CaloLat::Process(){
496    //      //  
497    if ( !L2 ){    if ( !L2 ){
# Line 327  void CaloLat::Process(){ Line 529  void CaloLat::Process(){
529    //    //
530    for (Int_t i=0; i<L2->GetCaloLevel1()->istrip ; i++){    for (Int_t i=0; i<L2->GetCaloLevel1()->istrip ; i++){
531      mip1 = L2->GetCaloLevel1()->DecodeEstrip(i,view1,plane1,strip1);      mip1 = L2->GetCaloLevel1()->DecodeEstrip(i,view1,plane1,strip1);
532        //
533        if ( !usepl18x && view1==0 && plane1==18 ) mip1 = 0.;
534        //
535      estrip[view1][plane1][strip1] = mip1;      estrip[view1][plane1][strip1] = mip1;
536    };    };
537    //    //
# Line 335  void CaloLat::Process(){ Line 540  void CaloLat::Process(){
540    //    //
541  };  };
542    
543  void CaloLong::Process(){  void Calo2D::Process(){
544    //      //  
545    if ( !L2 ){    if ( !L2 ){
546      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");
# Line 364  void CaloLong::Process(){ Line 569  void CaloLong::Process(){
569    //    //
570    // let's start    // let's start
571    //    //
572      Float_t es[2][22][96];
573      memset(es,0, 4224*sizeof(Float_t));
574      memset(estrip,0, 4224*sizeof(Float_t));
575      Float_t mip1 = 0.;
576      Int_t view1 = 0;
577      Int_t plane1 = 0;
578      Int_t strip1 = 0;
579      //
580      for (Int_t i=0; i<L2->GetCaloLevel1()->istrip ; i++){
581        mip1 = L2->GetCaloLevel1()->DecodeEstrip(i,view1,plane1,strip1);
582        //
583        if ( !usepl18x && view1==0 && plane1==18 ) mip1 = 0.;
584        //
585        es[view1][plane1][strip1] = mip1;
586      };
587      //
588      Int_t plane2 = 0;
589      Float_t emax[23];
590      memset(emax,0,sizeof(Float_t)*23);
591      memset(smax,0,sizeof(Int_t)*23);
592      memset(smay,0,sizeof(Int_t)*23);
593      //
594      //
595      for (Int_t p=0; p < 23 ; p++){
596        //
597        plane1 = p-1;
598        plane2 = p;
599        //
600        if ( p == 0 ){
601          plane1 = -1;
602          plane2 = 0;
603        };
604        if ( p == 22 ){
605          plane1 = 21;
606          plane2 = -1;
607        };
608        //
609        for (Int_t s=0; s < 96 ; s++){ // x
610          for (Int_t ss=0; ss < 96 ; ss++){ // y
611            if ( p == 0 ){
612              estrip[p][s][ss] += es[1][plane2][ss];
613              if ( (es[1][plane2][ss]) > emax[p] ){
614                smax[p] = 45;
615                smay[p] = ss;
616                emax[p] = es[1][plane2][ss] ;
617              };
618            };
619            if ( p > 0 && p < 22 ){
620              estrip[p][s][ss] += es[0][plane1][s] + es[1][plane2][ss];
621              if ( (es[0][plane1][s] + es[1][plane2][ss]) > emax[p] ){
622                smax[p] = s;
623                smay[p] = ss;
624                emax[p] = es[0][plane1][s] + es[1][plane2][ss] ;
625              };
626            };      
627            if ( p == 22 ){
628              estrip[p][s][ss] += es[0][plane1][s];
629              if ( (es[1][plane2][s]) > emax[p] ){
630                smax[p] = s;
631                smay[p] = 45;
632                emax[p] = es[1][plane2][s] ;
633              };
634            };
635          };
636        };
637        //
638      };
639      //
640      if ( debug ) this->Print();
641      if ( debug ) printf(" exit \n");
642      //
643    };
644    
645    
646    /**
647     * Default constructor
648     */
649    CaloLong::CaloLong(){
650      Clear();
651    };
652    
653    CaloLong::CaloLong(PamLevel2 *l2p){  
654      //
655      Clear();
656      //
657      L2 = l2p;
658      //
659      if ( !L2->IsORB() ) printf(" WARNING: OrbitalInfo Tree is needed, the plugin could not work properly without it \n");
660      //
661      OBT = 0;
662      PKT = 0;
663      atime = 0;
664      //
665      clp = L2->GetCaloLevel2();
666      //
667      sel = true;
668      cont = false;
669      N = 0;
670      NC = 22;
671      mask18b = -1;
672      //
673      no18x = true;
674      usepl18x = !no18x;
675      debug = false;
676      maskXE = false;
677      maskXO = false;
678      maskYE = false;
679      maskYO = false;
680      //
681      lmax = 0.;
682      umax = 100.;
683      slmax = "";
684      sumax = "";
685      suf = "";
686      xyaverage = true;
687      //
688      letmax = -1;
689      heavytail = false;
690      //  lmipth = 100.;
691      lmipth = 0.;
692      //
693      lfit = 0;
694    };
695    
696    TF1 *CaloLong::GetFit(){
697    //  TString fnam=Form("lfit%s",suf.Data());
698    //  TF1 *lfit  = dynamic_cast<TF1*>(gDirectory->FindObject(fnam));
699      if ( lfit ) return lfit;
700      return NULL;
701    }
702    
703    void CaloLong::MaskSection(TString sec){
704      sec.ToUpper();
705      if ( sec.Contains("XO") ) maskXO = true;
706      if ( sec.Contains("YO") ) maskYO = true;
707      if ( sec.Contains("XE") ) maskXE = true;
708      if ( sec.Contains("YE") ) maskYE = true;
709    }
710    
711    void CaloLong::UnMaskSections(){
712      this->UnMaskSection("XEXOYEYO");
713    }
714    
715    void CaloLong::UnMaskSection(TString sec){
716      sec.ToUpper();
717      if ( sec.Contains("XO") ) maskXO = false;
718      if ( sec.Contains("YO") ) maskYO = false;
719      if ( sec.Contains("XE") ) maskXE = false;
720      if ( sec.Contains("YE") ) maskYE = false;
721    }
722    
723    void CaloLong::Clear(){
724      //
725      if ( debug ) printf(" Clear called \n");
726      //
727    memset(eplane,0, 2*22*sizeof(Float_t));    memset(eplane,0, 2*22*sizeof(Float_t));
728    //    //
729      chi2 = 0.;
730      ndf = 0.;
731      E0 = 0.;
732      defE0 = 0.;
733      a = 0.;
734      b = 0.;
735      errE0 = 0.;
736      erra = 0.;
737      errb = 0.;
738      etmax = 0.;
739      asymm = 0.;
740      fitresult = 0;
741      //
742      X0pl = 0.76;
743      //
744    };
745    
746    void CaloLong::Print(){
747      //
748      Fit();
749      //
750      printf("==================== Calorimeter Longitudinal Profile =======================\n");
751      printf(" OBT: %u PKT: %u ATIME: %u suf %s \n",OBT,PKT,atime,suf.Data());
752      printf(" fitresult:.. %i\n",fitresult);
753      printf(" chi2     :.. %f\n",chi2);
754      printf(" ndf      :.. %f\n",ndf);
755      printf(" nchi2    :.. %f\n",chi2/ndf);
756      printf(" E0       :.. %f\n",E0);
757      printf(" E0/260.  :.. %f\n",E0/260.);
758      printf(" defE0    :.. %f\n",defE0);
759      printf(" lower    :.. %f\n",lmax);
760      printf(" upper    :.. %f\n",umax);
761      printf(" s lower  :.. %s\n",slmax.Data());
762      printf(" s upper  :.. %s\n",sumax.Data());
763      printf(" a        :.. %f\n",a);
764      printf(" b        :.. %f\n",b);
765      printf(" errE0    :.. %f\n",errE0);
766      printf(" erra     :.. %f\n",erra);
767      printf(" errb     :.. %f\n",errb);
768      printf(" asymm    :.. %f\n",asymm);
769      printf(" tmax     :.. %f\n",((a-1.)/b));
770      printf(" etmax    :.. %f\n",etmax);
771      printf(" X0pl     :.. %f\n",X0pl);
772      printf("========================================================================\n");
773      //
774    };
775    
776    void CaloLong::SetNoWpreSampler(Int_t n){
777      //
778      if ( NC+n <= 22 && NC+n >= 0 ){
779        N = n;
780      } else {
781        printf(" ERROR! Calorimeter is made of 22 W planes\n");
782        printf(" you are giving N presampler = %i and N calo = %i \n",n,NC);
783        printf(" WARNING: using default values NWpre = 0, NWcalo = 22\n");
784        NC = 22;
785        N = 0;
786      };
787    }
788    
789    void CaloLong::SetNoWcalo(Int_t n){
790      if ( N+n <= 22 && N+n >= 0 ){
791        NC = n;
792      } else {
793        printf(" ERROR! Calorimeter is made of 22 W planes\n");
794        printf(" you are giving N W presampler = %i and N W calo = %i \n",N,n);
795        printf(" WARNING: using default values NWpre = 0, NWcalo = 22\n");
796        NC = 22;
797        N = 0;
798      };
799    }
800    
801    void CaloLong::SplitInto(Int_t NoWpreSampler, Int_t NoWcalo){
802      this->SetNoWpreSampler(0);
803      this->SetNoWcalo(0);
804      if ( NoWpreSampler < NoWcalo ){
805              this->SetNoWpreSampler(NoWpreSampler);
806              this->SetNoWcalo(NoWcalo);    
807      } else {
808              this->SetNoWcalo(NoWcalo);    
809              this->SetNoWpreSampler(NoWpreSampler);
810      };
811    }
812    
813    void CaloLong::Process(){
814      //  
815      if ( !L2 ){
816        printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");
817        printf(" ERROR: CaloHough variables not filled \n");
818        return;
819      };
820      //
821      Bool_t newentry = false;
822      //
823      if ( L2->IsORB() ){
824        if ( L2->GetOrbitalInfo()->pkt_num != PKT || L2->GetOrbitalInfo()->OBT != OBT || L2->GetOrbitalInfo()->absTime != atime ){
825          newentry = true;
826          OBT = L2->GetOrbitalInfo()->OBT;
827          PKT = L2->GetOrbitalInfo()->pkt_num;
828          atime = L2->GetOrbitalInfo()->absTime;
829        };
830      } else {
831        newentry = true;
832      };
833      //
834      if ( !newentry ) return;
835      //
836      if ( debug ) printf(" Start processing event at OBT %u PKT %u time %u \n",OBT,PKT,atime);
837      //
838      Clear();
839      //
840      // let's start
841      //
842      if ( cont ){
843        for (Int_t i=0; i<22; i++){
844          if ( i == (18+N) ){
845            mask18b =  18 + N;
846            break;
847          };
848        };
849      };  
850      //  
851      if ( sel ){
852        for (Int_t i=0; i<22; i++){
853          if ( i == (18-N) ){
854            mask18b =  18 - N;
855            break;
856          };
857        };
858      };
859      //
860      //  if ( mask18b == 18 ) mask18b = -1;
861      //
862      Int_t view = 0;
863      Int_t plane = 0;
864      Int_t strip = 0;
865      Float_t mip = 0.;
866      Bool_t gof = true;
867      for (Int_t i=0; i < L2->GetCaloLevel1()->istrip; i++){
868        mip = L2->GetCaloLevel1()->DecodeEstrip(i,view,plane,strip);
869        //
870        if ( !usepl18x && view==0 && plane==18 ) mip = 0.;
871        //
872        gof = true;
873        if ( maskXE && (plane%2)==0 && view==1 ) gof = false;
874        if ( maskXO && (plane%2)!=0 && view==1 ) gof = false;
875        if ( maskYE && (plane%2)!=0 && view==0 ) gof = false;
876        if ( maskYO && (plane%2)==0 && view==0 ) gof = false;
877        if ( gof ) eplane[view][plane] += mip;
878      };
879      //
880      // inclination factor (stolen from Daniele's code)
881      //
882      Float_t ytgx = 0.;
883      Float_t ytgy = 0.;
884      ytgx = 0.76 * clp->tanx[0];
885      ytgy = 0.76 * clp->tany[0];  
886      X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) );
887      //
888      // Find experimental plane of maximum
889      //
890      Int_t pmax = 0;
891      Int_t vmax = 0;
892      Float_t emax = 0.;
893    for (Int_t v=0; v<2; v++){    for (Int_t v=0; v<2; v++){
894     for (Int_t i=0; i<22; i++){     for (Int_t i=0; i<22; i++){
895       eplane[v][i] = L2->GetCaloLevel1()->qtotpl(v,i);       if ( eplane[v][i] > emax ){
896           emax = eplane[v][i];
897           vmax = v;
898           pmax = i;
899         };
900     };     };
901    };    };
902    //    //
903      //
904      //
905      if ( vmax == 0 ) pmax++;
906      etmax = pmax * X0pl;
907      //
908    if ( debug ) this->Print();    if ( debug ) this->Print();
909    if ( debug ) printf(" exit \n");    if ( debug ) printf(" exit \n");
910    //    //
911  };  };
912    
913    void CaloLong::SetEnergies(Float_t myene[][22]){
914      //  
915      if ( !L2 ){
916        printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");
917        printf(" ERROR: CaloHough variables not filled \n");
918        return;
919      };
920      //
921      Bool_t newentry = false;
922      //
923      if ( L2->IsORB() ){
924        if ( L2->GetOrbitalInfo()->pkt_num != PKT || L2->GetOrbitalInfo()->OBT != OBT || L2->GetOrbitalInfo()->absTime != atime ){
925          newentry = true;
926          OBT = L2->GetOrbitalInfo()->OBT;
927          PKT = L2->GetOrbitalInfo()->pkt_num;
928          atime = L2->GetOrbitalInfo()->absTime;
929        };
930      } else {
931        newentry = true;
932      };
933      //
934      if ( !newentry ) return;
935      //
936      if ( debug ) printf(" SET ENERGIES Start processing event at OBT %u PKT %u time %u \n",OBT,PKT,atime);
937      //
938      Clear();
939      //
940      // let's start
941      //
942      if ( cont ){
943        for (Int_t i=0; i<22; i++){
944          if ( i == (18+N) ){
945            mask18b =  18 + N;
946            break;
947          };
948        };
949      };  
950      //  
951      if ( sel ){
952        for (Int_t i=0; i<22; i++){
953          if ( i == (18-N) ){
954            mask18b =  18 - N;
955            break;
956          };
957        };
958      };
959      //
960      //  if ( mask18b == 18 ) mask18b = -1;
961      //
962      Int_t view = 0;
963      Int_t plane = 0;
964      Bool_t gof = true;
965      for (view=0; view < 2; view++){
966        for (plane=0; plane < 22; plane++){      
967          gof = true;
968          if ( maskXE && (plane%2)==0 && view==1 ) gof = false;
969          if ( maskXO && (plane%2)!=0 && view==1 ) gof = false;
970          if ( maskYE && (plane%2)!=0 && view==0 ) gof = false;
971          if ( maskYO && (plane%2)==0 && view==0 ) gof = false;
972          if ( gof ) eplane[view][plane] = myene[view][plane];
973          if ( debug ) printf(" XE %i XO %i YE %i YO %i eplane %i %i = %f ........... myene %i %i = %f\n", maskXE, maskXO, maskYE, maskYO,view,plane,eplane[view][plane],view,plane,myene[view][plane]);
974        };
975      };
976      //
977      if ( !clp && debug ) printf(" ERROR!! no CaloLevel2 obj!!\n");
978      //
979      // inclination factor (stolen from Daniele's code)
980      //
981      Float_t ytgx = 0.;
982      Float_t ytgy = 0.;
983      ytgx = 0.76 * clp->tanx[0];
984      ytgy = 0.76 * clp->tany[0];  
985      X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) );
986      if ( debug ) printf(" X0pl %f tanx %f tany %f \n",X0pl,ytgx,ytgy);
987      //
988      // Find experimental plane of maximum
989      //
990      Int_t pmax = 0;
991      Int_t vmax = 0;
992      Float_t emax = 0.;
993      for (Int_t v=0; v<2; v++){
994       for (Int_t i=0; i<22; i++){
995         if ( eplane[v][i] > emax ){
996           emax = eplane[v][i];
997           vmax = v;
998           pmax = i;
999         };
1000       };
1001      };
1002      //
1003      //
1004      //
1005      if ( vmax == 0 ) pmax++;
1006      etmax = pmax * X0pl;
1007      letmax = etmax;
1008      //
1009      if ( debug ) this->Print();
1010      if ( debug ) printf(" exit \n");
1011      //
1012    };
1013    
1014    Double_t ccurve(Double_t *ti,Double_t *par){
1015      //
1016      Double_t t = *ti;
1017      Double_t cE0 = par[0];
1018      Double_t ca = par[1];
1019      Double_t cb = par[2];
1020      Double_t gammaa = TMath::Gamma(ca);
1021      //
1022      Double_t value = cE0 * cb * ( pow((cb * t),(ca - 1)) * exp( -cb * t ) ) / gammaa;
1023      //
1024      return value;
1025      //
1026    }
1027    
1028    void CaloLong::Fit(){
1029      this->Fit(false);
1030    };
1031    
1032    Float_t CaloLong::Evaluate(TString s, Float_t tmax, Float_t X0pl){
1033      /* SAMPLE OUTPUT:
1034         Enter Infix Expression : A + B + C / (E - F)
1035         Postfix Expression is : A B + C E F - / +  
1036      */  
1037      if ( !s.Contains("tmax") && !s.Contains("X0pl") ){
1038        printf(" ERROR, the input formula must contain \"t\"\n");
1039        return 0.;
1040      };
1041      if ( tmax != tmax ){
1042        printf(" ERROR, tmax is nan! \n");
1043        return 0.;
1044      };
1045      TString g=Form("%f",tmax);
1046      TString h=Form("%f",X0pl);
1047      TString *ts= new TString("");
1048      ts->Prepend(s.Data());
1049      ts->ReplaceAll("tmax",4,g.Data(),g.Capacity());
1050      ts->ReplaceAll("X0pl",4,h.Data(),h.Capacity());
1051      ts->Prepend("(");
1052      ts->Append(")");
1053      if ( debug )  printf(" ts %s tssize %i capac %i s %s g %s \n",ts->Data(),ts->Sizeof(),ts->Capacity(),s.Data(),g.Data());
1054      char in[50],post[50];  
1055      strcpy(&in[0],"                                                 ");
1056      strcpy(&in[0],ts->Data());
1057      strcpy(&post[0],"                                                 ");
1058      if ( debug ) printf("Infix Expression is : ##%s##\n",&in[0]);
1059      infix2postfix(&in[0],&post[0],1);
1060      if ( debug ) printf("Postfix Expression is : ##%s##\n",&post[0]);
1061      /* SAMPLE OUTPUT:
1062         Enter Postfix Expression : 3 5 + 2 /
1063         3 5 + 2 / EQUALS 4
1064      */
1065      Float_t res = evaluate(&post[0]);
1066      if ( debug ) printf("%s EQUALS %f\n",post,res);
1067      return res;
1068    }
1069    
1070    void CaloLong::Fit(Bool_t draw){
1071      //
1072      Process();
1073      //
1074      //  
1075      if ( !L2 ){
1076        printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");
1077        printf(" ERROR: CaloHough variables not filled \n");
1078        return;
1079      };
1080      //
1081      Bool_t newentry = false;
1082      //
1083      if ( L2->IsORB() ){
1084        if ( L2->GetOrbitalInfo()->pkt_num != fPKT || L2->GetOrbitalInfo()->OBT != fOBT || L2->GetOrbitalInfo()->absTime != fatime ){
1085          newentry = true;
1086          fOBT = L2->GetOrbitalInfo()->OBT;
1087          fPKT = L2->GetOrbitalInfo()->pkt_num;
1088          fatime = L2->GetOrbitalInfo()->absTime;
1089        };
1090      } else {
1091        newentry = true;
1092      };
1093      //
1094      if ( !newentry ) return;
1095      //
1096      if ( debug ) printf(" Start fitting event at OBT %u PKT %u time %u \n",fOBT,fPKT,fatime);
1097      //
1098      if ( draw ){
1099        gStyle->SetLabelSize(0.04);
1100        gStyle->SetNdivisions(510,"XY");
1101      };
1102      //
1103      TString hid = Form("clongfit%s",suf.Data());  
1104      TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
1105      //  if ( tc ) tc->Delete();
1106      //  if ( tc ) tc->Close();
1107      if ( !tc && draw ){
1108        tc = new TCanvas(hid,hid);
1109      } else {
1110        if ( tc ) tc->cd();
1111      };
1112      //
1113      TString thid = Form("hlongfit%s",suf.Data());
1114      TH2F *th  = dynamic_cast<TH2F*>(gDirectory->FindObject(thid));
1115      if ( th ) th->Delete();
1116      //
1117      TString ghid = Form("glongfit%s",suf.Data());
1118      TGraphErrors *gh  = dynamic_cast<TGraphErrors*>(gDirectory->FindObject(ghid));
1119      if ( gh ) gh->Delete();
1120      //
1121      Float_t xpos = 0.;
1122      Float_t enemip = 0.;
1123      Float_t xmax = NC * X0pl + 0.2;
1124      //
1125      Double_t xxx[50];
1126      Double_t exx[50];
1127      Double_t yyy[50];
1128      Double_t eyy[50];
1129      Int_t numpo = 0;
1130      memset(xxx,0,50*sizeof(Double_t));
1131      memset(exx,0,50*sizeof(Double_t));
1132      memset(yyy,0,50*sizeof(Double_t));
1133      memset(eyy,0,50*sizeof(Double_t));
1134      //
1135      // AGH, BUG!
1136      //
1137      Int_t mmin = 0;
1138      Int_t mmax = 0;
1139      if ( cont ){
1140        mmin = N;
1141        mmax = NC+N;
1142      } else {
1143        mmin = 0;
1144        mmax = NC;
1145      };
1146      //
1147      Float_t emax = 0.;
1148      Float_t qtotparz = 0.;
1149      for (Int_t st=mmin;st<mmax+1;st++){
1150        enemip = 0.;
1151        xpos = (st - mmin) * X0pl;
1152        //
1153        if ( xyaverage ){
1154          //
1155          if ( st > mmin && st < mmax ){      
1156            if ( no18x && ( st == 18+1 || st == mask18b+1 )){
1157              if ( !maskYO ){
1158                enemip = 2. * eplane[1][st];
1159              } else {
1160                enemip = eplane[1][st];
1161              };
1162            } else {
1163              enemip = eplane[0][st-1] + eplane[1][st];
1164            };
1165          } else {
1166            if ( st == mmin ){
1167              if ( !maskYE ){
1168                enemip = 2. * eplane[1][st];
1169              } else {
1170                enemip = eplane[1][st];
1171              };
1172            };
1173            if ( st == mmax ){
1174              if ( !maskXE ){
1175                enemip = 2. * eplane[0][st-1];
1176              } else {
1177                enemip = eplane[0][st-1];
1178              };
1179            };
1180          };
1181          //
1182          qtotparz += enemip;
1183          if ( enemip > emax ) emax = enemip;
1184          if ( enemip > 0. ){
1185            xxx[numpo] = xpos;
1186            exx[numpo] = 0.1;
1187            yyy[numpo] = enemip;
1188            eyy[numpo] = sqrt(enemip*3.)+sqrt(5.);
1189            //if ( xpos > letmax && enemip > lmipth && heavytail) eyy[numpo] = (sqrt(enemip*3.)+sqrt(5.))/numpo;
1190            if ( xpos > letmax && enemip > lmipth && heavytail) eyy[numpo] = sqrt(enemip)/5.;
1191            if ( xpos > letmax-1 && xpos < letmax+1 && heavytail ) eyy[numpo] /= 5.;
1192            //if ( xpos > letmax-2 && xpos < letmax+1 && heavytail ) eyy[numpo] /= 7.;
1193            if ( xpos < 3. && heavytail ) eyy[numpo] /= 5.;
1194            numpo++;
1195            //      th->Fill(xpos,enemip);
1196            if ( debug ) printf(" AVE Filling: st %i xpos %f energy %f error %f letmax %f lmpith %f \n",st,xpos,enemip,eyy[numpo-1],letmax,lmipth);
1197          };
1198        } else {
1199          for (Int_t jj=0; jj<2; jj++){
1200            if ( st > mmin && st < mmax ){
1201              if ( jj == 0 && no18x && ( st == 18+1 || st == mask18b+1 )){
1202                enemip = eplane[1][st];
1203              } else {
1204                if ( jj == 0 ){
1205                  enemip = eplane[jj][st-1];
1206                } else {
1207                  enemip = eplane[jj][st];
1208                };        
1209              };
1210            } else {
1211              if ( st == mmin && jj == 1 ){
1212                enemip = eplane[jj][st];
1213              };
1214              if ( st == mmax && jj == 0){
1215                enemip = eplane[jj][st-1];
1216              };
1217            };
1218            //
1219            qtotparz += enemip;
1220          if ( enemip > emax ) emax = enemip;
1221            if ( enemip > 0. ){
1222              xxx[numpo] = xpos;
1223              exx[numpo] = 0.1;
1224              yyy[numpo] = enemip;
1225              eyy[numpo] = sqrt(enemip*3.)+sqrt(5.);
1226              //if ( xpos > letmax && enemip > lmipth && heavytail) eyy[numpo] = (sqrt(enemip*3.)+sqrt(5.))/numpo;
1227              //if ( xpos > letmax && enemip > lmipth && heavytail ) eyy[numpo] = (sqrt(enemip*3.)+sqrt(5.))/numpo;
1228              if ( xpos > letmax && enemip > lmipth && heavytail) eyy[numpo] = sqrt(enemip)/5.;
1229              if ( xpos > letmax-1 && xpos < letmax+1 && heavytail ) eyy[numpo] /= 5.;
1230            //if ( xpos > letmax-2 && xpos < letmax+1 && heavytail ) eyy[numpo] /= 7.;
1231              if ( xpos < 3. && heavytail ) eyy[numpo] /= 5.;        
1232              //      eyy[numpo] = sqrt(enemip)/(st*0.95);
1233              numpo++;
1234              //      th->Fill(xpos,enemip);
1235              if ( debug ) printf(" XY Filling: st %i xpos %f energy %f \n",st,xpos,enemip);
1236            };
1237          };        
1238        };
1239    
1240        //
1241        //    for (Int_t v=1; v>=0;v--)// {
1242    //       //
1243    //       if ( v == 1 ){
1244    //      xpos = (st - N) * X0pl;
1245    //       } else {
1246    //      xpos = (st + 1 - N) * X0pl;
1247    //       };
1248    //       //
1249    //       if ( no18x && st == 18 && v == 0 ){
1250    //      // skip plane 18x
1251    //       } else {
1252    //      if ( v == 1 && st == mask18b ){
1253    //        // emulate plane 18x
1254    //      } else {
1255    //        if ( eplane[v][st] > 0. ){
1256    //          th->Fill(xpos,eplane[v][st]);
1257    //          if ( debug ) printf(" Filling: st %i v %i xpos %f energy %f \n",st,v,xpos,eplane[v][st]);
1258    //        };
1259    //      };
1260    //       };
1261    //       //
1262    //     };
1263      };
1264      //
1265      //  th = new TH2F(thid,thid,int(NC*1.5),-0.2,xmax);
1266      th = new TH2F(thid,thid,1000,-0.2,xmax,1000,0.,emax*1.2);
1267      gh = new TGraphErrors(numpo,xxx,yyy,exx,eyy);
1268      TString fnam=Form("lfit%s",suf.Data());
1269    //  TF1 *lfit  = dynamic_cast<TF1*>(gDirectory->FindObject(fnam));
1270      if ( lfit ) lfit->Delete();
1271      lfit = new TF1(fnam,ccurve,0.,xmax,3);
1272    //   if ( !lfit ){
1273    //     lfit = new TF1(fnam,ccurve,0.,xmax,3);
1274    //   };
1275    //   lfit->Clear();
1276    //   lfit->SetName(fnam);
1277    //   lfit->SetRange(0.,xmax);
1278    //  lfit->Print();
1279      //
1280      if ( debug ) printf("qtot %f qtotparz %f \n",clp->qtot,qtotparz);
1281      E0 = qtotparz;
1282      //  E0 = clp->qtot;
1283      a = 5.;
1284      b = 0.5;
1285      if ( debug ) printf(" STARTING PARAMETERS: E0 %f a %f b %f \n",E0,a,b);
1286      lfit->SetParameters(E0,a,b);
1287      //  lfit->SetParLimits(0,0.,1000.);
1288      //  lfit->SetParLimits(1,-1.,80.);
1289      //  lfit->SetParLimits(2,-1.,10.);
1290      // TString optio = "ROW"; // "RO"
1291      TString optio = "RO"; // "RO"
1292      if ( debug ){
1293        optio += "NV";
1294        if ( draw ) optio += "V";
1295      } else {
1296        optio += "NQ";
1297        if ( draw ) optio += "Q";
1298      };
1299      //
1300      if ( debug ) printf(" OK, start the fitting procedure... options: %s \n",optio.Data());
1301      //
1302      //  fitresult = th->Fit("lfit",optio);
1303      //  lfit->Update();
1304      fitresult = gh->Fit("lfit",optio);
1305      //  fitresult = gh->Fit("lfit","ROQW");
1306      //
1307      if ( debug ) printf(" the fit is done! result: %i \n",fitresult);
1308      //
1309      E0 = lfit->GetParameter(0);
1310      a = lfit->GetParameter(1);
1311      b = lfit->GetParameter(2);
1312      errE0 = lfit->GetParError(0);
1313      erra = lfit->GetParError(1);
1314      errb = lfit->GetParError(2);
1315      chi2 = lfit->GetChisquare();
1316      ndf = lfit->GetNDF();
1317      Float_t tmax = 0.;
1318      if ( debug ) printf(" Parameters are retrieved \n");
1319      if ( b != 0 ) tmax = (a - 1.)/b;
1320      //
1321      if ( fitresult != 0 ){
1322        if ( debug ) printf(" The fit failed, no integrals calculation and asymm is set to -1. \n");
1323        asymm = -1.;
1324        defE0 = -1.;
1325      } else {
1326        if ( slmax.MaybeRegexp() ) lmax = this->Evaluate(slmax,tmax,X0pl);
1327        if ( sumax.MaybeRegexp() ) umax = this->Evaluate(sumax,tmax,X0pl);
1328        Int_t npp = 1000;
1329        double *xpp=new double[npp];
1330        double *wpp=new double[npp];
1331        lfit->CalcGaussLegendreSamplingPoints(npp,xpp,wpp,1e-12);
1332        defE0 = lfit->IntegralFast(npp,xpp,wpp,lmax,umax);
1333        //
1334        double *xp=new double[npp];
1335        double *wp=new double[npp];
1336        lfit->CalcGaussLegendreSamplingPoints(npp,xp,wp,1e-12);
1337        Float_t imax = lfit->IntegralFast(npp,xp,wp,0.,tmax);
1338        //    Float_t imax = lfit->Integral(0.,tmax);
1339        if ( debug ) printf(" Integral till maximum (%f): %f \n",tmax,imax);
1340        Int_t np = 1000;
1341        double *x=new double[np];
1342        double *w=new double[np];
1343        lfit->CalcGaussLegendreSamplingPoints(np,x,w,1e-12);
1344        Float_t i10max = lfit->IntegralFast(np,x,w,0.,10.*tmax);
1345        delete x;
1346        delete w;
1347        delete xp;
1348        delete wp;
1349        delete xpp;
1350        delete wpp;
1351        //    Float_t i10max = lfit->Integral(0.,10.*tmax);
1352        if ( debug ) printf(" Integral: %f \n",i10max);
1353        //
1354        if ( i10max != imax ){
1355          asymm = imax / (i10max-imax);
1356        } else {
1357          if ( debug ) printf(" i10max == imax, asymm undefined\n");
1358          asymm = -2.;
1359        };
1360        if ( asymm != asymm ){
1361          if ( debug ) printf(" asymm is nan \n");      
1362          asymm = -3.;
1363        };
1364        //lfit->Integral(0.,tmax)/(lfit->Integral(0.,10.*tmax)-lfit->Integral(0.,tmax));
1365        if ( debug ) printf(" Asymmetry has been calculated \n");
1366      };
1367      //
1368      if ( asymm < 0. || ndf <= 0. || chi2 < 0. || tmax < 0. ){
1369        if ( debug ) printf(" Funny asymm||ndf||chi2||tmax values, fit failed \n");
1370        fitresult = 100;
1371      };
1372      //
1373      if ( draw ){
1374        //
1375        tc->cd();    
1376        //    gStyle->SetOptStat(11111);
1377        tc->SetTitle();
1378        th->SetTitle("");
1379        th->SetName("");
1380        th->SetMarkerStyle(20);
1381        // axis titles
1382        th->SetXTitle("Depth [X0]");
1383        th->SetYTitle("Energy [MIP]");
1384        //    th->DrawCopy("Perror");
1385        th->DrawCopy();
1386        gh->SetMarkerSize(1.);
1387        gh->SetMarkerStyle(21);
1388        //    gh->SetMarkerColor(kRed);
1389        gh->Draw("Psame");
1390        lfit->Draw("same");
1391        tc->Modified();
1392        tc->Update();
1393        //
1394        gStyle->SetLabelSize(0);
1395        gStyle->SetNdivisions(1,"XY");
1396        //
1397      } else {
1398        if ( th ) th->Delete();
1399      };
1400      //
1401    //   gMinuit->SetPrintLevel(1);
1402    //   gMinuit->mnamin();
1403    //   //  gMinuit->mnfree(0);
1404    //   gMinuit->mncler();
1405    //   gMinuit->mnrset(1);
1406      //  gMinuit->Delete();
1407      //ROOT::Math::Minimizer::Clear();
1408      //  gMinuit->mnhelp("*");
1409      //
1410      //  delete lfit;
1411      //
1412    };
1413    
1414    void CaloLong::Draw(){
1415      //
1416      Process();
1417      Draw(-1);
1418    };
1419    
1420    void CaloLong::Draw(Int_t view){
1421      //
1422      Int_t minv = 0;
1423      Int_t maxv = 0;
1424      //
1425      if ( view == -1 ){
1426        maxv = -1;
1427      } else {
1428        minv = view;
1429        maxv = view+1;
1430      };
1431      //
1432      Process();
1433      //
1434      gStyle->SetLabelSize(0.04);
1435      gStyle->SetNdivisions(510,"XY");
1436      //
1437      if ( maxv != -1 ){
1438        for (Int_t v=minv; v<maxv;v++){
1439          TString hid = Form("clongv%i%s",v,suf.Data());    
1440          TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
1441          if ( tc ){
1442            //       tc->Clear();
1443          } else {
1444            tc = new TCanvas(hid,hid);
1445          };
1446          //
1447          TString thid = Form("hlongv%i%s",v,suf.Data());  
1448          TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));
1449          if ( th ) th->Delete();
1450          //         th->Clear();
1451          //         th->Reset();
1452          //        } else {
1453          th = new TH1F(thid,thid,22,-0.5,21.5);
1454          //        };
1455          tc->cd();
1456          //
1457          for (Int_t st=0;st<22;st++){
1458            th->Fill(st,eplane[v][st]);
1459          };
1460          th->Draw();
1461          tc->Modified();
1462          tc->Update();
1463        };
1464      } else {
1465        //
1466        TString hid = Form("clongvyvx%s",suf.Data());      
1467        TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
1468        if ( tc ){
1469        } else {
1470          tc = new TCanvas(hid,hid);
1471        };
1472        //
1473        TString thid = Form("hlongvyvx%s",suf.Data());      
1474        TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));
1475        if ( th ) th->Delete();
1476        th = new TH1F(thid,thid,44,-0.5,43.5);
1477        tc->cd();
1478        Int_t pp=0;
1479        for (Int_t st=0;st<22;st++){
1480          for (Int_t v=1; v>=0;v--){
1481            //
1482            th->Fill(pp,eplane[v][st]);
1483            //
1484            pp++;
1485          };
1486        };
1487        th->Draw();
1488        tc->Modified();
1489        tc->Update();
1490      };
1491      //
1492      gStyle->SetLabelSize(0);
1493      gStyle->SetNdivisions(1,"XY");
1494      //
1495    };
1496    
1497    void CaloLong::Delete(){
1498      Clear();
1499      //delete this;
1500    };
1501    
1502    
1503    
1504    
1505    
1506    
1507    
1508    
1509    
1510    
1511    
1512    
1513    
1514    
1515    
1516    
1517    
1518    
1519    
1520    
Legend:
Removed from v.1.1.1.1  
changed lines
  Added in v.1.19
  ViewVC Help
Powered by ViewVC 1.1.23