C++ TString::Data方法代码示例

void plotxsec(TString xsectype/*= "vm"*/, bool ploty/*=kFALSE*/,bool sim/*=kFALSE*/,bool comp/*=kFALSE*/){
  if (setup(xsectype)==kFALSE) return;

  gStyle->SetOptStat(0);

  TString methodB = TString::Format("Yield - Method B {%s}",_topNames[4].Data());
  TString methodA = TString::Format("Yield - Method A");

  TH1F* hYW_Dir[5];
  TH1F* hYW_5D[5];
  TH1F* hYW_TH[5];
  TH1F* hYWnorm_Dir[5];
  TH1F* hYWnorm_5D[5];
  TLegend* l_Dir = new TLegend(0.1,0.7,0.5,0.9);
  TLegend* l_5D  = new TLegend(0.1,0.7,0.5,0.9);
  TString hname_hYW_Dir = TString::Format("hYW_Dir/Varset1/hYW_ACC_CORR");
  TString hname_hYW_5D  = TString::Format("hYW/Varset1/hYW");
  TString hname_hYW_TH  = TString::Format("hYW_Dir/Varset1/hYW_TH");

  TCanvas* c_Dir = new TCanvas(TString::Format("Dir-%s:%s", xsectype.Data(), _q2w_bng.Data()), 
                               TString::Format("Dir-%s:%s", xsectype.Data(), _q2w_bng.Data()) );
  c_Dir->SetGrid(); 
  TCanvas* c_5D = new TCanvas(TString::Format("5D-%s:%s", xsectype.Data(), _q2w_bng.Data()), 
                              TString::Format("5D-%s:%s", xsectype.Data(), _q2w_bng.Data()) );
  c_5D->SetGrid();
  //! For each Top: hYW_Dir,hYW_5D --> hYWnorm_Dir,hYWnorm_5D & then draw
  TFile* fy[5];
  if (sim) memcpy(fy,_fysim,sizeof(_fysim));
  else     memcpy(fy,_fyexp,sizeof(_fyexp));
  for (int iTop=0;iTop<5;iTop++){
    if (fy[iTop]==NULL) continue;
    cout << "Top = " << iTop << endl;
    hYW_Dir[iTop] = (TH1F*)fy[iTop]->Get(hname_hYW_Dir);
    hYW_5D[iTop]  = (TH1F*)fy[iTop]->Get(hname_hYW_5D);
    if(sim) hYW_TH[iTop]  = (TH1F*)fy[iTop]->Get(hname_hYW_TH);
    if (hYW_Dir[iTop]==NULL || hYW_5D[iTop]==NULL) continue;
    hYW_Dir[iTop]->SetMarkerStyle(20+iTop);
    hYW_Dir[iTop]->SetMarkerColor(2+iTop);
    //hYW_Dir[iTop]->SetMaximum(150000);
    hYW_5D[iTop]->SetMarkerStyle(20+iTop); //+5
    hYW_5D[iTop]->SetMarkerColor(2+iTop);  //+5
    //hYW_5D[iTop]->SetMaximum(500000);
    /*if (sim) hYW_TH[iTop]->SetMarkerColor(kFullStar);  
    if (sim) hYW_TH[iTop]->SetMarkerColor(kBlack);*/
    if (sim) hYW_TH[iTop]->SetFillColor(kRed);  
    if (sim) hYW_TH[iTop]->SetFillStyle(3001);
    //Normalize yields to LUM*vgflux*binw_w*binw_q2
    hYWnorm_Dir[iTop] = normalizeYield(hYW_Dir[iTop]);
    hYWnorm_Dir[iTop]->SetMaximum(20);
    hYWnorm_5D[iTop]  = normalizeYield(hYW_5D[iTop]);
    hYWnorm_5D[iTop]->SetMaximum(20);

   if (iTop==0){
     c_5D->cd();
     l_5D->SetHeader( TString::Format("%s - 5D Method", XSECTITLE.Data()) );
     l_5D->AddEntry( hYWnorm_5D[iTop], TString::Format("%s", _topNames[iTop].Data()) );
     hYWnorm_5D[iTop]->SetMinimum(0);
     ploty?hYW_5D[iTop]->Draw("p"):hYWnorm_5D[iTop]->Draw("p");
     if (sim)
     {
      hYW_TH[iTop]->Draw("sames");
      l_5D->AddEntry( hYW_TH[iTop], TString::Format("Thrown") );
     } 

     c_Dir->cd();
     l_Dir->SetHeader( TString::Format("%s - Direct Method", XSECTITLE.Data()) );
     l_Dir->AddEntry( hYWnorm_Dir[iTop], TString::Format("%s", _topNames[iTop].Data()) );
     hYWnorm_Dir[iTop]->SetMinimum(0);
     ploty?hYW_Dir[iTop]->Draw("p"):hYWnorm_Dir[iTop]->Draw("p");
     if (sim)
     {
      hYW_TH[iTop]->Draw("bar sames");
      l_Dir->AddEntry( hYW_TH[iTop], TString::Format("Thrown") );
     }    
   }else{
     c_5D->cd();
     l_5D->AddEntry( hYWnorm_5D[iTop], TString::Format("%s", _topNames[iTop].Data()) );
     ploty?hYW_5D[iTop]->Draw("p sames"):hYWnorm_5D[iTop]->Draw("p sames");
     if (sim)
     {
      hYW_TH[iTop]->Draw("bar sames");
     }

     c_Dir->cd();
     l_Dir->AddEntry( hYWnorm_Dir[iTop], TString::Format("%s", _topNames[iTop].Data()) );
     ploty?hYW_Dir[iTop]->Draw("p sames"):hYWnorm_Dir[iTop]->Draw("p sames");
     if (sim)
     {
      hYW_TH[iTop]->Draw("bar sames");
     }
   }
   c_5D->cd();
   l_5D->Draw("same");
   c_Dir->cd();
   l_Dir->Draw("same"); 
  }
  //new - comp 5D & Dir for Top2 & 5
  if (comp){
    TCanvas* c_comp = new TCanvas(TString::Format("5D_Dir_Cpmp-%s", xsectype.Data()), TString::Format("5D & Dir Comp %s", xsectype.Data()));
    TLegend* l_comp = new TLegend(0.1,0.7,0.5,0.9);
//.........这里部分代码省略.........

Int_t EventCount(TString inputlist, TString outfile, Int_t nev=-1)
{  

  HLoop* loop = new HLoop   (kTRUE    );  // kTRUE : create Hades  (needed to work with standard eventstructure)
         loop ->addMultFiles(inputlist);

//-----Inicialize-classes-we-need-to-use-----------------
  PidParticle p;
  dEdx        d;
  HMultCorr   cCorr;
//-------------------------------------------------------

//--------------------------------------------------------------------------------------------------------
#include "/u/parfenov/anaflow/EventCountVar.h"     //include all constants and variables we need          
#include "/u/parfenov/anaflow/outsigma.C"          //Declaration of the 3-sigma cuts for our pt & rapidity
TH2F*  hProtPtVsY_Eff;
TFile* foutfile       = new TFile(outfile.Data(), "recreate"); 
TFile* FileProtEff    = new TFile("EffCorrForProtPtVsY_WithWideBetaPcuts_ShieldGosia.root", "read");
       hProtPtVsY_Eff = (TH2F*) FileProtEff->Get("hProtPtVsY_EffSm");
       foutfile->cd();

#include "/u/parfenov/anaflow/EventCountHisto.h"   //include all histograms we need                       

FFlow pfy0pt[NC][NRAPI][NPTRN]; //--flow--9-bins-Yo--12-bins-in-Pt-in-3-centrality-sets--
  for (Int_t i=0; i<NC;i++){
   for(Int_t j=0; j<NRAPI; j++){
     for(Int_t k=0; k<NPTRN; k++){
       //-proton------------------------------------------
       pfy0pt[i][j][k].THDeclare("pCent",i,"Yo",j,"Pt",k);
     }
   }
  }
//--------------------------------------------------------------------------------------------------------

  if(!loop->setInput("-*,+HParticleCand,+HParticleEvtInfo,+HStart2Hit,+HStart2Cal,+HWallHit,+HRpcCluster,+HTofHit,+HWallEventPlane")) { exit(1); }
  loop->printCategories();
  
  TIterator* iterWallHit = 0;
  if (loop->getCategory("HWallHit")  ) iterWallHit   = loop->getCategory("HWallHit")->MakeIterator();
  
  HParticleEvtInfo* evtinfo;
  
  HParticleCand* pParticleCand; HCategory* candCat = (HCategory*)HCategoryManager::getCategory(catParticleCand); if(!candCat){exit(1);}
  HWallHit*      pWallHit;      HCategory* wallCat = (HCategory*)HCategoryManager::getCategory(catWallHit);      if(!wallCat){exit(1);}

  HEnergyLossCorrPar *momCorPar = new HEnergyLossCorrPar();
  momCorPar->setDefaultPar("apr12");
      
  /////////////////////////////////////////////////////////
  //-----Loop-over-events----------------------------------
  /////////////////////////////////////////////////////////
  
  loop->addMultFiles(inputlist);

  Int_t events = loop->getEntries();
  if(nev < events && nev >= 0 ) events = nev;
  
  for (Int_t i=1; i<events;i++){
    
    if (loop->nextEvent(i)<=0){
      cout << "End recieved with IF exit" << endl;
      break;
    }
    if (i%1000 == 0){
      cout << "    event " << i << endl;
    }
    
    if( loop->isNewFile(currFName) ){ 
        currFName = currFName( currFName.Last('/')+1,currFName.Length()-currFName.Last('/')-1 );
        currBeName = currFName(0, 13); //-simply-the-file-name-w/o-EBchunk--
        currFDay = currFName( 4, 3 ); //-we-cut-out-day-number-position-starting-from-4th-digit-with-length-of-3digits--
        DAY_NUM = currFDay.Atoi();
        currTime = currFName( 7, 2 ); //-we-cut-out-day-hour--position-starting-from-4th-digit-with-length-of-3digits--
        HR=currTime.Atoi();
        currTime = currFName( 9, 2 ); //-we-cut-out-day-minute-position-starting-from-4th-digit-with-length-of-3digits--
        MN=currTime.Atoi();
        currTime = currFName(11, 2 ); //-we-cut-out-day-second-position-starting-from-4th-digit-with-length-of-3digits--
        SC=currTime.Atoi();
        AbsMinute = MN + HR*60 + DAY_NUM*24*60;

        #include "/u/parfenov/anaflow/anaKaons_params_gen7_Auto_RC_RW_updated_gen108RminEqZero_11MultRecent.cc" //-for recentering---

        //--Now-we-read-multiplicity-parameters-for-this-beXXXXXXXXXX-file--------------------//
        mulVal = cCorr.getLineValuesAsVectFromCalibFile( currBeName );
        cout<<"mulVal "<< mulVal[0] <<" "<< mulVal[1] <<" "<< mulVal[2] <<" "<< mulVal[3] <<endl;
        fAverageMTS[0]=mulVal[4];  //multiplicity in tracking sector 1
        fAverageMTS[1]=mulVal[5];  //multiplicity in tracking sector 2
        fAverageMTS[2]=mulVal[6];  //...
        fAverageMTS[3]=mulVal[7];
        fAverageMTS[4]=mulVal[8];
        fAverageMTS[5]=mulVal[9];  //multiplicity in tracking sector 6
        printf("fAverageMTS[0,1,2,3,4,5]=[%6.3f,%6.3f,%6.3f,%6.3f,%6.3f,%6.3f,]\n", fAverageMTS[0], fAverageMTS[1], fAverageMTS[2], fAverageMTS[3], fAverageMTS[4], fAverageMTS[5]);
        
        //fTrkMultDay108
        if(mulVal[2]>0.0){
          //-if-correction-is-found-or-if-it-was-not-zero-defined-(in-case-of-not-working-MDCs)---
          fTrkMultScaler = fTrkMultDay108/mulVal[2];
          printf("==>Correction multiplier nTrkMult for this file is: (%f/%f)=%f\n", fTrkMultDay108, mulVal[2], fTrkMultScaler );
        }else{
          //--if-correction-entry-is-not-found-(may-be-not-created-once)------------------
//.........这里部分代码省略.........

void plot(TChain *tree_data, TChain *tree_mc, const TString &varRD, const TString &varMC, const TString& hName, const TString& xtitle, const TString& ytitle, int nbin, float xlow, float xmax, TCut cutRD, TCut cutMC, TCut cutjet, string name, bool print){

  TCanvas* c = new TCanvas(Form("c_%s_%s",hName.Data(),name),Form("c_%s_%s",hName.Data(),name),1); 
  TLegend *l;

  if(hName.Contains("deposit") || hName.Contains("pT")){ 
    gPad->SetLogy();
    l = new TLegend(0.68,0.78,0.90,0.90);
  }else{
    l = new TLegend(0.18,0.67,0.42,0.82);
  }
  //gPad->SetGrid(1,1);
  
  TH1 *h_data = new TH1D(Form("h_data_%s_%s",hName.Data(),name),"h_data",nbin,xlow,xmax);
  TH1 *h_mc = new TH1D(Form("h_mc_%s_%s",hName.Data(),name),"h_mc",nbin,xlow,xmax);
  h_data->Sumw2();

  tree_mc->Draw(Form("%s>>h_mc_%s_%s",varMC.Data(),hName.Data(),name),cutMC,"");
  tree_data->Draw(Form("%s>>h_data_%s_%s",varRD.Data(),hName.Data(),name),cutRD && cutjet,"Esame");

  double num_data = tree_data->GetEntries();
  double num_mc = tree_mc->GetEntries();
  double scale = num_data/num_mc;

  double pivotal_mc = 0;
  double pivotal_data = 0;
  TH1 *h_mc_eta= new TH1F("h_mc_eta","h_mc_eta",60,-3,3);
  TH1 *h_data_eta= new TH1F("h_data_eta","h_data_eta",60,-3,3);
  tree_mc->Project("h_mc_eta","muon_eta","muon_eta < 2.1 && muon_eta > -2.1");
  tree_data->Project("h_data_eta","random_eta",cutjet);
  pivotal_mc = h_mc_eta->GetEntries();
  pivotal_data = h_data_eta->GetEntries();
  h_mc_eta->Delete();
  h_data_eta->Delete();
  cout << "Data= " << num_data << endl;
  cout << "MC= " << num_mc << endl;
  cout << "CH pivotal data= " << pivotal_data << endl;
  cout << "CH pivotal mc= " << pivotal_mc << endl;

  TLatex *label= new TLatex;
  label->SetNDC();
  label->SetTextSize(0.05);
  label->DrawLatex(0.18,0.87,"CMS Preliminary 2010");

  h_data->SetMarkerStyle(kFullCircle);
  h_data->SetMarkerSize(0.8);//1
  h_data->SetMarkerStyle(8);
  h_data->SetStats(0);
  h_data->GetXaxis()->SetTitle(xtitle.Data());
  h_data->GetYaxis()->SetTitle(ytitle.Data());

  h_mc->SetLineColor(2);
  h_mc->SetLineWidth(2);
  h_mc->SetFillStyle(1001);
  h_mc->SetFillColor(kRed-8);
  h_mc->SetStats(0);
  h_mc->GetXaxis()->SetTitle(xtitle.Data());
  h_mc->GetYaxis()->SetTitle(ytitle.Data());

  h_data->Scale(1/pivotal_data/h_data->GetBinWidth(1));
  h_mc->Scale(1/pivotal_mc/h_mc->GetBinWidth(1));

  h_mc->SetMaximum(4); // 2.5

  //l->AddEntry(h_data,"random cone for MinBias","p");
  l->AddEntry(h_data,"random cone for W #rightarrow #mu#nu","p");
  l->AddEntry(h_mc,"muon for W #rightarrow #mu#nu","F");
  l->SetTextSize(0.04);
  l->SetFillColor(0);
  l->SetLineColor(0);
  l->Draw();
  
  if(print){ 
    c->Print(Form("c_%s_%s_Jun1.pdf",hName.Data(),name));
    c->Print(Form("c_%s_%s_Jun1.eps",hName.Data(),name));
    c->Print(Form("c_%s_%s_Jun1.png",hName.Data(),name));
  }
}

void AnalysisSparse(Bool_t save_output = kFALSE)
{
  gStyle->SetGridColor(kGray);
  //  TString tmpstr(fname);
  //  if (tmpstr.Contains("data")) {
  //    Printf("!!! Real Data !!!");
  //    mc = kFALSE;
  //  }
  TString gtitle = Form("Monte Carlo, %s", graph_name.Data());
  grapht = graph_name.Data();
  Double_t grx[999], gry[999], gry2[999], gry3[999], gry4[999],
    gry_eff[999], gry_fix[999], grxE[999];
  Double_t gry22[999], gry22E[999], grx22E[999];
  Double_t gry_true[999], gry_true_eff[999], gry_true_effE[999];
  TH1::AddDirectory(kFALSE);
  TFile::SetCacheFileDir(gSystem->HomeDirectory());
  TFile *f = TFile::Open(fname.Data(), "CACHEREAD");
  if (!f) return;
  TList *l; f->GetObject(lname.Data(), l);
  if (!l) return;
  Int_t bf[999], bl[999];
  Int_t nn = FindExactRange(((THnSparse *)(l->FindObject(s1name.Data())))->
                            Projection(1), del_step, bf, bl);
  //  Int_t nn = FindRange5(bf, bl);
  Bool_t binhaluska = kFALSE;
  if (binAnders) {
    nn = 8;
    bf[0] = 6;bf[1] =  9;bf[2] = 11;bf[3] = 16;bf[4] = 21;bf[5] = 26;
    bl[0] = 8;bl[1] = 10;bl[2] = 15;bl[3] = 20;bl[4] = 25;bl[5] = 30;

    bf[6] = 31;bf[7] = 41;
    bl[6] = 40;bl[7] = 50;
  }
  Printf("number of intervals = %d =>", nn);

  Int_t count = 0;
  Double_t ptmean = 0, value = 0;
  Int_t fitStatus = -1;
  gStyle->SetOptStat(0);
  TCanvas *c = new TCanvas("c", "Signal & Background");
  c->Divide(5, 5); c->Modified(); c->Draw();
  TCanvas *c2 = (TCanvas *)c->DrawClone("c2");
  c2->SetTitle("Phi mesons (raw)"); c2->Modified(); c2->Draw();
  TCanvas *c3, *c4;
  if (mc) {
    c3 = (TCanvas *)c->DrawClone("c3");
    c3->SetTitle("Phi mesons (gen)"); c3->Modified(); c3->Draw();
    c4 = (TCanvas *)c->DrawClone("c4");
    c4->SetTitle("Phi mesons (true)"); c4->Modified(); c4->Draw();
  }

  for (Int_t i = 0; i < nn; i++) {
    c->cd(count + 1)->SetGrid();
    h1 = (TH1D *)PullHisto(l, s1name.Data(), bf[i], bl[i], ptmean);
    h1->SetLineColor(kRed);
    h1->GetXaxis()->SetTitle("inv. mass, GeV/c^2");
    h1->Draw("hist");

    h3_p = (TH1D *)PullHisto(l, s3name_p.Data(), bf[i], bl[i], ptmean);
    h3_m = (TH1D *)PullHisto(l, s3name_m.Data(), bf[i], bl[i], ptmean);
    // !!!!!!!!!!!!!!!!!!!!!!!!
    if (count==0) h3_p = h1;
    // !!!!!!!!!!!!!!!!!!!!!!!!
    else {
      h3_p->Add(h3_m);
      //      h3_p->Add((TH1D *)PullHisto(l, smix.Data(), bf[i], bl[i], ptmean));
      //      h3_p->Add((TH1D *)PullHisto(l, smixpp.Data(), bf[i], bl[i], ptmean));
      //      h3_p->Add((TH1D *)PullHisto(l, smixmm.Data(), bf[i], bl[i], ptmean));
      Norm(h1, h3_p, norm[0], norm[1]);
    }
    h3_p->SetLineColor(kBlue);
    h3_p->Draw("hist, same");

    if (mc) {
      c3->cd(count + 1)->SetGrid();
      Printf("%s", s1namegen.Data());
      hg = (TH1D *)PullHisto(l, s1namegen.Data(), bf[i], bl[i], ptmean);
      hg->SetLineColor(kMagenta);
      hg->GetXaxis()->SetTitle("inv. mass, GeV/c^2");
      hg->Draw("hist");
      c4->cd(count + 1)->SetGrid();
      ht = (TH1D *)PullHisto(l, s1nametrue.Data(), bf[i], bl[i], ptmean);
      ht->SetLineColor(kMagenta-5);
      ht->GetXaxis()->SetTitle("inv. mass, GeV/c^2");
      ht->Draw("hist");
    }
    c2->cd(count + 1)->SetGrid();
    TH1 *hh = (TH1 *)h1->Clone("hh");
    hh->SetLineColor(kRed+1);
    hh->Add(h3_p, -1);
    /// !!!!!!!!!!!!!!!!!!!!!!
    //////////    if ((ilist == 3) && (count < 2)) hh->Reset();
    // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    hh->Draw("hist");

    // !!!!!!!!!!!!!!!!!!
    ff->SetParameters(0.1, 1.02, 0.004, -25000., 0., 0., 0.);
    ff->SetLineColor(hh->GetLineColor());
    ff->SetLineWidth(1);
    //    ff->SetLineStyle(kDashed);
//.........这里部分代码省略.........

void TestSPlot()
{
//This tutorial illustrates the use of class TSPlot and of the sPlots method
//
//It is an example of analysis of charmless B decays, performed for BABAR.
//One is dealing with a data sample in which two species are present:
//the first is termed signal and the second background. 
//A maximum Likelihood fit is performed to obtain the two yields N1 and N2
//The fit relies on two discriminating variables collectively denoted y, 
//which are chosen within three possible variables denoted Mes, dE and F.
//The variable which is not incorporated in y, is used as the control variable x.
//The distributions of discriminating variables and more details about the method 
//can be found in the TSPlot class description
//
// NOTE: This script requires a data file "TestSPlot_toyMC.dat".
//       This data file is not distributed in the standard ROOT binary tar file.
//       You can download it from ftp://root.cern.ch/root/TestSPlot_toyMC.dat
//
//Authors: Anna Kreshuk, Muriel Pivc

   TString dir = gSystem->UnixPathName(gInterpreter->GetCurrentMacroName());
   dir.ReplaceAll("TestSPlot.C","");
   dir.ReplaceAll("/./","/");
   TString dataFile = Form("%sTestSPlot_toyMC.dat",dir.Data());
   
   //Read the data and initialize a TSPlot object
   TTree *datatree = new TTree("datatree", "datatree");
   datatree->ReadFile(dataFile, 
                      "Mes/D:dE/D:F/D:MesSignal/D:MesBackground/D:dESignal/D:dEBackground/D:FSignal/D:FBackground/D",' ');

   TSPlot *splot = new TSPlot(0, 3, 5420, 2, datatree);

   //Set the selection for data tree
   //Note the order of the variables: 
   //first the control variables (not presented in this example),
   //then the 3 discriminating variables, then their probability distribution 
   //functions for the first species(signal) and then their pdfs for the 
   //second species(background)
   splot->SetTreeSelection(
      "Mes:dE:F:MesSignal:dESignal:FSignal:MesBackground:"
      "dEBackground:FBackground");

   //Set the initial estimates of the number of events in each species 
   //- used as initial parameter values for the Minuit likelihood fit
   Int_t ne[2];
   ne[0]=500; ne[1]=5000;
   splot->SetInitialNumbersOfSpecies(ne);

   //Compute the weights
   splot->MakeSPlot();

   //Fill the sPlots
   splot->FillSWeightsHists(25);

   //Now let's look at the sPlots
   //The first two histograms are sPlots for the Mes variable signal and 
   //background. dE and F were chosen as discriminating variables to determine 
   //N1 and N2, through a maximum Likelihood fit, and thus the sPlots for the 
   //control variable Mes, unknown to the fit, was contructed.
   //One can see that the sPlot for signal reproduces the PDF correctly, 
   //even when the latter vanishes.
   //
   //The lower two histograms are sPlots for the F variables signal and 
   //background. dE and Mes were chosen as discriminating variables to 
   //determine N1 and N2, through a maximum Likelihood fit, and thus the 
   //sPlots for the control variable F, unknown to the fit, was contructed.

   TCanvas *myc = new TCanvas("myc", 
   "sPlots of Mes and F signal and background", 800, 600);
   myc->SetFillColor(40);

   TPaveText *pt = new TPaveText(0.02,0.85,0.98,0.98);
   pt->SetFillColor(18);
   pt->SetTextFont(20);
   pt->SetTextColor(4);
   pt->AddText("sPlots of Mes and F signal and background,");
   pt->AddText("obtained by the tutorial TestSPlot.C on BABAR MC "
               "data (sPlot_toyMC.fit)");
   TText *t3=pt->AddText(
      "M. Pivk and F. R. Le Diberder, Nucl.Inst.Meth.A, physics/0402083");
   t3->SetTextColor(1);
   t3->SetTextFont(30);
   pt->Draw();

   TPad* pad1 = new TPad("pad1","Mes signal",0.02,0.43,0.48,0.83,33);
   TPad* pad2 = new TPad("pad2","Mes background",0.5,0.43,0.98,0.83,33);
   TPad* pad3 = new TPad("pad3", "F signal", 0.02, 0.02, 0.48, 0.41,33);
   TPad* pad4 = new TPad("pad4", "F background", 0.5, 0.02, 0.98, 0.41,33);
   pad1->Draw();
   pad2->Draw();
   pad3->Draw();
   pad4->Draw();

   pad1->cd();
   pad1->SetGrid();
   TH1D *sweight00 = splot->GetSWeightsHist(-1, 0, 0);
   sweight00->SetTitle("Mes signal");
   sweight00->SetStats(kFALSE);
   sweight00->Draw("e");
   sweight00->SetMarkerStyle(21);
//.........这里部分代码省略.........

void Plot_AM_events_07Sep_Susy_2_auto() {
 TString cutNameBefore = Form("Data/");

 //                            cut_variable 
 TString cutNameAfter  = Form("_2_6_Tot_temp");

 gROOT->LoadMacro("PlotVHqqHggH.C+");
 gInterpreter->ExecuteMacro("LatinoStyle2.C");

 TCanvas* c1 = new TCanvas("events","events",500,600);
 TFile* f = new TFile("~/Cern/Code/VBF/qqHWW/AnalysisPackage_qqHWWlnulnu/out_test_Latinos_07Sep2012_2000_Run2012AB_8TeV_SUSY.root");


 PlotVHqqHggH* hs = new PlotVHqqHggH();

 hs->setLumi(5.063);
 hs->setLabel("event");
 hs->addLabel("    #sqrt{s} = 8 TeV");

 TString name;

 std::vector<int> vectColourBkg;        
 std::vector<double> vectSystBkg;       
 std::vector<std::string> vectNameBkg;  
 std::vector<TH1F*> vectTHBkg;          

 std::vector<int> vectColourSig;      
 std::vector<double> vectSystSig;       
 std::vector<std::string> vectNameSig;  
 std::vector<TH1F*> vectTHSig;          

 ///==== signal (begin) ====

 name = Form("%sT2tt-350-70%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHSig.push_back ( (TH1F*) f->Get(name) );
 vectNameSig.push_back ("T2tt-stop350-lsp70");
 vectColourSig.push_back(6);

 name = Form("%sT2tt-500-70%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHSig.push_back ( (TH1F*) f->Get(name) );
 vectNameSig.push_back ("T2tt-stop500-lsp70");
 vectColourSig.push_back(97);

 name = Form("%sT2tt-350-100%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHSig.push_back ( (TH1F*) f->Get(name) );
 vectNameSig.push_back ("T2tt-stop350-lsp100");
 vectColourSig.push_back(70);

 name = Form("%sT2tt-500-100%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHSig.push_back ( (TH1F*) f->Get(name) );
 vectNameSig.push_back ("T2tt-stop500-lsp100");
 vectColourSig.push_back(65);

 name = Form("%sT2tt-500-200%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHSig.push_back ( (TH1F*) f->Get(name) );
 vectNameSig.push_back ("T2tt-stop500-lsp200");
 vectColourSig.push_back(5);

 name = Form("%sT2tt-200-0%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHSig.push_back ( (TH1F*) f->Get(name) );
 vectNameSig.push_back ("T2tt-stop200-lsp0");
 vectColourSig.push_back(7);

 ///==== signal (end)  ====

 name = Form("%sDATA%s",cutNameBefore.Data(),cutNameAfter.Data());
 hs->setDataHist ((TH1F*)f->Get(name));



 hs->setBlindBinSx(0);
 hs->setBlindBinDx(0);

 hs->setCutSx(-999,">");
 hs->setCutDx(-999,"<");


 ///==== background (begin)  ====

 name = Form("%sH-125%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHBkg.push_back ( (TH1F*) f->Get(name) );
 vectNameBkg.push_back ("H-125");
 vectColourBkg.push_back(633);
 vectSystBkg.push_back(0.00);


 name = Form("%sVV%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHBkg.push_back ( (TH1F*) f->Get(name) );
 vectNameBkg.push_back ("WZ/ZZ");
 vectColourBkg.push_back(858);
 vectSystBkg.push_back(0.00);

 name = Form("%sWJets%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHBkg.push_back ( (TH1F*) f->Get(name) );
 vectNameBkg.push_back ("W+jets");
 vectColourBkg.push_back(921);
 vectSystBkg.push_back(0.36);

 name = Form("%sTop%s",cutNameBefore.Data(),cutNameAfter.Data());
 vectTHBkg.push_back ( (TH1F*) f->Get(name) );
//.........这里部分代码省略.........

void Converts2x2Matrix(TCanvas* c, TCanvas* &c_output,Bool_t transpose){
    
    
    Double_t LimitLeft =0.15; //limit from left margin (15% of size)
    Double_t LimitBottom =0.10; //limit from bottom margin (15% of size)
    Double_t LimitRight =0.08; //limit from right margin (5% of size)
    Double_t LimitTop =0.05; //limit from top margin (5% of size)
    
    if(!c){
        cout << "Cannot read canvas!" << endl; return;
    }
    if(!c_output){
        cout << "output canvas not defined... quitting" << endl; return;
    }
    
    TString inputcanvasname = c->GetName();
    
    TPad * pd1 = (TPad*)c->FindObject(Form("%s_1",inputcanvasname.Data()));
    TPad * pd2 = (TPad*)c->FindObject(Form("%s_2",inputcanvasname.Data()));
    TPad * pd3 = (TPad*)c->FindObject(Form("%s_3",inputcanvasname.Data()));
    TPad * pd4 = (TPad*)c->FindObject(Form("%s_4",inputcanvasname.Data()));
    
    Double_t horizontalwidth=1.*(1-LimitLeft-LimitRight+LimitRight*LimitLeft)/(2-LimitLeft-LimitRight);
    Double_t verticalwidth=1.*(1-LimitBottom-LimitTop+LimitTop*LimitBottom)/(2-LimitBottom-LimitTop);
    
    //Double_t verticalwidth = (1-LimitBottom)*(1-LimitTop);
    
    Double_t x1 = 0;
    Double_t x2 = horizontalwidth/(1-LimitLeft);
    Double_t x3 = x2+horizontalwidth/(1-LimitRight);
    Double_t x4 = -1;
    
    Double_t y1 = 0;
    Double_t y2 = verticalwidth/(1-LimitBottom);
    Double_t y3 = y2+ verticalwidth/(1-LimitTop);
    Double_t y4 = -1;
    
    cout << "x =  " << x1 << "," << x2 <<" ,"<< x3 <<endl;
    cout << "y =  " << y1 << "," << y2 <<" ,"<< y3 <<endl;
    
  
    pd1->SetPad(x1,y3,x2,y2);
    pd2->SetPad(x2,y3,x3,y2);
    pd3->SetPad(x1,y2,x2,y1);
    pd4->SetPad(x2,y2,x3,y1);
    
    pd1->SetMargin(LimitLeft,0,0,LimitTop);
    pd2->SetMargin(0,LimitRight,0,LimitTop);
    pd3->SetMargin(LimitLeft,0,LimitBottom,0);
    pd4->SetMargin(0,LimitRight,LimitBottom,0);
    
    pd1->SetTicks();
    pd2->SetTicks();
    pd3->SetTicks();
    pd4->SetTicks();
    
    c_output->cd();
    pd1->Draw();
    pd2->Draw();
    pd3->Draw();
    pd4->Draw();
    
}

///////////////////////////////////////////////////////////////////////
//         TOOL BOX
//////////////////////////////////////////////////////////////////////
void makeMultiPanelCanvas(TCanvas*& canv,
                          const Int_t columns,
                          const Int_t rows,
                          const Float_t leftOffset,
                          const Float_t bottomOffset,
                          const Float_t leftMargin,
                          const Float_t bottomMargin,
                          const Float_t edge) {
   if (canv==0) {
      Error("makeMultiPanelCanvas","Got null canvas.");
      return;
   }
   canv->Clear();

   TPad* pad[columns][rows];
   Float_t Xlow[columns];
   Float_t Xup[columns];
   Float_t Ylow[rows];
   Float_t Yup[rows];
   Float_t PadWidth =
   (1.0-leftOffset)/((1.0/(1.0-leftMargin)) +
   (1.0/(1.0-edge))+(Float_t)columns-2.0);
   Float_t PadHeight =
   (1.0-bottomOffset)/((1.0/(1.0-bottomMargin)) +
   (1.0/(1.0-edge))+(Float_t)rows-2.0);
   Xlow[0] = leftOffset;
   Xup[0] = leftOffset + PadWidth/(1.0-leftMargin);
   Xup[columns-1] = 1;
   Xlow[columns-1] = 1.0-PadWidth/(1.0-edge);

   Yup[0] = 1;
   Ylow[0] = 1.0-PadHeight/(1.0-edge);
   Ylow[rows-1] = bottomOffset;
   Yup[rows-1] = bottomOffset + PadHeight/(1.0-bottomMargin);

   for(Int_t i=1;i<columns-1;i++) {
      Xlow[i] = Xup[0] + (i-1)*PadWidth;
      Xup[i] = Xup[0] + (i)*PadWidth;
   }
   Int_t ct = 0;
   for(Int_t i=rows-2;i>0;i--) {
      Ylow[i] = Yup[rows-1] + ct*PadHeight;
      Yup[i] = Yup[rows-1] + (ct+1)*PadHeight;
      ct++;
   }
   TString padName;
   for(Int_t i=0;i<columns;i++) {
      for(Int_t j=0;j<rows;j++) {
         canv->cd();
         padName = Form("p_%d_%d",i,j);
         pad[i][j] = new TPad(padName.Data(),padName.Data(),
            Xlow[i],Ylow[j],Xup[i],Yup[j]);
         if(i==0) pad[i][j]->SetLeftMargin(leftMargin);
         else pad[i][j]->SetLeftMargin(0);

         if(i==(columns-1)) pad[i][j]->SetRightMargin(edge);
         else pad[i][j]->SetRightMargin(0);

         if(j==0) pad[i][j]->SetTopMargin(edge);
         else pad[i][j]->SetTopMargin(0);

         if(j==(rows-1)) pad[i][j]->SetBottomMargin(bottomMargin);
         else pad[i][j]->SetBottomMargin(0);

         pad[i][j]->Draw();
         pad[i][j]->cd();
         pad[i][j]->SetNumber(columns*j+i+1);
      }
   }
}

void TriggerInputsForMuonEventCuts ( TString runListFilename, TString selectedInputs="", TString defaultStorage = "raw://" )
{
  AliCDBManager::Instance()->SetDefaultStorage(defaultStorage.Data());
  TObjArray inputsList;
  inputsList.SetOwner();
  
  TObjArray* selectedInputsList = selectedInputs.Tokenize(",");

  // Read input run list
  ifstream inFile(runListFilename.Data());
  TString srun = "";
  if ( inFile.is_open() ) {
    while ( ! inFile.eof() ) {
      srun.ReadLine(inFile,kFALSE);
      if ( ! srun.IsDigit() ) continue;
      
      // For each run, read trigger inputs from OCDB
      Int_t runNumber = srun.Atoi();
      AliCDBManager::Instance()->SetRun(runNumber);
      
      // Get trigger class configuration
      AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/CTP/Config");
      if ( ! entry ) continue;
      
      THashList* runInputs = new THashList();
      runInputs->SetOwner();
      runInputs->SetUniqueID((UInt_t)runNumber);
      AliTriggerConfiguration* trigConf = (AliTriggerConfiguration*)entry->GetObject();
      const TObjArray& trigInputsArray = trigConf->GetInputs();
      AliTriggerInput* trigInput = 0x0;
      TIter next(&trigInputsArray);
      while ( ( trigInput = static_cast<AliTriggerInput*>(next()) ) ) {
        if ( selectedInputsList->GetEntriesFast() > 0 && ! selectedInputsList->FindObject(trigInput->GetName()) ) continue;
        Int_t inputId = (Int_t)TMath::Log2(trigInput->GetMask());
        TObjString* currInput = new TObjString(trigInput->GetName());
        currInput->SetUniqueID(inputId);
        runInputs->Add(currInput);
      }
      inputsList.Add(runInputs);
    }
    inFile.close();
  }
  delete selectedInputsList;
  
  // Loop on the trigger inputs
  // and group runs with an equal list of inputs
  Int_t nentries = inputsList.GetEntries();
  TArrayI checkMask(nentries);
  checkMask.Reset(1);
  for ( Int_t irun=0; irun<nentries; irun++ ) {
    if ( checkMask[irun] == 0 ) continue;
    THashList* currList = static_cast<THashList*>(inputsList.At(irun));
    TString runRange = Form("Run range: %u", currList->GetUniqueID());
    for ( Int_t jrun=irun+1; jrun<nentries; jrun++ ) {
      if ( checkMask[jrun] == 0 ) continue;
      THashList* checkList = static_cast<THashList*>(inputsList.At(jrun));
      Bool_t isDifferent = kFALSE;
      for ( Int_t itrig=0; itrig<currList->GetEntries(); itrig++ ) {
        TObjString* currInput = static_cast<TObjString*>(currList->At(itrig));
        TObject* checkInput = checkList->FindObject(currInput->GetName());
        if ( ! checkInput || checkInput->GetUniqueID() != currInput->GetUniqueID() ) {
          isDifferent = kTRUE;
          break;
        }
      } // loop on trigger inputs
      if ( isDifferent ) continue;
      checkMask[jrun] = 0;
      runRange += Form(",%u", checkList->GetUniqueID());
    } // loop on runs
    
    TString outString = "\nSetTrigInputsMap(\"";
    for ( Int_t itrig=0; itrig<currList->GetEntries(); itrig++ ) {
      TObjString* currInput = static_cast<TObjString*>(currList->At(itrig));
      outString += Form("%s:%u,",currInput->GetString().Data(), currInput->GetUniqueID());
    }
    outString.Append("\");\n");
    outString.ReplaceAll(",\"","\"");
    outString += runRange;
    printf("%s\n", outString.Data());
  } // loop on runs
}

void run(Int_t runNumber)
{
    TStopwatch timer;
    timer.Start();
    
    TString strRunNumber = "run";
    strRunNumber += runNumber;
    const Int_t nev = -1;                                // number of events to read, -1 - untill CTRL+C
    const Int_t trigger = 1;                             // 1 - onspill, 2 - offspill. -1 - all
    TString inDir = "/Users/kresan/data/s438b/lmd/";     // directory with lmd files
    TString outDir = "/Users/kresan/data/s438b/data/";   // output directory

    TString filename = inDir + strRunNumber + "_*.lmd";
    TString outputFileName = outDir + strRunNumber + "_raw.root";                  // name of output file
    TString parFileName = outDir + "params_" + strRunNumber + "_raw.root";         // name of parameter file

    const Int_t updateRate = 150000;
    const Int_t minStats = 10000;                           // minimum number of entries for TCAL calibration
    const Int_t nModules = 800;                             // number of photomultipliers (for TCAL calibration)

    // Create source with unpackers ----------------------------------------------
    TString ntuple_options = "UNPACK:EVENTNO,UNPACK:TRIGGER,RAW";
    TString ucesb_dir = getenv("UCESB_DIR");
    TString ucesb_path = ucesb_dir + "/../upexps/s438b/s438b";

    EXT_STR_h101 ucesb_struct;
    R3BUcesbSource* source = new R3BUcesbSource(filename, ntuple_options,
                                                ucesb_path, &ucesb_struct, sizeof(ucesb_struct));
    source->SetMaxEvents(nev);
    source->AddReader(new R3BUnpackReader((EXT_STR_h101_unpack*)&ucesb_struct.unpack, offsetof(EXT_STR_h101, unpack)));
    source->AddReader(new R3BNeulandTacquilaReader((EXT_STR_h101_raw_nnp*)&ucesb_struct.nnp, offsetof(EXT_STR_h101, nnp)));
    source->AddReader(new R3BLosReader((EXT_STR_h101_LOS*)&ucesb_struct.los, offsetof(EXT_STR_h101, los)));
    // ---------------------------------------------------------------------------

    // Create online run ---------------------------------------------------------
    FairRunOnline* run = new FairRunOnline(source);
    run->SetRunId(runNumber);
    run->SetOutputFile(outputFileName.Data());
    //run->ActivateHttpServer(100000);
    // ---------------------------------------------------------------------------

    // Create ALADIN field map ---------------------------------------------------
    R3BAladinFieldMap* magField = new R3BAladinFieldMap("AladinMaps");
    Double_t fMeasCurrent = 2500.; // I_current [A]
    magField->SetCurrent(fMeasCurrent);
    magField->SetScale(1.);
    run->SetField(magField);
    // ---------------------------------------------------------------------------

    // TCAL ----------------------------------------------------------------------
    R3BNeulandMapped2CalPar* tcalFill = new R3BNeulandMapped2CalPar("TcalFill");
    tcalFill->SetUpdateRate(updateRate);
    tcalFill->SetMinStats(minStats);
    tcalFill->SetTrigger(trigger);
    tcalFill->SetNofModules(nModules);
    run->AddTask(tcalFill);

    R3BLosMapped2CalPar* losTcalFill = new R3BLosMapped2CalPar("LosTcalFill");
    losTcalFill->SetUpdateRate(updateRate);
    losTcalFill->SetMinStats(minStats);
    losTcalFill->SetNofModules(20, 4);
    run->AddTask(losTcalFill);
    // ---------------------------------------------------------------------------

    // Add analysis task ---------------------------------------------------------
    R3BNeulandMappedHist* ana = new R3BNeulandMappedHist("LandRawAna", 1);
    run->AddTask(ana);
    // ---------------------------------------------------------------------------

    // Initialize ----------------------------------------------------------------
    run->Init();
    //((TTree*)gFile->Get("cbmsim"))->SetMaxTreeSize(maxSize);
    FairLogger::GetLogger()->SetLogScreenLevel("INFO");
    // ---------------------------------------------------------------------------

    // Runtime data base ---------------------------------------------------------
    FairRuntimeDb* rtdb = run->GetRuntimeDb();
    R3BFieldPar* fieldPar = (R3BFieldPar*)rtdb->getContainer("R3BFieldPar");
    fieldPar->SetParameters(magField);
    fieldPar->setChanged();
    Bool_t kParameterMerged = kTRUE;
    FairParRootFileIo* parOut = new FairParRootFileIo(kParameterMerged);
    parOut->open(parFileName);
    rtdb->setOutput(parOut);
    rtdb->print();
    // ---------------------------------------------------------------------------

    // Run -----------------------------------------------------------------------
    if(nev < 0)
    {
        run->Run(nev, 0);
    }
    else
    {
        run->Run(0, nev);
    }
    rtdb->saveOutput();
    delete run;
    // ---------------------------------------------------------------------------

//.........这里部分代码省略.........

void computeAccGenWm(const TString conf,             // input file
                     const TString outputDir,        // output directory
		     const Int_t   charge            // 0 = inclusive, +1 = W+, -1 = W-
) {
  gBenchmark->Start("computeAccGenWm");

  //--------------------------------------------------------------------------------------------------------------
  // Settings 
  //============================================================================================================== 

  const Double_t PT_CUT     = 25;
  const Double_t ETA_CUT    = 2.1;
  const Double_t ETA_BARREL = 1.2;
  const Double_t ETA_ENDCAP = 1.2;

  const Int_t BOSON_ID  = 24;
  Int_t LEPTON_ID = 13;
  if (charge==1) LEPTON_ID=-13;

  //--------------------------------------------------------------------------------------------------------------
  // Main analysis code 
  //==============================================================================================================  

  vector<TString> fnamev;  // file name per input file
  vector<TString> labelv;  // TLegend label per input file
  vector<Int_t>   colorv;  // plot color per input file
  vector<Int_t>   linev;   // plot line style per input file

  //
  // parse .conf file
  //
  ifstream ifs;
  ifs.open(conf.Data());
  assert(ifs.is_open());
  string line;
  while(getline(ifs,line)) {
    if(line[0]=='#') continue;
    
    string fname;
    Int_t color, linesty;
    stringstream ss(line);
    ss >> fname >> color >> linesty;
    string label = line.substr(line.find('@')+1);
    fnamev.push_back(fname);
    labelv.push_back(label);
    colorv.push_back(color);
    linev.push_back(linesty);
  }
  ifs.close();

  // Create output directory
  gSystem->mkdir(outputDir,kTRUE);
  
  // Data structures to store info from TTrees
  baconhep::TGenEventInfo *gen = new baconhep::TGenEventInfo();
  TClonesArray *genPartArr     = new TClonesArray("baconhep::TGenParticle");
  
  TFile *infile=0;
  TTree *eventTree=0;

  // Variables to store acceptances and uncertainties (per input file)
  vector<Double_t> nEvtsv, nSelv, nSelBv, nSelEv;
  vector<Double_t> accv, accBv, accEv;
  vector<Double_t> accErrv, accErrBv, accErrEv;

  //
  // loop through files
  //
  for(UInt_t ifile=0; ifile<fnamev.size(); ifile++) {  

    // Read input file and get the TTrees
    cout << "Processing " << fnamev[ifile] << " ..." << endl;
    infile = TFile::Open(fnamev[ifile]); 
    assert(infile);
  
    eventTree = (TTree*)infile->Get("Events"); 
    assert(eventTree);
    eventTree->SetBranchAddress("GenEvtInfo", &gen); TBranch *genBr = eventTree->GetBranch("GenEvtInfo"); 
    eventTree->SetBranchAddress("GenParticle", &genPartArr); TBranch *partBr = eventTree->GetBranch("GenParticle");

    nEvtsv.push_back(0);
    nSelv.push_back(0);
    nSelBv.push_back(0);
    nSelEv.push_back(0);
    
    //
    // loop over events
    //    
    for(UInt_t ientry=0; ientry<eventTree->GetEntries(); ientry++) {
    //for(UInt_t ientry=0; ientry<10000000; ientry++) {
      genBr->GetEntry(ientry);
      genPartArr->Clear(); partBr->GetEntry(ientry);

      TLorentzVector *vec=0, *lep1=0, *lep2=0;
      if (charge==-1 && toolbox::flavor(genPartArr, -BOSON_ID, vec, lep1, lep2)!=LEPTON_ID) {
        continue;
      }
      if (charge==1 && toolbox::flavor(genPartArr, BOSON_ID, vec, lep1, lep2)!=LEPTON_ID) {
        continue;
      }
//.........这里部分代码省略.........

void plotphi(int top, bool phi=kFALSE){
  if (setup("vm","","pol__")==kFALSE) return;
  gStyle->SetOptFit(1111); //pcev = 1111
  //gStyle->SetTitleW(1.5);

  int itop = top-1;
  TFile* fy = _fyexp[itop];

  int nq2wbins = fy->GetNkeys();
  TIter nextkey(fy->GetListOfKeys());
  TKey *key;
  int iq2wbin = 0;
  while (key = (TKey*)nextkey()) {
    TString Q2Wdirname = key->GetName();
    if(Q2Wdirname.EqualTo("hYW_Dir") || Q2Wdirname.EqualTo("hYW"))continue;
    cout << "Q2Wdirname = " << Q2Wdirname << endl;

    THStack* hs[10];
    TF1 *fphi = new TF1("fphi", "([0] + [1]*cos(x*TMath::DegToRad()) + [2]*cos(2*x*TMath::DegToRad()) + [3]*sin(x*TMath::DegToRad()))",0, 360);
    fphi->SetParameter(0,1);
    fphi->SetParameter(1,10);
    fphi->SetParameter(2,20);
    fphi->SetParameter(3,100);
    fphi->SetParName(0, "A");
    fphi->SetParName(1, "B");
    fphi->SetParName(2, "C");
    fphi->SetParName(3, "hPR");
    for (int i = 0; i < 10; i++)
    {
      TString hname_unp = TString::Format("%s/hPhi/Varset1/theta/hphi_%02d",Q2Wdirname.Data(),i+1);
      TString hname_pos = TString::Format("%s/hPhi_POS/Varset1/theta/hphi_%02d",Q2Wdirname.Data(),i+1);
      TString hname_neg = TString::Format("%s/hPhi_NEG/Varset1/theta/hphi_%02d",Q2Wdirname.Data(),i+1);
      cout << "hname_unp = " << hname_unp << endl;
      cout << "hname_pos = " << hname_pos << endl;
      cout << "hname_neg = " << hname_neg << endl;
      TH1D* hunp = (TH1D*)fy->Get(hname_unp);
      TH1D* hpos = (TH1D*)fy->Get(hname_pos);
      TH1D* hneg = (TH1D*)fy->Get(hname_neg);
      if (hunp==NULL || hpos==NULL || hneg==NULL) cout << "histogram not found" << endl;
      hunp->SetLineColor(kBlack);
      hpos->SetLineColor(kRed);
      hneg->SetLineColor(kBlue);
      hunp->Fit(fphi);
      hpos->Fit(fphi);
      hneg->Fit(fphi);

      //! Modify titles
      TObjArray *tarr;
      char t[100];

      TPaveText *ptunp = new TPaveText(0.0, 0.9, 1.0, 1, "NDC");
      TString tunp = hunp->GetTitle();
      tarr = tunp.Tokenize("|");
      sprintf(t, "%s:%s:%s", tarr->At(1)->GetName(), tarr->At(2)->GetName(), tarr->At(3)->GetName());
      ptunp->AddText(tarr->At(0)->GetName());
      ptunp->AddText(t);
      hunp->SetTitle("");

      TPaveText *ptpos = new TPaveText(0.0, 0.9, 1.0, 1, "NDC");
      TString tpos = hpos->GetTitle();
      tarr = tpos.Tokenize("|");
      sprintf(t, "%s:%s:%s", tarr->At(1)->GetName(), tarr->At(2)->GetName(), tarr->At(3)->GetName());
      ptpos->AddText(tarr->At(0)->GetName());
      ptpos->AddText(t);
      hpos->SetTitle("");

      TPaveText *ptneg = new TPaveText(0.0, 0.9, 1.0, 1, "NDC");
      TString tneg = hneg->GetTitle();
      tarr = tneg.Tokenize("|");
      sprintf(t, "%s:%s:%s", tarr->At(1)->GetName(), tarr->At(2)->GetName(), tarr->At(3)->GetName());
      ptneg->AddText(tarr->At(0)->GetName());
      ptneg->AddText(t);
      hneg->SetTitle("");

      //! Draw histograms
      TCanvas *c = new TCanvas(hpos->GetName(),hpos->GetName(),900, 600);
      c->Divide(3,1);
      c->cd(1);
      hunp->Draw();
      ptunp->Draw();
      c->cd(2);
      hpos->Draw();
      ptpos->Draw();
      c->cd(3);
      hneg->Draw();
      ptneg->Draw();
      TString dir = TString::Format("./polobs/top%d/%s/Varset1/theta", top, Q2Wdirname.Data());
      gSystem->mkdir(dir,1);
      c->SaveAs(TString::Format("%s/%s.png", dir.Data(), c->GetName()));
      c->Close();
      c->Delete();
    }
  }
}

void plotcompxsec(TString xsectype/*= "vm"*/, bool ploty/*=kFALSE*/,bool comp/*=kFALSE*/){
  TString methodB = TString::Format("Yield - Method B {%s}",_topNames[4].Data());
  TString methodA = TString::Format("Yield - Method A");

  TH1F* hYW_Dir[5];
  TH1F* hYW_5D[5];
  TH1F* hYWnorm_Dir[5];
  TH1F* hYWnorm_5D[5];
  TLegend* l_Dir = new TLegend(0.1,0.7,0.5,0.9);
  TLegend* l_5D  = new TLegend(0.1,0.7,0.5,0.9);
  TString hname_hYW_Dir  = TString::Format("hYW_Dir/Varset1/hYW_ACC_CORR");
  TString hname_hYW_5D   = TString::Format("hYW/Varset1/hYW");

  TCanvas* c_Dir = new TCanvas(TString::Format("Dir-%s:%s", xsectype.Data(), _q2w_bng.Data()), 
                               TString::Format("Dir-%s:%s", xsectype.Data(), _q2w_bng.Data()) );
  c_Dir->SetGrid(); 
  TCanvas* c_5D = new TCanvas(TString::Format("5D-%s:%s", xsectype.Data(), _q2w_bng.Data()), 
                              TString::Format("5D-%s:%s", xsectype.Data(), _q2w_bng.Data()) );
  c_5D->SetGrid();

  if (setup(xsectype,Q2W_BNG1)==kFALSE) return;

  //! For each Top: hYW_Dir,hYW_5D --> hYWnorm_Dir,hYWnorm_5D & then draw
  for (int iTop=0;iTop<5;iTop++){
    if (iTop!=1 && iTop!=4) continue;
    if (_fyexp[iTop]==NULL) continue;
    cout << "Top = " << iTop << endl;
    hYW_Dir[iTop] = (TH1F*)_fyexp[iTop]->Get(hname_hYW_Dir);
    hYW_5D[iTop]  = (TH1F*)_fyexp[iTop]->Get(hname_hYW_5D);
    if (hYW_Dir[iTop]==NULL || hYW_5D[iTop]==NULL) continue;
    hYW_Dir[iTop]->SetMarkerStyle(20+iTop);
    hYW_Dir[iTop]->SetMarkerColor(2+iTop);
    hYW_Dir[iTop]->SetMaximum(150000);
    hYW_5D[iTop]->SetMarkerStyle(20+iTop); //+5
    hYW_5D[iTop]->SetMarkerColor(2+iTop);  //+5
    hYW_5D[iTop]->SetMaximum(500000);
    //Normalize yields to LUM*vgflux*binw_w*binw_q2
    hYWnorm_Dir[iTop] = normalizeYield(hYW_Dir[iTop]);
    hYWnorm_Dir[iTop]->SetMaximum(10);
    hYWnorm_5D[iTop]  = normalizeYield(hYW_5D[iTop]);
    hYWnorm_5D[iTop]->SetMaximum(10);

   if (iTop==1){
     c_5D->cd();
     l_5D->SetHeader( TString::Format("%s - 5D", XSECTITLE.Data()) );
     l_5D->AddEntry( hYWnorm_5D[iTop], TString::Format("5D-%s", _topNames[iTop].Data()) );
     hYWnorm_5D[iTop]->SetMinimum(0);
     ploty?hYW_5D[iTop]->Draw("p"):hYWnorm_5D[iTop]->Draw("p");

     c_Dir->cd();
     l_Dir->SetHeader( TString::Format("%s - 5D", XSECTITLE.Data()) );
     l_Dir->AddEntry( hYWnorm_Dir[iTop], TString::Format("Direct-%s", _topNames[iTop].Data()) );
     hYWnorm_Dir[iTop]->SetMinimum(0);
     ploty?hYW_Dir[iTop]->Draw("p"):hYWnorm_Dir[iTop]->Draw("p");
   }else{
     c_5D->cd();
     l_5D->AddEntry( hYWnorm_5D[iTop], TString::Format("5D-%s", _topNames[iTop].Data()) );
     ploty?hYW_5D[iTop]->Draw("p sames"):hYWnorm_5D[iTop]->Draw("p sames");

     c_Dir->cd();
     l_Dir->AddEntry( hYWnorm_Dir[iTop], TString::Format("Direct-%s", _topNames[iTop].Data()) );
     ploty?hYW_Dir[iTop]->Draw("p sames"):hYWnorm_Dir[iTop]->Draw("p sames");
   }
   c_5D->cd();
   l_5D->Draw("same");
   c_Dir->cd();
   l_Dir->Draw("same"); 
  }

  if (setup(xsectype,Q2W_BNG2)==kFALSE) return;

  //! For each Top: hYW_Dir,hYW_5D --> hYWnorm_Dir,hYWnorm_5D & then draw
  for (int iTop=0;iTop<5;iTop++){
    if (iTop!=1 && iTop!=4) continue;
    if (_fyexp[iTop]==NULL) continue;
    cout << "Top = " << iTop << endl;
    hYW_Dir[iTop] = (TH1F*)_fyexp[iTop]->Get(hname_hYW_Dir);
    hYW_5D[iTop]  = (TH1F*)_fyexp[iTop]->Get(hname_hYW_5D);
    if (hYW_Dir[iTop]==NULL || hYW_5D[iTop]==NULL) continue;
    hYW_Dir[iTop]->SetMarkerStyle(20+iTop);
    hYW_Dir[iTop]->SetMarkerColor(2+iTop);
    hYW_Dir[iTop]->SetMaximum(150000);
    hYW_5D[iTop]->SetMarkerStyle(20+iTop); //+5
    hYW_5D[iTop]->SetMarkerColor(2+iTop);  //+5
    hYW_5D[iTop]->SetMaximum(500000);
    //Normalize yields to LUM*vgflux*binw_w*binw_q2
    hYWnorm_Dir[iTop] = normalizeYield(hYW_Dir[iTop]);
    hYWnorm_Dir[iTop]->SetMaximum(10);
    hYWnorm_5D[iTop]  = normalizeYield(hYW_5D[iTop]);
    hYWnorm_5D[iTop]->SetMaximum(10);

   /*if (iTop==0){
     c_5D->cd();
     l_5D->SetHeader( TString::Format("%s - 5D", XSECTITLE.Data()) );
     l_5D->AddEntry( hYWnorm_5D[iTop], TString::Format("5D-%s", _topNames[iTop].Data()) );
     hYWnorm_5D[iTop]->SetMinimum(0);
     ploty?hYW_5D[iTop]->Draw("p"):hYWnorm_5D[iTop]->Draw("p");

     c_Dir->cd();
     l_Dir->SetHeader( TString::Format("%s - 5D", XSECTITLE.Data()) );
//.........这里部分代码省略.........

Bool_t setup(TString xsectype/*= "vm"*/, TString q2w_bng/*=""*/, bool pol/*=false*/){
  Bool_t ret = kTRUE;

  //gStyle->SetOptStat(0);

  TString pwd = gSystem->pwd();
  cout <<"pwd = "<<pwd << endl;

  TString path;
  if(q2w_bng=="")
  {
    if      (pwd.Contains("1.4-1.5__1.6-1.8"))  _q2w_bng = Q2W_BNG1;    
    else if (pwd.Contains("1.9-2.5__1.3-1.9"))  _q2w_bng = Q2W_BNG2;

    if      (pwd.Contains("experiment"))      _dtype = "exp";
    else if (pwd.Contains("simulation"))      _dtype = "sim";
    path=".";
  }else{
    _q2w_bng=q2w_bng;
    printf("q2wbng = %s\n", _q2w_bng.Data());
    if(_q2w_bng.Contains(Q2W_BNG1))
    {
      path="/data/trivedia/e1f/ana2pi/experiment/e1f.goldenruns.qskim_071513/Q2W__1.4-1.5__1.6-1.8";
    }else if (_q2w_bng.Contains(Q2W_BNG2)){
      path="/data/trivedia/e1f/ana2pi/experiment/e1f.goldenruns.qskim_071513/Q2W__1.9-2.5__1.3-1.9";
    }else 
    {
      printf("Could not determine path\n");
      return kFALSE; 
    }
    _dtype="exp";
  }

  if (pol) _dtype = "pol__exp";
  

  cout << "xsectype = " << xsectype << endl;
  if (xsectype=="at"){
    copy(AT_TOPNAMES,AT_TOPNAMES+5, _topNames);
    _fyexp[0] = TFile::Open(TString::Format("%s/1__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fyexp[1] = TFile::Open(TString::Format("%s/2__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fyexp[2] = TFile::Open(TString::Format("%s/3__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fyexp[3] = TFile::Open(TString::Format("%s/4__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fyexp[4] = TFile::Open(TString::Format("%s/1:2:3:4__%s__%s.root", path.Data(), _q2w_bng.Data(), _dtype.Data()));
    _fysim[0] = TFile::Open(TString::Format("%s/simdir/1__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
    _fysim[1] = TFile::Open(TString::Format("%s/simdir/2__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
    _fysim[2] = TFile::Open(TString::Format("%s/simdir/3__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
    _fysim[3] = TFile::Open(TString::Format("%s/simdir/4__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
    _fysim[4] = TFile::Open(TString::Format("%s/simdir/1:2:3:4__%s__%s.root", _q2w_bng.Data(), "sim"));
  }else if (xsectype=="vm"){
    copy(VM_TOPNAMES,VM_TOPNAMES+5, _topNames);
    _fyexp[0] = TFile::Open(TString::Format("%s/1__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fyexp[1] = TFile::Open(TString::Format("%s/2:1__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fyexp[2] = TFile::Open(TString::Format("%s/3:1__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fyexp[3] = TFile::Open(TString::Format("%s/4:1__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fyexp[4] = TFile::Open(TString::Format("%s/1:2:3:4__%s__%s.root", path.Data(),_q2w_bng.Data(), _dtype.Data()));
    _fysim[0] = TFile::Open(TString::Format("%s/simdir/1__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
    _fysim[1] = TFile::Open(TString::Format("%s/simdir/2:1__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
    _fysim[2] = TFile::Open(TString::Format("%s/simdir/3:1__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
    _fysim[3] = TFile::Open(TString::Format("%s/simdir/4:1__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
    _fysim[4] = TFile::Open(TString::Format("%s/simdir/1:2:3:4__%s__%s.root", path.Data(),_q2w_bng.Data(), "sim"));
  }else {
    printf("Could not determine xsectype\n");
    ret=kFALSE;
  }

  TString q2bng = _q2w_bng.Tokenize("__")->At(0)->GetName();
  _q2min = ((TString)q2bng.Tokenize("-")->At(1)->GetName()).Atof();
  _q2max = ((TString)q2bng.Tokenize("-")->At(2)->GetName()).Atof();
  _dq2 = _q2max - _q2min;
  
  return ret;
}

Int_t SETUP(TList *inputList = NULL) {

  TString aliRootDir, aliPhysicsDir;

  if (gProof && !gProof->IsMaster()) {

    //
    // On client
    //

    gMessTag = "client";
    aliRootDir = gSystem->Getenv("ALICE_ROOT");  // NULL --> ""
    aliPhysicsDir = gSystem->Getenv("ALICE_PHYSICS");  // NULL --> ""

    if (aliRootDir.IsNull()) {
      // ALICE_ROOT is mandatory; ALICE_PHYSICS is optional
      ::Error(gMessTag.Data(), "ALICE_ROOT environment variable not defined on client");
      return -1;
    }

    if (aliPhysicsDir.EqualTo(aliRootDir)) {
      // AliPhysics does not exist if it points to the same directory of AliRoot
      aliPhysicsDir = "";
    }

    ::Info(gMessTag.Data(), "Enabling ALICE software on client");

  }
  else {

    //
    // On master/workers
    //

    gMessTag = gSystem->HostName();
    if (gProof && gProof->IsMaster()) {
      gMessTag.Append("(master)");
    }
    else {
      gMessTag.Append("(worker)");
    }

    // Here we have two working modes: AAF and VAF.
    //  - AAF mode: look for [email protected](AliRoot|AliPhysics)::<version>
    //  - VAF mode: any other PARfile name

    TString buf;
    buf = gSystem->BaseName(gSystem->pwd());
    TPMERegexp reAaf("^[email protected](AliRoot|AliPhysics)::(.*)$");

    if (reAaf.Match(buf) == 3) {

      // AliRoot enabled from a metaparfile whose name matches
      // [email protected]::<version>: set up ALICE_ROOT environment variable
      // accordingly from there.
      // Note: this is the AAF case.

      TString swName = reAaf[1].Data();
      TString swVer = reAaf[2].Data();

      // Get ALICE_ROOT and ALICE_PHYSICS from the env set by Modules
      buf.Form( ". /cvmfs/alice.cern.ch/etc/login.sh && "
        "eval `alienv printenv [email protected]%s::%s` && "
        "echo \"$ALICE_ROOT\" ; "
        "echo \"$ALICE_PHYSICS\"", swName.Data(), swVer.Data() );
      buf = gSystem->GetFromPipe( buf.Data() );

      // Output on two lines
      TString tok;
      Ssiz_t from = 0;
      UInt_t count = 0;
      while ( buf.Tokenize(tok, from, "\n") ) {
        if (count == 0) {
          aliRootDir = tok;
        }
        else if (count == 1) {
          aliPhysicsDir = tok;
        }
        else {
          break;
        }
        count++;
      }

      // Set (or override) environment for AliRoot and AliPhysics.
      // LD_LIBRARY_PATH: current working directory always has precedence.
      // Note: supports both current $ALICE_ROOT/lib and legacy
      //       $ALICE_ROOT/lib/tgt_<arch> format.
      gSystem->Setenv("ALICE_ROOT", aliRootDir.Data());
      gSystem->SetDynamicPath(
        Form(".:%s/lib:%s/lib/tgt_%s:%s", aliRootDir.Data(), aliRootDir.Data(),
          gSystem->GetBuildArch(), gSystem->GetDynamicPath()) );

      ::Info(gMessTag.Data(),
        "Enabling %s %s on PROOF (AAF mode)...", swName.Data(), swVer.Data());

      if (!aliPhysicsDir.IsNull()) {
        gSystem->Setenv("ALICE_PHYSICS", aliPhysicsDir.Data());
        gSystem->SetDynamicPath(
          Form("%s/lib:%s", aliPhysicsDir.Data(), gSystem->GetDynamicPath()) );
//.........这里部分代码省略.........