/[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.13 by mocchiut, Wed Aug 19 07:25:34 2009 UTC revision 1.19 by mocchiut, Tue May 18 04:03:29 2010 UTC
# Line 283  CaloLat::CaloLat(PamLevel2 *l2p){   Line 283  CaloLat::CaloLat(PamLevel2 *l2p){  
283    //    //
284    suf = "";    suf = "";
285    debug = false;    debug = false;
286      usepl18x = false;
287    //    //
288  };  };
289    
# Line 300  Calo2D::Calo2D(PamLevel2 *l2p){   Line 301  Calo2D::Calo2D(PamLevel2 *l2p){  
301    //    //
302    suf = "";    suf = "";
303    debug = false;    debug = false;
304      usepl18x = false;
305    //    //
306  };  };
307    
# Line 527  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 574  void Calo2D::Process(){ Line 579  void Calo2D::Process(){
579    //    //
580    for (Int_t i=0; i<L2->GetCaloLevel1()->istrip ; i++){    for (Int_t i=0; i<L2->GetCaloLevel1()->istrip ; i++){
581      mip1 = L2->GetCaloLevel1()->DecodeEstrip(i,view1,plane1,strip1);      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;      es[view1][plane1][strip1] = mip1;
586    };    };
587    //    //
# Line 663  CaloLong::CaloLong(PamLevel2 *l2p){   Line 671  CaloLong::CaloLong(PamLevel2 *l2p){  
671    mask18b = -1;    mask18b = -1;
672    //    //
673    no18x = true;    no18x = true;
674      usepl18x = !no18x;
675    debug = false;    debug = false;
676    maskXE = false;    maskXE = false;
677    maskXO = false;    maskXO = false;
# Line 676  CaloLong::CaloLong(PamLevel2 *l2p){   Line 685  CaloLong::CaloLong(PamLevel2 *l2p){  
685    suf = "";    suf = "";
686    xyaverage = true;    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){  void CaloLong::MaskSection(TString sec){
704    sec.ToUpper();    sec.ToUpper();
705    if ( sec.Contains("XO") ) maskXO = true;    if ( sec.Contains("XO") ) maskXO = true;
# Line 700  void CaloLong::UnMaskSection(TString sec Line 722  void CaloLong::UnMaskSection(TString sec
722    
723  void CaloLong::Clear(){  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.;    chi2 = 0.;
# Line 724  void CaloLong::Print(){ Line 748  void CaloLong::Print(){
748    Fit();    Fit();
749    //    //
750    printf("==================== Calorimeter Longitudinal Profile =======================\n");    printf("==================== Calorimeter Longitudinal Profile =======================\n");
751    printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);    printf(" OBT: %u PKT: %u ATIME: %u suf %s \n",OBT,PKT,atime,suf.Data());
752    printf(" fitresult:.. %i\n",fitresult);    printf(" fitresult:.. %i\n",fitresult);
753    printf(" chi2     :.. %f\n",chi2);    printf(" chi2     :.. %f\n",chi2);
754    printf(" ndf      :.. %f\n",ndf);    printf(" ndf      :.. %f\n",ndf);
# Line 842  void CaloLong::Process(){ Line 866  void CaloLong::Process(){
866    Bool_t gof = true;    Bool_t gof = true;
867    for (Int_t i=0; i < L2->GetCaloLevel1()->istrip; i++){    for (Int_t i=0; i < L2->GetCaloLevel1()->istrip; i++){
868      mip = L2->GetCaloLevel1()->DecodeEstrip(i,view,plane,strip);      mip = L2->GetCaloLevel1()->DecodeEstrip(i,view,plane,strip);
869        //
870        if ( !usepl18x && view==0 && plane==18 ) mip = 0.;
871        //
872      gof = true;      gof = true;
873      if ( maskXE && (plane%2)==0 && view==1 ) gof = false;      if ( maskXE && (plane%2)==0 && view==1 ) gof = false;
874      if ( maskXO && (plane%2)!=0 && view==1 ) gof = false;      if ( maskXO && (plane%2)!=0 && view==1 ) gof = false;
# Line 852  void CaloLong::Process(){ Line 879  void CaloLong::Process(){
879    //    //
880    // inclination factor (stolen from Daniele's code)    // inclination factor (stolen from Daniele's code)
881    //    //
882    Float_t ytgx = 0;    Float_t ytgx = 0.;
883    Float_t ytgy = 0;    Float_t ytgy = 0.;
884    ytgx = 0.76 * clp->tanx[0];    ytgx = 0.76 * clp->tanx[0];
885    ytgy = 0.76 * clp->tany[0];      ytgy = 0.76 * clp->tany[0];  
886    X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) );    X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) );
# Line 947  void CaloLong::SetEnergies(Float_t myene Line 974  void CaloLong::SetEnergies(Float_t myene
974      };      };
975    };    };
976    //    //
977      if ( !clp && debug ) printf(" ERROR!! no CaloLevel2 obj!!\n");
978      //
979    // inclination factor (stolen from Daniele's code)    // inclination factor (stolen from Daniele's code)
980    //    //
981    Float_t ytgx = 0;    Float_t ytgx = 0.;
982    Float_t ytgy = 0;    Float_t ytgy = 0.;
983    ytgx = 0.76 * clp->tanx[0];    ytgx = 0.76 * clp->tanx[0];
984    ytgy = 0.76 * clp->tany[0];      ytgy = 0.76 * clp->tany[0];  
985    X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) );    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    // Find experimental plane of maximum
989    //    //
# Line 974  void CaloLong::SetEnergies(Float_t myene Line 1004  void CaloLong::SetEnergies(Float_t myene
1004    //    //
1005    if ( vmax == 0 ) pmax++;    if ( vmax == 0 ) pmax++;
1006    etmax = pmax * X0pl;    etmax = pmax * X0pl;
1007      letmax = etmax;
1008    //    //
1009    if ( debug ) this->Print();    if ( debug ) this->Print();
1010    if ( debug ) printf(" exit \n");    if ( debug ) printf(" exit \n");
# Line 1155  void CaloLong::Fit(Bool_t draw){ Line 1186  void CaloLong::Fit(Bool_t draw){
1186          exx[numpo] = 0.1;          exx[numpo] = 0.1;
1187          yyy[numpo] = enemip;          yyy[numpo] = enemip;
1188          eyy[numpo] = sqrt(enemip*3.)+sqrt(5.);          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++;          numpo++;
1195          //      th->Fill(xpos,enemip);          //      th->Fill(xpos,enemip);
1196          if ( debug ) printf(" AVE Filling: st %i xpos %f energy %f \n",st,xpos,enemip);          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 {      } else {
1199        for (Int_t jj=0; jj<2; jj++){        for (Int_t jj=0; jj<2; jj++){
# Line 1187  void CaloLong::Fit(Bool_t draw){ Line 1223  void CaloLong::Fit(Bool_t draw){
1223            exx[numpo] = 0.1;            exx[numpo] = 0.1;
1224            yyy[numpo] = enemip;            yyy[numpo] = enemip;
1225            eyy[numpo] = sqrt(enemip*3.)+sqrt(5.);            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);            //      eyy[numpo] = sqrt(enemip)/(st*0.95);
1233            numpo++;            numpo++;
1234            //      th->Fill(xpos,enemip);            //      th->Fill(xpos,enemip);
# Line 1224  void CaloLong::Fit(Bool_t draw){ Line 1266  void CaloLong::Fit(Bool_t draw){
1266    th = new TH2F(thid,thid,1000,-0.2,xmax,1000,0.,emax*1.2);    th = new TH2F(thid,thid,1000,-0.2,xmax,1000,0.,emax*1.2);
1267    gh = new TGraphErrors(numpo,xxx,yyy,exx,eyy);    gh = new TGraphErrors(numpo,xxx,yyy,exx,eyy);
1268    TString fnam=Form("lfit%s",suf.Data());    TString fnam=Form("lfit%s",suf.Data());
1269    TF1 *lfit  = dynamic_cast<TF1*>(gDirectory->FindObject(fnam));  //  TF1 *lfit  = dynamic_cast<TF1*>(gDirectory->FindObject(fnam));
1270    if ( lfit ) lfit->Delete();    if ( lfit ) lfit->Delete();
1271    lfit = new TF1(fnam,ccurve,0.,xmax,3);    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);    if ( debug ) printf("qtot %f qtotparz %f \n",clp->qtot,qtotparz);
1281    E0 = qtotparz;    E0 = qtotparz;
1282    //  E0 = clp->qtot;    //  E0 = clp->qtot;
# Line 1237  void CaloLong::Fit(Bool_t draw){ Line 1287  void CaloLong::Fit(Bool_t draw){
1287    //  lfit->SetParLimits(0,0.,1000.);    //  lfit->SetParLimits(0,0.,1000.);
1288    //  lfit->SetParLimits(1,-1.,80.);    //  lfit->SetParLimits(1,-1.,80.);
1289    //  lfit->SetParLimits(2,-1.,10.);    //  lfit->SetParLimits(2,-1.,10.);
1290    //  TString optio = "ROW"; // "RO"    // TString optio = "ROW"; // "RO"
1291    TString optio = "RO"; // "RO"    TString optio = "RO"; // "RO"
1292    if ( debug ){    if ( debug ){
1293      optio += "NV";      optio += "NV";
# Line 1247  void CaloLong::Fit(Bool_t draw){ Line 1297  void CaloLong::Fit(Bool_t draw){
1297      if ( draw ) optio += "Q";      if ( draw ) optio += "Q";
1298    };    };
1299    //    //
1300    if ( debug ) printf(" OK, start the fitting procedure...\n");    if ( debug ) printf(" OK, start the fitting procedure... options: %s \n",optio.Data());
1301    //    //
1302    //  fitresult = th->Fit("lfit",optio);    //  fitresult = th->Fit("lfit",optio);
1303      //  lfit->Update();
1304    fitresult = gh->Fit("lfit",optio);    fitresult = gh->Fit("lfit",optio);
1305      //  fitresult = gh->Fit("lfit","ROQW");
1306    //    //
1307    if ( debug ) printf(" the fit is done! result: %i \n",fitresult);    if ( debug ) printf(" the fit is done! result: %i \n",fitresult);
1308    //    //
# Line 1294  void CaloLong::Fit(Bool_t draw){ Line 1346  void CaloLong::Fit(Bool_t draw){
1346      delete w;      delete w;
1347      delete xp;      delete xp;
1348      delete wp;      delete wp;
1349        delete xpp;
1350        delete wpp;
1351      //    Float_t i10max = lfit->Integral(0.,10.*tmax);      //    Float_t i10max = lfit->Integral(0.,10.*tmax);
1352      if ( debug ) printf(" Integral: %f \n",i10max);      if ( debug ) printf(" Integral: %f \n",i10max);
1353      //      //
# Line 1344  void CaloLong::Fit(Bool_t draw){ Line 1398  void CaloLong::Fit(Bool_t draw){
1398      if ( th ) th->Delete();      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;    //  delete lfit;
1411    //    //
1412  };  };

Legend:
Removed from v.1.13  
changed lines
  Added in v.1.19

  ViewVC Help
Powered by ViewVC 1.1.23