/[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 by mocchiut, Fri Nov 9 09:11:24 2007 UTC revision 1.2 by mocchiut, Mon Sep 22 20:08:25 2008 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  /**  /**
# Line 60  CaloLat::CaloLat(){ Line 64  CaloLat::CaloLat(){
64  /**  /**
65   * Default constructor   * Default constructor
66   */   */
67  CaloLong::CaloLong(){  Calo2D::Calo2D(){
68    Clear();    Clear();
69  };  };
70    
# Line 80  CaloLat::CaloLat(PamLevel2 *l2p){   Line 84  CaloLat::CaloLat(PamLevel2 *l2p){  
84    //    //
85  };  };
86    
87  CaloLong::CaloLong(PamLevel2 *l2p){    Calo2D::Calo2D(PamLevel2 *l2p){  
88    //    //
89    Clear();    Clear();
90    //    //
# Line 100  void CaloLat::Clear(){ Line 104  void CaloLat::Clear(){
104    //    //
105  };  };
106    
107  void CaloLong::Clear(){  void Calo2D::Clear(){
108    //    //
109  };  };
110    
# Line 116  void CaloLat::Print(){ Line 120  void CaloLat::Print(){
120    //    //
121  };  };
122    
123  void CaloLong::Print(){  void Calo2D::Print(){
124    //    //
125    Process();    Process();
126    //    //
127    printf("==================== Calorimeter Lateral Profile =======================\n");    printf("==================== Calorimeter 2D Profile =======================\n");
128    printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);    printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);
129  //  printf(" nx [number of X combination]:.. %i\n",nx);  //  printf(" nx [number of X combination]:.. %i\n",nx);
130  //  printf(" ny [number of Y combination]:.. %i\n",ny);  //  printf(" ny [number of Y combination]:.. %i\n",ny);
# Line 134  void CaloLat::Draw(){ Line 138  void CaloLat::Draw(){
138    Draw(-1,-1);    Draw(-1,-1);
139  };  };
140    
141  void CaloLong::Draw(){  void Calo2D::Draw(){
142    //    //
143    Process();    Process();
144    Draw(-1);    Draw(-1);
145  };  };
146    
   
147  void CaloLat::Draw(Int_t view,Int_t plane){  void CaloLat::Draw(Int_t view,Int_t plane){
148    //    //
149    Int_t minv = 0;    Int_t minv = 0;
# Line 203  void CaloLat::Draw(Int_t view,Int_t plan Line 206  void CaloLat::Draw(Int_t view,Int_t plan
206    //    //
207  };  };
208    
209  void CaloLong::Draw(Int_t view){  void Calo2D::Draw(Int_t plane){
210    //    //
211    Int_t minv = 0;    Int_t minp = 0;
212    Int_t maxv = 0;    Int_t maxp = 0;
213    //    //
214    if ( view == -1 ){    if ( plane == -1 ){
215          maxv = -1;          minp = 0;
216            maxp = 23;
217    } else {    } else {
218          minv = view;          minp = plane;
219          maxv = view+1;          maxp = plane+1;
220    };    };
221    //    //
222    Process();    Process();
# Line 220  void CaloLong::Draw(Int_t view){ Line 224  void CaloLong::Draw(Int_t view){
224    gStyle->SetLabelSize(0.04);    gStyle->SetLabelSize(0.04);
225    gStyle->SetNdivisions(510,"XY");    gStyle->SetNdivisions(510,"XY");
226    //    //
227    if ( maxv != -1 ){    for (Int_t p=minp; p<maxp;p++){
228    for (Int_t v=minv; v<maxv;v++){      TString hid = Form("c2dp%i",p);    
229          TString hid = Form("clongv%i",v);            TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
230          TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));      if ( tc ){
231          if ( tc ){        //         tc->Clear();
232  //       tc->Clear();      } else {
233          } else {        tc = new TCanvas(hid,hid);
234           tc = new TCanvas(hid,hid);      };
235          };      //
236          //      TString thid = Form("h2dp%i",p);    
237          TString thid = Form("hlongv%i",v);            TH2F *th  = dynamic_cast<TH2F*>(gDirectory->FindObject(thid));
238          TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));      if ( th ) th->Delete();
239          if ( th ) th->Delete();      //   th->Clear();
240  //       th->Clear();      //   th->Reset();
241  //       th->Reset();      //  } else {
242  //      } else {      Int_t minx = smax[p] - 10;
243           th = new TH1F(thid,thid,22,-0.5,21.5);      if ( minx < 0 ) minx = 0;
244  //      };      Int_t maxx = minx + 20;
245          tc->cd();      if ( maxx > 95 ){
246          //        maxx = 95;
247          for (Int_t st=0;st<22;st++){        minx = 75;
248           th->Fill(st,eplane[v][st]);      };
249          };      Int_t miny = smay[p] - 10;
250          th->Draw();      if ( miny < 0 ) miny = 0;
251          tc->Modified();      Int_t maxy = miny + 20;
252          tc->Update();      if ( maxy > 95 ){
253    };        maxy = 95;
254    } else {        miny = 75;
255          //      };
256          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);
257          TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));      //    th = new TH2F(thid,thid,96,-0.5,95.5,96,-0.5,95.5);
258          if ( tc ){      //  };
259          } else {      tc->cd();
260           tc = new TCanvas(hid,hid);      //
261          };      for (Int_t stx=minx;stx<maxx+1;stx++){
262          //        for (Int_t sty=miny;sty<maxy+1;sty++){
263          TString thid = Form("hlongvyvx");                th->Fill(stx,sty,estrip[p][stx][sty]);
264          TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));        };
265          if ( th ) th->Delete();      };
266          th = new TH1F(thid,thid,44,-0.5,43.5);      gStyle->SetPalette(1);
267          tc->cd();      //    tc->SetLogz();    
268          Int_t pp=0;      //    th->Draw("colbox");
269          for (Int_t st=0;st<22;st++){      th->Draw("cont4");
270            for (Int_t v=1; v>=0;v--){      tc->Modified();
271                  //      tc->Update();
                 th->Fill(pp,eplane[v][st]);  
                 //  
                 pp++;  
           };  
         };  
         th->Draw();  
         tc->Modified();  
         tc->Update();  
272    };    };
273    //    //
274    gStyle->SetLabelSize(0);    gStyle->SetLabelSize(0);
# Line 285  void CaloLat::Delete(){ Line 281  void CaloLat::Delete(){
281    //delete this;    //delete this;
282  };  };
283    
284  void CaloLong::Delete(){  void Calo2D::Delete(){
285    Clear();    Clear();
286    //delete this;    //delete this;
287  };  };
288    
289    
290  void CaloLat::Process(){  void CaloLat::Process(){
291    //      //  
292    if ( !L2 ){    if ( !L2 ){
# Line 335  void CaloLat::Process(){ Line 332  void CaloLat::Process(){
332    //    //
333  };  };
334    
335  void CaloLong::Process(){  void Calo2D::Process(){
336    //      //  
337    if ( !L2 ){    if ( !L2 ){
338      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 361  void CaloLong::Process(){
361    //    //
362    // let's start    // let's start
363    //    //
364      Float_t es[2][22][96];
365      memset(es,0, 4224*sizeof(Float_t));
366      memset(estrip,0, 4224*sizeof(Float_t));
367      Float_t mip1 = 0.;
368      Int_t view1 = 0;
369      Int_t plane1 = 0;
370      Int_t strip1 = 0;
371      //
372      for (Int_t i=0; i<L2->GetCaloLevel1()->istrip ; i++){
373        mip1 = L2->GetCaloLevel1()->DecodeEstrip(i,view1,plane1,strip1);
374        es[view1][plane1][strip1] = mip1;
375      };
376      //
377      Int_t plane2 = 0;
378      Float_t emax[23];
379      memset(emax,0,sizeof(Float_t)*23);
380      memset(smax,0,sizeof(Int_t)*23);
381      memset(smay,0,sizeof(Int_t)*23);
382      //
383      //
384      for (Int_t p=0; p < 23 ; p++){
385        //
386        plane1 = p-1;
387        plane2 = p;
388        //
389        if ( p == 0 ){
390          plane1 = -1;
391          plane2 = 0;
392        };
393        if ( p == 22 ){
394          plane1 = 21;
395          plane2 = -1;
396        };
397        //
398        for (Int_t s=0; s < 96 ; s++){ // x
399          for (Int_t ss=0; ss < 96 ; ss++){ // y
400            if ( p == 0 ){
401              estrip[p][s][ss] += es[1][plane2][ss];
402              if ( (es[1][plane2][ss]) > emax[p] ){
403                smax[p] = 45;
404                smay[p] = ss;
405                emax[p] = es[1][plane2][ss] ;
406              };
407            };
408            if ( p > 0 && p < 22 ){
409              estrip[p][s][ss] += es[0][plane1][s] + es[1][plane2][ss];
410              if ( (es[0][plane1][s] + es[1][plane2][ss]) > emax[p] ){
411                smax[p] = s;
412                smay[p] = ss;
413                emax[p] = es[0][plane1][s] + es[1][plane2][ss] ;
414              };
415            };      
416            if ( p == 22 ){
417              estrip[p][s][ss] += es[0][plane1][s];
418              if ( (es[1][plane2][s]) > emax[p] ){
419                smax[p] = s;
420                smay[p] = 45;
421                emax[p] = es[1][plane2][s] ;
422              };
423            };
424          };
425        };
426        //
427      };
428      //
429      if ( debug ) this->Print();
430      if ( debug ) printf(" exit \n");
431      //
432    };
433    
434    
435    /**
436     * Default constructor
437     */
438    CaloLong::CaloLong(){
439      Clear();
440    };
441    
442    CaloLong::CaloLong(PamLevel2 *l2p){  
443      //
444      Clear();
445      //
446      L2 = l2p;
447      //
448      if ( !L2->IsORB() ) printf(" WARNING: OrbitalInfo Tree is needed, the plugin could not work properly without it \n");
449      //
450      OBT = 0;
451      PKT = 0;
452      atime = 0;
453      //
454      sel = true;
455      cont = false;
456      N = 0;
457      NC = 22;
458      mask18b = -1;
459      //
460      no18x = true;
461      debug = false;
462      //
463    };
464    
465    void CaloLong::Clear(){
466      //
467    memset(eplane,0, 2*22*sizeof(Float_t));    memset(eplane,0, 2*22*sizeof(Float_t));
468    //    //
469      chi2 = 0.;
470      ndf = 0.;
471      E0 = 0.;
472      a = 0.;
473      b = 0.;
474      errE0 = 0.;
475      erra = 0.;
476      errb = 0.;
477      etmax = 0.;
478      asymm = 0.;
479      fitresult = 0;
480      //
481      X0pl = 0.76;
482      //
483    };
484    
485    void CaloLong::Print(){
486      //
487      Process();
488      //
489      printf("==================== Calorimeter Longitudinal Profile =======================\n");
490      printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);
491      printf(" fitresult:.. %i\n",fitresult);
492      printf(" chi2     :.. %f\n",chi2);
493      printf(" ndf      :.. %f\n",ndf);
494      printf(" nchi2    :.. %f\n",chi2/ndf);
495      printf(" E0       :.. %f\n",E0);
496      printf(" E0/260.  :.. %f\n",E0/260.);
497      printf(" a        :.. %f\n",a);
498      printf(" b        :.. %f\n",b);
499      printf(" errE0    :.. %f\n",errE0);
500      printf(" erra     :.. %f\n",erra);
501      printf(" errb     :.. %f\n",errb);
502      printf(" asymm    :.. %f\n",asymm);
503      printf(" tmax     :.. %f\n",((a-1.)/b));
504      printf(" etmax    :.. %f\n",etmax);
505      printf(" X0pl     :.. %f\n",X0pl);
506      printf("========================================================================\n");
507      //
508    };
509    
510    void CaloLong::SetNoWpreSampler(Int_t n){
511      //
512      if ( NC+n <= 22 && NC+n >= 0 ){
513        N = n;
514      } else {
515        printf(" ERROR! Calorimeter is made of 22 W planes\n");
516        printf(" you are giving N presampler = %i and N calo = %i \n",n,NC);
517        printf(" WARNING: using default values NWpre = 0, NWcalo = 22\n");
518        NC = 22;
519        N = 0;
520      };
521    }
522    
523    void CaloLong::SetNoWcalo(Int_t n){
524      if ( N+n <= 22 && N+n >= 0 ){
525        NC = n;
526      } else {
527        printf(" ERROR! Calorimeter is made of 22 W planes\n");
528        printf(" you are giving N W presampler = %i and N W calo = %i \n",N,n);
529        printf(" WARNING: using default values NWpre = 0, NWcalo = 22\n");
530        NC = 22;
531        N = 0;
532      };
533    }
534    
535    void CaloLong::Process(){
536      //  
537      if ( !L2 ){
538        printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");
539        printf(" ERROR: CaloHough variables not filled \n");
540        return;
541      };
542      //
543      Bool_t newentry = false;
544      //
545      if ( L2->IsORB() ){
546        if ( L2->GetOrbitalInfo()->pkt_num != PKT || L2->GetOrbitalInfo()->OBT != OBT || L2->GetOrbitalInfo()->absTime != atime ){
547          newentry = true;
548          OBT = L2->GetOrbitalInfo()->OBT;
549          PKT = L2->GetOrbitalInfo()->pkt_num;
550          atime = L2->GetOrbitalInfo()->absTime;
551        };
552      } else {
553        newentry = true;
554      };
555      //
556      if ( !newentry ) return;
557      //
558      if ( debug ) printf(" Start processing event at OBT %u PKT %u time %u \n",OBT,PKT,atime);
559      //
560      Clear();
561      //
562      // let's start
563      //
564      if ( cont ){
565        for (Int_t i=0; i<22; i++){
566          if ( i == (18+N) ){
567            mask18b =  18 + N;
568            break;
569          };
570        };
571      };  
572      //  
573      if ( sel ){
574        for (Int_t i=0; i<22; i++){
575          if ( i == (18-N) ){
576            mask18b =  18 - N;
577            break;
578          };
579        };
580      };
581      //
582      //  if ( mask18b == 18 ) mask18b = -1;
583      //
584      Int_t view = 0;
585      Int_t plane = 0;
586      Int_t strip = 0;
587      Float_t mip = 0.;
588      for (Int_t i=0; i < L2->GetCaloLevel1()->istrip; i++){
589        mip = L2->GetCaloLevel1()->DecodeEstrip(i,view,plane,strip);
590        eplane[view][plane] += mip;
591      };
592      //
593      // inclination factor (steal from Daniele's code)
594      //
595      Float_t ytgx = 0;
596      Float_t ytgy = 0;
597      ytgx = 0.76 * L2->GetCaloLevel2()->tanx[0];
598      ytgy = 0.76 * L2->GetCaloLevel2()->tany[0];  
599      X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) );
600      //
601      // Find experimental plane of maximum
602      //
603      Int_t pmax = 0;
604      Int_t vmax = 0;
605      Float_t emax = 0.;
606    for (Int_t v=0; v<2; v++){    for (Int_t v=0; v<2; v++){
607     for (Int_t i=0; i<22; i++){     for (Int_t i=0; i<22; i++){
608       eplane[v][i] = L2->GetCaloLevel1()->qtotpl(v,i);       if ( eplane[v][i] > emax ){
609           emax = eplane[v][i];
610           vmax = v;
611           pmax = i;
612         };
613     };     };
614    };    };
615    //    //
616      //
617      //
618      if ( vmax == 0 ) pmax++;
619      etmax = pmax * X0pl;
620      //
621    if ( debug ) this->Print();    if ( debug ) this->Print();
622    if ( debug ) printf(" exit \n");    if ( debug ) printf(" exit \n");
623    //    //
624  };  };
625    
626    
627    Double_t ccurve(Double_t *ti,Double_t *par){
628      //
629      Double_t t = *ti;
630      Double_t cE0 = par[0];
631      Double_t ca = par[1];
632      Double_t cb = par[2];
633      Double_t gammaa = TMath::Gamma(ca);
634      //
635      Double_t value = cE0 * cb * ( pow((cb * t),(ca - 1)) * exp( -cb * t ) ) / gammaa;
636      //
637      return value;
638      //
639    }
640    
641    void CaloLong::Fit(){
642      this->Fit(false);
643    };
644    
645    void CaloLong::Fit(Bool_t draw){
646      //
647      Process();
648      //
649      //  
650      if ( !L2 ){
651        printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");
652        printf(" ERROR: CaloHough variables not filled \n");
653        return;
654      };
655      //
656      Bool_t newentry = false;
657      //
658      if ( L2->IsORB() ){
659        if ( L2->GetOrbitalInfo()->pkt_num != fPKT || L2->GetOrbitalInfo()->OBT != fOBT || L2->GetOrbitalInfo()->absTime != fatime ){
660          newentry = true;
661          fOBT = L2->GetOrbitalInfo()->OBT;
662          fPKT = L2->GetOrbitalInfo()->pkt_num;
663          fatime = L2->GetOrbitalInfo()->absTime;
664        };
665      } else {
666        newentry = true;
667      };
668      //
669      if ( !newentry ) return;
670      //
671      if ( debug ) printf(" Start fitting event at OBT %u PKT %u time %u \n",fOBT,fPKT,fatime);
672      //
673      if ( draw ){
674        gStyle->SetLabelSize(0.04);
675        gStyle->SetNdivisions(510,"XY");
676      };
677      //
678      TString hid = Form("clongfit");      
679      TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
680      //  if ( tc ) tc->Delete();
681      //  if ( tc ) tc->Close();
682      if ( !tc && draw ){
683        tc = new TCanvas(hid,hid);
684      } else {
685        if ( tc ) tc->cd();
686      };
687      //
688      TString thid = Form("hlongfit");      
689      TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));
690      if ( th ) th->Delete();
691      Float_t xpos = 0.;
692      Float_t enemip = 0.;
693      Float_t xmax = NC * X0pl + 0.2;
694      //  th = new TH1F(thid,thid,int(NC*1.5),-0.2,xmax);
695      th = new TH1F(thid,thid,100,-0.2,xmax);
696      //
697      for (Int_t st=N;st<(N+NC);st++){
698        enemip = 0.;
699        xpos = (st - N) * X0pl;
700        if ( st > N && st < (N+NC-1) ){      
701          if ( no18x && ( st == 18+1 || st == mask18b+1 )){
702            enemip = 2. * eplane[1][st];
703          } else {
704            enemip = eplane[0][st-1] + eplane[1][st];
705          };
706        } else {
707          if ( st == N ) enemip = 2. * eplane[1][st];
708          if ( st == (N+NC-1) ) enemip = 2. * eplane[0][st];
709        };
710        //
711        if ( enemip > 0. ){
712          th->Fill(xpos,enemip);
713          if ( debug ) printf(" Filling: st %i xpos %f energy %f \n",st,xpos,enemip);
714        };
715        //
716        //    for (Int_t v=1; v>=0;v--)// {
717    //       //
718    //       if ( v == 1 ){
719    //      xpos = (st - N) * X0pl;
720    //       } else {
721    //      xpos = (st + 1 - N) * X0pl;
722    //       };
723    //       //
724    //       if ( no18x && st == 18 && v == 0 ){
725    //      // skip plane 18x
726    //       } else {
727    //      if ( v == 1 && st == mask18b ){
728    //        // emulate plane 18x
729    //      } else {
730    //        if ( eplane[v][st] > 0. ){
731    //          th->Fill(xpos,eplane[v][st]);
732    //          if ( debug ) printf(" Filling: st %i v %i xpos %f energy %f \n",st,v,xpos,eplane[v][st]);
733    //        };
734    //      };
735    //       };
736    //       //
737    //     };
738      };
739      //
740      TF1 *lfit = new TF1("lfit",ccurve,0.,xmax,3);
741      E0 = L2->GetCaloLevel2()->qtot;
742      a = 5.;
743      b = 0.5;
744      if ( debug ) printf(" STARTING PARAMETERS: E0 %f a %f b %f \n",E0,a,b);
745      lfit->SetParameters(E0,a,b);
746      //  lfit->SetParLimits(0,0.,1000.);
747      //  lfit->SetParLimits(1,-1.,80.);
748      //  lfit->SetParLimits(2,-1.,10.);
749      TString optio;
750      if ( debug ){
751        //    optio = "MERBOV";
752        //    optio = "MEROV";
753        //    optio = "EROV";
754        optio = "RNOV";
755        if ( draw ) optio = "ROV";
756      } else {
757        //    optio = "MERNOQ";
758        //    optio = "ERNOQ";
759        optio = "RNOQ";
760        if ( draw ) optio = "ROQ";
761      };
762      //
763      if ( debug ) printf(" OK, start the fitting procedure...\n");
764      //
765      fitresult = th->Fit("lfit",optio);
766      //
767      if ( debug ) printf(" the fit is done! result: %i \n",fitresult);
768      //
769      E0 = lfit->GetParameter(0);
770      a = lfit->GetParameter(1);
771      b = lfit->GetParameter(2);
772      errE0 = lfit->GetParError(0);
773      erra = lfit->GetParError(1);
774      errb = lfit->GetParError(2);
775      chi2 = lfit->GetChisquare();
776      ndf = lfit->GetNDF();
777      Float_t tmax = 0.;
778      if ( debug ) printf(" Parameters are retrieved \n");
779      if ( b != 0 ) tmax = (a - 1.)/b;
780      //
781      if ( fitresult != 0 ){
782        if ( debug ) printf(" The fit failed, no integrals calculation and asymm is set to -1. \n");
783        asymm = -1.;
784      } else {
785        Int_t npp = 1000;
786        double *xp=new double[npp];
787        double *wp=new double[npp];
788        lfit->CalcGaussLegendreSamplingPoints(npp,xp,wp,1e-12);
789        Float_t imax = lfit->IntegralFast(npp,xp,wp,0.,tmax);
790        //    Float_t imax = lfit->Integral(0.,tmax);
791        if ( debug ) printf(" Integral till maximum (%f): %f \n",tmax,imax);
792        Int_t np = 1000;
793        double *x=new double[np];
794        double *w=new double[np];
795        lfit->CalcGaussLegendreSamplingPoints(np,x,w,1e-12);
796        Float_t i10max = lfit->IntegralFast(np,x,w,0.,10.*tmax);
797        delete x;
798        delete w;
799        delete xp;
800        delete wp;
801        //    Float_t i10max = lfit->Integral(0.,10.*tmax);
802        if ( debug ) printf(" Integral: %f \n",i10max);
803        //
804        if ( i10max != imax ){
805          asymm = imax / (i10max-imax);
806        } else {
807          if ( debug ) printf(" i10max == imax, asymm undefined\n");
808          asymm = -2.;
809        };
810        //lfit->Integral(0.,tmax)/(lfit->Integral(0.,10.*tmax)-lfit->Integral(0.,tmax));
811        if ( debug ) printf(" Asymmetry has been calculated \n");
812      };
813      //
814      if ( draw ){
815        //
816        tc->cd();    
817        //    gStyle->SetOptStat(11111);
818        tc->SetTitle();
819        th->SetTitle("");
820        th->SetName("");
821        th->SetMarkerStyle(20);
822        // axis titles
823        th->SetXTitle("Depth [X0]");
824        th->SetYTitle("Energy [MIP]");
825        th->DrawCopy("Perror");
826        lfit->Draw("same");
827        tc->Modified();
828        tc->Update();
829        //
830        gStyle->SetLabelSize(0);
831        gStyle->SetNdivisions(1,"XY");
832        //
833      };
834      //
835      delete lfit;
836      //
837    };
838    
839    void CaloLong::Draw(){
840      //
841      Process();
842      Draw(-1);
843    };
844    
845    void CaloLong::Draw(Int_t view){
846      //
847      Int_t minv = 0;
848      Int_t maxv = 0;
849      //
850      if ( view == -1 ){
851        maxv = -1;
852      } else {
853        minv = view;
854        maxv = view+1;
855      };
856      //
857      Process();
858      //
859      gStyle->SetLabelSize(0.04);
860      gStyle->SetNdivisions(510,"XY");
861      //
862      if ( maxv != -1 ){
863        for (Int_t v=minv; v<maxv;v++){
864          TString hid = Form("clongv%i",v);
865          TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
866          if ( tc ){
867            //       tc->Clear();
868          } else {
869            tc = new TCanvas(hid,hid);
870          };
871          //
872          TString thid = Form("hlongv%i",v);        
873          TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));
874          if ( th ) th->Delete();
875          //         th->Clear();
876          //         th->Reset();
877          //        } else {
878          th = new TH1F(thid,thid,22,-0.5,21.5);
879          //        };
880          tc->cd();
881          //
882          for (Int_t st=0;st<22;st++){
883            th->Fill(st,eplane[v][st]);
884          };
885          th->Draw();
886          tc->Modified();
887          tc->Update();
888        };
889      } else {
890        //
891        TString hid = Form("clongvyvx");    
892        TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
893        if ( tc ){
894        } else {
895          tc = new TCanvas(hid,hid);
896        };
897        //
898        TString thid = Form("hlongvyvx");  
899        TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));
900        if ( th ) th->Delete();
901        th = new TH1F(thid,thid,44,-0.5,43.5);
902        tc->cd();
903        Int_t pp=0;
904        for (Int_t st=0;st<22;st++){
905          for (Int_t v=1; v>=0;v--){
906            //
907            th->Fill(pp,eplane[v][st]);
908            //
909            pp++;
910          };
911        };
912        th->Draw();
913        tc->Modified();
914        tc->Update();
915      };
916      //
917      gStyle->SetLabelSize(0);
918      gStyle->SetNdivisions(1,"XY");
919      //
920    };
921    
922    void CaloLong::Delete(){
923      Clear();
924      //delete this;
925    };
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.2
  ViewVC Help
Powered by ViewVC 1.1.23