/[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.12 by mocchiut, Mon Aug 17 15:00:36 2009 UTC revision 1.18 by mocchiut, Thu May 13 13:55:37 2010 UTC
# Line 281  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    
# Line 297  Calo2D::Calo2D(PamLevel2 *l2p){   Line 299  Calo2D::Calo2D(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 375  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 383  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 426  void Calo2D::Draw(Int_t plane){ Line 430  void Calo2D::Draw(Int_t plane){
430    gStyle->SetNdivisions(510,"XY");    gStyle->SetNdivisions(510,"XY");
431    //    //
432    for (Int_t p=minp; p<maxp;p++){    for (Int_t p=minp; p<maxp;p++){
433      TString hid = Form("c2dp%i",p);          TString hid = Form("c2dp%i%s",p,suf.Data());        
434      TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));      TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
435      if ( tc ){      if ( tc ){
436        //         tc->Clear();        //         tc->Clear();
# Line 434  void Calo2D::Draw(Int_t plane){ Line 438  void Calo2D::Draw(Int_t plane){
438        tc = new TCanvas(hid,hid);        tc = new TCanvas(hid,hid);
439      };      };
440      //      //
441      TString thid = Form("h2dp%i",p);          TString thid = Form("h2dp%i%s",p,suf.Data());      
442      TH2F *th  = dynamic_cast<TH2F*>(gDirectory->FindObject(thid));      TH2F *th  = dynamic_cast<TH2F*>(gDirectory->FindObject(thid));
443      if ( th ) th->Delete();      if ( th ) th->Delete();
444      //   th->Clear();      //   th->Clear();
# Line 525  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 572  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 661  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 671  CaloLong::CaloLong(PamLevel2 *l2p){   Line 682  CaloLong::CaloLong(PamLevel2 *l2p){  
682    umax = 100.;    umax = 100.;
683    slmax = "";    slmax = "";
684    sumax = "";    sumax = "";
685      suf = "";
686    xyaverage = true;    xyaverage = true;
687    //    //
688      letmax = -1;
689      heavytail = false;
690      //  lmipth = 100.;
691      lmipth = 0.;
692      //
693  };  };
694    
695    TF1 *CaloLong::GetFit(){
696      TString fnam=Form("lfit%s",suf.Data());
697      TF1 *lfit  = dynamic_cast<TF1*>(gDirectory->FindObject(fnam));
698      if ( lfit ) return lfit;
699      return NULL;
700    }
701    
702  void CaloLong::MaskSection(TString sec){  void CaloLong::MaskSection(TString sec){
703    sec.ToUpper();    sec.ToUpper();
704    if ( sec.Contains("XO") ) maskXO = true;    if ( sec.Contains("XO") ) maskXO = true;
# Line 697  void CaloLong::UnMaskSection(TString sec Line 721  void CaloLong::UnMaskSection(TString sec
721    
722  void CaloLong::Clear(){  void CaloLong::Clear(){
723    //    //
724      if ( debug ) printf(" Clear called \n");
725      //
726    memset(eplane,0, 2*22*sizeof(Float_t));    memset(eplane,0, 2*22*sizeof(Float_t));
727    //    //
728    chi2 = 0.;    chi2 = 0.;
# Line 839  void CaloLong::Process(){ Line 865  void CaloLong::Process(){
865    Bool_t gof = true;    Bool_t gof = true;
866    for (Int_t i=0; i < L2->GetCaloLevel1()->istrip; i++){    for (Int_t i=0; i < L2->GetCaloLevel1()->istrip; i++){
867      mip = L2->GetCaloLevel1()->DecodeEstrip(i,view,plane,strip);      mip = L2->GetCaloLevel1()->DecodeEstrip(i,view,plane,strip);
868        //
869        if ( !usepl18x && view==0 && plane==18 ) mip = 0.;
870        //
871      gof = true;      gof = true;
872      if ( maskXE && (plane%2)==0 && view==1 ) gof = false;      if ( maskXE && (plane%2)==0 && view==1 ) gof = false;
873      if ( maskXO && (plane%2)!=0 && view==1 ) gof = false;      if ( maskXO && (plane%2)!=0 && view==1 ) gof = false;
# Line 849  void CaloLong::Process(){ Line 878  void CaloLong::Process(){
878    //    //
879    // inclination factor (stolen from Daniele's code)    // inclination factor (stolen from Daniele's code)
880    //    //
881    Float_t ytgx = 0;    Float_t ytgx = 0.;
882    Float_t ytgy = 0;    Float_t ytgy = 0.;
883    ytgx = 0.76 * clp->tanx[0];    ytgx = 0.76 * clp->tanx[0];
884    ytgy = 0.76 * clp->tany[0];      ytgy = 0.76 * clp->tany[0];  
885    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 944  void CaloLong::SetEnergies(Float_t myene Line 973  void CaloLong::SetEnergies(Float_t myene
973      };      };
974    };    };
975    //    //
976      if ( !clp && debug ) printf(" ERROR!! no CaloLevel2 obj!!\n");
977      //
978    // inclination factor (stolen from Daniele's code)    // inclination factor (stolen from Daniele's code)
979    //    //
980    Float_t ytgx = 0;    Float_t ytgx = 0.;
981    Float_t ytgy = 0;    Float_t ytgy = 0.;
982    ytgx = 0.76 * clp->tanx[0];    ytgx = 0.76 * clp->tanx[0];
983    ytgy = 0.76 * clp->tany[0];      ytgy = 0.76 * clp->tany[0];  
984    X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) );    X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) );
985      if ( debug ) printf(" X0pl %f tanx %f tany %f \n",X0pl,ytgx,ytgy);
986    //    //
987    // Find experimental plane of maximum    // Find experimental plane of maximum
988    //    //
# Line 971  void CaloLong::SetEnergies(Float_t myene Line 1003  void CaloLong::SetEnergies(Float_t myene
1003    //    //
1004    if ( vmax == 0 ) pmax++;    if ( vmax == 0 ) pmax++;
1005    etmax = pmax * X0pl;    etmax = pmax * X0pl;
1006      letmax = etmax;
1007    //    //
1008    if ( debug ) this->Print();    if ( debug ) this->Print();
1009    if ( debug ) printf(" exit \n");    if ( debug ) printf(" exit \n");
# Line 1066  void CaloLong::Fit(Bool_t draw){ Line 1099  void CaloLong::Fit(Bool_t draw){
1099      gStyle->SetNdivisions(510,"XY");      gStyle->SetNdivisions(510,"XY");
1100    };    };
1101    //    //
1102    TString hid = Form("clongfit");          TString hid = Form("clongfit%s",suf.Data());  
1103    TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));    TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
1104    //  if ( tc ) tc->Delete();    //  if ( tc ) tc->Delete();
1105    //  if ( tc ) tc->Close();    //  if ( tc ) tc->Close();
# Line 1076  void CaloLong::Fit(Bool_t draw){ Line 1109  void CaloLong::Fit(Bool_t draw){
1109      if ( tc ) tc->cd();      if ( tc ) tc->cd();
1110    };    };
1111    //    //
1112    TString thid = Form("hlongfit");          TString thid = Form("hlongfit%s",suf.Data());
1113    TH2F *th  = dynamic_cast<TH2F*>(gDirectory->FindObject(thid));    TH2F *th  = dynamic_cast<TH2F*>(gDirectory->FindObject(thid));
1114    if ( th ) th->Delete();    if ( th ) th->Delete();
1115    //    //
1116    TString ghid = Form("glongfit");          TString ghid = Form("glongfit%s",suf.Data());
1117    TGraphErrors *gh  = dynamic_cast<TGraphErrors*>(gDirectory->FindObject(ghid));    TGraphErrors *gh  = dynamic_cast<TGraphErrors*>(gDirectory->FindObject(ghid));
1118    if ( gh ) gh->Delete();    if ( gh ) gh->Delete();
1119    //    //
# Line 1152  void CaloLong::Fit(Bool_t draw){ Line 1185  void CaloLong::Fit(Bool_t draw){
1185          exx[numpo] = 0.1;          exx[numpo] = 0.1;
1186          yyy[numpo] = enemip;          yyy[numpo] = enemip;
1187          eyy[numpo] = sqrt(enemip*3.)+sqrt(5.);          eyy[numpo] = sqrt(enemip*3.)+sqrt(5.);
1188            //if ( xpos > letmax && enemip > lmipth && heavytail) eyy[numpo] = (sqrt(enemip*3.)+sqrt(5.))/numpo;
1189            if ( xpos > letmax && enemip > lmipth && heavytail) eyy[numpo] = sqrt(enemip)/5.;
1190            if ( xpos > letmax-1 && xpos < letmax+1 && heavytail ) eyy[numpo] /= 5.;
1191            //if ( xpos > letmax-2 && xpos < letmax+1 && heavytail ) eyy[numpo] /= 7.;
1192            if ( xpos < 3. && heavytail ) eyy[numpo] /= 5.;
1193          numpo++;          numpo++;
1194          //      th->Fill(xpos,enemip);          //      th->Fill(xpos,enemip);
1195          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);
1196        };        };
1197      } else {      } else {
1198        for (Int_t jj=0; jj<2; jj++){        for (Int_t jj=0; jj<2; jj++){
# Line 1184  void CaloLong::Fit(Bool_t draw){ Line 1222  void CaloLong::Fit(Bool_t draw){
1222            exx[numpo] = 0.1;            exx[numpo] = 0.1;
1223            yyy[numpo] = enemip;            yyy[numpo] = enemip;
1224            eyy[numpo] = sqrt(enemip*3.)+sqrt(5.);            eyy[numpo] = sqrt(enemip*3.)+sqrt(5.);
1225              //if ( xpos > letmax && enemip > lmipth && heavytail) eyy[numpo] = (sqrt(enemip*3.)+sqrt(5.))/numpo;
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)/5.;
1228              if ( xpos > letmax-1 && xpos < letmax+1 && heavytail ) eyy[numpo] /= 5.;
1229            //if ( xpos > letmax-2 && xpos < letmax+1 && heavytail ) eyy[numpo] /= 7.;
1230              if ( xpos < 3. && heavytail ) eyy[numpo] /= 5.;        
1231            //      eyy[numpo] = sqrt(enemip)/(st*0.95);            //      eyy[numpo] = sqrt(enemip)/(st*0.95);
1232            numpo++;            numpo++;
1233            //      th->Fill(xpos,enemip);            //      th->Fill(xpos,enemip);
# Line 1220  void CaloLong::Fit(Bool_t draw){ Line 1264  void CaloLong::Fit(Bool_t draw){
1264    //  th = new TH2F(thid,thid,int(NC*1.5),-0.2,xmax);    //  th = new TH2F(thid,thid,int(NC*1.5),-0.2,xmax);
1265    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);
1266    gh = new TGraphErrors(numpo,xxx,yyy,exx,eyy);    gh = new TGraphErrors(numpo,xxx,yyy,exx,eyy);
1267    TF1 *lfit  = dynamic_cast<TF1*>(gDirectory->FindObject("lfit"));    TString fnam=Form("lfit%s",suf.Data());
1268      TF1 *lfit  = dynamic_cast<TF1*>(gDirectory->FindObject(fnam));
1269    if ( lfit ) lfit->Delete();    if ( lfit ) lfit->Delete();
1270    lfit = new TF1("lfit",ccurve,0.,xmax,3);    lfit = new TF1(fnam,ccurve,0.,xmax,3);
1271    //   if ( !lfit ){
1272    //     lfit = new TF1(fnam,ccurve,0.,xmax,3);
1273    //   };
1274    //   lfit->Clear();
1275    //   lfit->SetName(fnam);
1276    //   lfit->SetRange(0.,xmax);
1277    //  lfit->Print();
1278      //
1279    if ( debug ) printf("qtot %f qtotparz %f \n",clp->qtot,qtotparz);    if ( debug ) printf("qtot %f qtotparz %f \n",clp->qtot,qtotparz);
1280    E0 = qtotparz;    E0 = qtotparz;
1281    //  E0 = clp->qtot;    //  E0 = clp->qtot;
# Line 1233  void CaloLong::Fit(Bool_t draw){ Line 1286  void CaloLong::Fit(Bool_t draw){
1286    //  lfit->SetParLimits(0,0.,1000.);    //  lfit->SetParLimits(0,0.,1000.);
1287    //  lfit->SetParLimits(1,-1.,80.);    //  lfit->SetParLimits(1,-1.,80.);
1288    //  lfit->SetParLimits(2,-1.,10.);    //  lfit->SetParLimits(2,-1.,10.);
1289    //  TString optio = "ROW"; // "RO"    // TString optio = "ROW"; // "RO"
1290    TString optio = "RO"; // "RO"    TString optio = "RO"; // "RO"
1291    if ( debug ){    if ( debug ){
1292      optio += "NV";      optio += "NV";
# Line 1243  void CaloLong::Fit(Bool_t draw){ Line 1296  void CaloLong::Fit(Bool_t draw){
1296      if ( draw ) optio += "Q";      if ( draw ) optio += "Q";
1297    };    };
1298    //    //
1299    if ( debug ) printf(" OK, start the fitting procedure...\n");    if ( debug ) printf(" OK, start the fitting procedure... options: %s \n",optio.Data());
1300    //    //
1301    //  fitresult = th->Fit("lfit",optio);    //  fitresult = th->Fit("lfit",optio);
1302      //  lfit->Update();
1303    fitresult = gh->Fit("lfit",optio);    fitresult = gh->Fit("lfit",optio);
1304      //  fitresult = gh->Fit("lfit","ROQW");
1305    //    //
1306    if ( debug ) printf(" the fit is done! result: %i \n",fitresult);    if ( debug ) printf(" the fit is done! result: %i \n",fitresult);
1307    //    //
# Line 1290  void CaloLong::Fit(Bool_t draw){ Line 1345  void CaloLong::Fit(Bool_t draw){
1345      delete w;      delete w;
1346      delete xp;      delete xp;
1347      delete wp;      delete wp;
1348        delete xpp;
1349        delete wpp;
1350      //    Float_t i10max = lfit->Integral(0.,10.*tmax);      //    Float_t i10max = lfit->Integral(0.,10.*tmax);
1351      if ( debug ) printf(" Integral: %f \n",i10max);      if ( debug ) printf(" Integral: %f \n",i10max);
1352      //      //
# Line 1340  void CaloLong::Fit(Bool_t draw){ Line 1397  void CaloLong::Fit(Bool_t draw){
1397      if ( th ) th->Delete();      if ( th ) th->Delete();
1398    };    };
1399    //    //
1400    //   gMinuit->SetPrintLevel(1);
1401    //   gMinuit->mnamin();
1402    //   //  gMinuit->mnfree(0);
1403    //   gMinuit->mncler();
1404    //   gMinuit->mnrset(1);
1405      //  gMinuit->Delete();
1406      //ROOT::Math::Minimizer::Clear();
1407      //  gMinuit->mnhelp("*");
1408      //
1409    //  delete lfit;    //  delete lfit;
1410    //    //
1411  };  };
# Line 1369  void CaloLong::Draw(Int_t view){ Line 1435  void CaloLong::Draw(Int_t view){
1435    //    //
1436    if ( maxv != -1 ){    if ( maxv != -1 ){
1437      for (Int_t v=minv; v<maxv;v++){      for (Int_t v=minv; v<maxv;v++){
1438        TString hid = Form("clongv%i",v);        TString hid = Form("clongv%i%s",v,suf.Data());    
1439        TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));        TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
1440        if ( tc ){        if ( tc ){
1441          //       tc->Clear();          //       tc->Clear();
# Line 1377  void CaloLong::Draw(Int_t view){ Line 1443  void CaloLong::Draw(Int_t view){
1443          tc = new TCanvas(hid,hid);          tc = new TCanvas(hid,hid);
1444        };        };
1445        //        //
1446        TString thid = Form("hlongv%i",v);                TString thid = Form("hlongv%i%s",v,suf.Data());  
1447        TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));        TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));
1448        if ( th ) th->Delete();        if ( th ) th->Delete();
1449        //         th->Clear();        //         th->Clear();
# Line 1396  void CaloLong::Draw(Int_t view){ Line 1462  void CaloLong::Draw(Int_t view){
1462      };      };
1463    } else {    } else {
1464      //      //
1465      TString hid = Form("clongvyvx");          TString hid = Form("clongvyvx%s",suf.Data());      
1466      TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));      TCanvas *tc  = dynamic_cast<TCanvas*>(gDirectory->FindObject(hid));
1467      if ( tc ){      if ( tc ){
1468      } else {      } else {
1469        tc = new TCanvas(hid,hid);        tc = new TCanvas(hid,hid);
1470      };      };
1471      //      //
1472      TString thid = Form("hlongvyvx");        TString thid = Form("hlongvyvx%s",suf.Data());      
1473      TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));      TH1F *th  = dynamic_cast<TH1F*>(gDirectory->FindObject(thid));
1474      if ( th ) th->Delete();      if ( th ) th->Delete();
1475      th = new TH1F(thid,thid,44,-0.5,43.5);      th = new TH1F(thid,thid,44,-0.5,43.5);
Legend:
Removed from v.1.12  
changed lines
  Added in v.1.18
  ViewVC Help
Powered by ViewVC 1.1.23