C++ TCanvas类代码示例

void VetoEnergy3(Int_t Lo, Int_t Hi)
{
  TCanvas* Window;
  TTimer* Refresh;
  TH1F* Energy;
  TF1* Fit;
  Char_t Buff[256];
  Char_t* Keyb;
  Double_t Pos[438];
  Double_t Lower;
  Double_t Param[3];
  Int_t OK;

  Window = new TCanvas();
  Window->SetCrosshair();
  Window->ToggleEventStatus();
  Refresh = new TTimer("Flush()", 50, kFALSE);
  for(Int_t ch=Lo; ch<Hi+1; ch++)
  {
    sprintf(Buff, "Veto_Energy%d", ch);
    Energy = (TH1F*)gROOT->FindObject(Buff);
    Energy->Draw();
    Window->Update();

    Fit = new TF1("F1", "gaus", 0.7, 2.0);
    Energy->Fit("F1", "RQ");
    Fit->Draw("same");
    Window->Update();
    Fit->GetParameters(Param);

    Refresh->TurnOn();
    Refresh->Reset();
    sprintf(Buff, "OK for %d? ", ch);
    Keyb = Getline(Buff);
    Refresh->TurnOff();
    OK = atoi(Keyb);
    if(OK) goto Done;

Refit:
    Refresh->TurnOn();
    Refresh->Reset();
    sprintf(Buff, "Lower bound for %d: ", ch);
    Keyb = Getline(Buff);
    Refresh->TurnOff();
    Lower = atof(Keyb);

    Fit = new TF1("F1", "gaus", Lower, 2.0);
    Energy->Fit("F1", "RQ");
    Fit->Draw("same");
    Window->Update();
    Fit->GetParameters(Param);

    Refresh->TurnOn();
    Refresh->Reset();
    sprintf(Buff, "OK for %d? ", ch);
    Keyb = Getline(Buff);
    Refresh->TurnOff();
    OK = atoi(Keyb);
    if(!OK) goto Refit;

Done:
    Pos[ch] = Param[1];
  }
  for(Int_t ch=Lo; ch<Hi+1; ch++)
    printf("%8.6f\n", Pos[ch]);
}

//.........这里部分代码省略.........
  for ( unsigned i = 0; i < v_h_esum.size(); i++ )
    {

      TH1F* h_raw = v_h_esum.at( i );
      TH1F* h_digi = (TH1F*)h_raw->Clone("h_digi");
      h_digi->Reset();

      for ( int idx1 = 0; idx1 < h_raw->GetNbinsX(); idx1++ )
	{
	  float eraw_i = h_raw->GetXaxis()->GetBinCenter( idx1 );
	  float craw_i = h_raw->GetBinContent( idx1 );

	  for ( int idx2 = 0; idx2 < craw_i; idx2++ )
	    {
	      h_digi->Fill( rand->Poisson( NperMeV * EtoMeV * eraw_i ) / ( NperMeV * EtoMeV ) );
	    }
	}

      v_h_esum_v2.push_back( h_digi );
  
    }
  
  */
  /*----------------Gaussian Fits-------------------*/

  vector< float > v_gauss_mean;
  vector< float > v_gauss_merr;
  vector< float > v_gauss_sigm;
  vector< float > v_gauss_serr;

    for ( unsigned i = 0; i < v_h_esum.size(); i++ )
    //for ( unsigned i = 0; i < v_h_esum_v2.size(); i++ )
    {
      TCanvas *cfit = new TCanvas();
      cfit->cd();
      
      //without AdHoc Digitization
               v_h_esum.at( i )->Fit("gaus","","");
	    
	    v_gauss_mean.push_back( v_h_esum.at( i )->GetFunction("gaus")->GetParameter( 1 )/200 );
	    v_gauss_merr.push_back( v_h_esum.at( i )->GetFunction("gaus")->GetParError(1)/200);
	    v_gauss_sigm.push_back( v_h_esum.at( i )->GetFunction("gaus")->GetParameter( 2 ) / v_h_esum.at( i )->GetFunction("gaus")->GetParameter( 1 ) );

	    v_gauss_serr.push_back( v_gauss_sigm.at( i )*(sqrt(((v_gauss_merr.at(i))/(v_gauss_mean.at(i))) + ((v_h_esum.at( i )->GetFunction("gaus")->GetParError(2))/(v_h_esum.at( i )->GetFunction("gaus")->GetParameter(2))))));

	    //v_gauss_sigm.push_back( v_h_esum.at( i )->GetFunction("gaus")->GetParameter( 2 ) / v_pbeam.at( i ));
	  
	    /*          
	    //with AdHoc Digitization
	    
	    v_h_esum_v2.at( i )->Fit("gaus","","");

	    v_gauss_mean.push_back( v_h_esum_v2.at( i )->GetFunction("gaus")->GetParameter( 1 ) );
	    v_gauss_merr.push_back( v_h_esum_v2.at( i )->GetFunction("gaus")->GetParError(1));
      
	    //v_gauss_sigm.push_back( v_h_esum_v2.at( i )->GetFunction("gaus")->GetParameter( 2 ) / v_pbeam.at( i ));

	    v_gauss_sigm.push_back( v_h_esum_v2.at( i )->GetFunction("gaus")->GetParameter( 2 ) / v_h_esum_v2.at( i )->GetFunction("gaus")->GetParameter( 1 ) );

	    //v_gauss_serr.push_back( v_h_esum_v2.at(i)->GetFunction("gaus")->GetParError(2));

            v_gauss_serr.push_back( v_gauss_sigm.at( i )*(sqrt(((v_gauss_merr.at(i))/(v_gauss_mean.at(i))) + ((v_h_esum_v2.at( i )->GetFunction("gaus")->GetParError(2))/(v_h_esum_v2.at( i )->GetFunction("gaus")->GetParameter(2))))));
	    */
	    
    }
	    

TF1* fit(Float_t ptmin, Float_t ptmax)
{
  TCanvas* c = new TCanvas(Form("c_%.0f_%.0f",ptmin,ptmax),"",600,600);
  TFile* infile = new TFile(Form("%s_%s_cent_%.0f_%.0f_pt_%.0f_%.0f.root",infname.Data(),collisionsystem.Data(),centmin,centmax,ptmin,ptmax));
  TH1D* h = (TH1D*)infile->Get("h");                    h->SetName(Form("h_%.0f_%.0f",ptmin,ptmax));
  TH1D* hMCSignal = (TH1D*)infile->Get("hMCSignal");    hMCSignal->SetName(Form("hMCSignal_%.0f_%.0f",ptmin,ptmax));
  TH1D* hMCSwapped = (TH1D*)infile->Get("hMCSwapped");  hMCSwapped->SetName(Form("hMCSwapped_%.0f_%.0f",ptmin,ptmax));
  TF1* f = new TF1(Form("f_%.0f_%.0f",ptmin,ptmax),"[0]*([7]*([9]*Gaus(x,[1],[2]*(1+[11]))/(sqrt(2*3.14159)*[2]*(1+[11]))+(1-[9])*Gaus(x,[1],[10]*(1+[11]))/(sqrt(2*3.14159)*[10]*(1+[11])))+(1-[7])*Gaus(x,[1],[8]*(1+[11]))/(sqrt(2*3.14159)*[8]*(1+[11])))+[3]+[4]*x+[5]*x*x+[6]*x*x*x", 1.7, 2.0);

  f->SetParLimits(4,-1000,1000);
  f->SetParLimits(10,0.001,0.05);
  f->SetParLimits(2,0.01,0.1);
  f->SetParLimits(8,0.02,0.2);
  f->SetParLimits(7,0,1);
  f->SetParLimits(9,0,1);
  
  f->SetParameter(0,setparam0);
  f->SetParameter(1,setparam1);
  f->SetParameter(2,setparam2);
  f->SetParameter(10,setparam10);
  f->SetParameter(9,setparam9);

  f->FixParameter(8,setparam8);
  f->FixParameter(7,1);
  f->FixParameter(1,fixparam1);
  f->FixParameter(3,0);
  f->FixParameter(4,0);
  f->FixParameter(5,0);
  f->FixParameter(6,0);
  f->FixParameter(11,0);
  h->GetEntries();
  
  hMCSignal->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"q","",minhisto,maxhisto);
  hMCSignal->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"q","",minhisto,maxhisto);
  f->ReleaseParameter(1);
  hMCSignal->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L q","",minhisto,maxhisto);
  hMCSignal->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L q","",minhisto,maxhisto);
  hMCSignal->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L m q","",minhisto,maxhisto);
  
  f->FixParameter(1,f->GetParameter(1));
  f->FixParameter(2,f->GetParameter(2));
  f->FixParameter(10,f->GetParameter(10));
  f->FixParameter(9,f->GetParameter(9));
  f->FixParameter(7,0);
  f->ReleaseParameter(8);
  f->SetParameter(8,setparam8);
  
  hMCSwapped->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L q","",minhisto,maxhisto);
  hMCSwapped->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L q","",minhisto,maxhisto);
  hMCSwapped->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L q","",minhisto,maxhisto);
  hMCSwapped->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L m q","",minhisto,maxhisto);
  
  f->FixParameter(7,hMCSignal->Integral(0,1000)/(hMCSwapped->Integral(0,1000)+hMCSignal->Integral(0,1000)));
  f->FixParameter(8,f->GetParameter(8));
  f->ReleaseParameter(3);
  f->ReleaseParameter(4);
  f->ReleaseParameter(5);
  f->ReleaseParameter(6);

  f->SetLineColor(kRed);
  
  h->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"q","",minhisto,maxhisto);
  h->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"q","",minhisto,maxhisto);
  f->ReleaseParameter(1);
  f->SetParLimits(1,1.85,1.90);
  f->ReleaseParameter(11);
  f->SetParLimits(11,-1,1);
  h->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L q","",minhisto,maxhisto);
  h->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L q","",minhisto,maxhisto);
  f->FixParameter(11,f->GetParameter(11)-f->GetParError(11));
  h->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L q","",minhisto,maxhisto);
  h->Fit(Form("f_%.0f_%.0f",ptmin,ptmax),"L m","",minhisto,maxhisto);
  
  TF1* background = new TF1(Form("background_%.0f_%.0f",ptmin,ptmax),"[0]+[1]*x+[2]*x*x+[3]*x*x*x");
  background->SetParameter(0,f->GetParameter(3));
  background->SetParameter(1,f->GetParameter(4));
  background->SetParameter(2,f->GetParameter(5));
  background->SetParameter(3,f->GetParameter(6));
  background->SetLineColor(4);
  background->SetRange(minhisto,maxhisto);
  background->SetLineStyle(2);
  
  TF1* mass = new TF1(Form("fmass_%.0f_%.0f",ptmin,ptmax),"[0]*([3]*([4]*Gaus(x,[1],[2]*(1+[6]))/(sqrt(2*3.14159)*[2]*(1+[6]))+(1-[4])*Gaus(x,[1],[5]*(1+[6]))/(sqrt(2*3.14159)*[5]*(1+[6]))))");
  mass->SetParameters(f->GetParameter(0),f->GetParameter(1),f->GetParameter(2),f->GetParameter(7),f->GetParameter(9),f->GetParameter(10),f->GetParameter(11));
  mass->SetParError(0,f->GetParError(0));
  mass->SetParError(1,f->GetParError(1));
  mass->SetParError(2,f->GetParError(2));
  mass->SetParError(3,f->GetParError(7));
  mass->SetParError(4,f->GetParError(9));
  mass->SetParError(5,f->GetParError(10));
  mass->SetFillColor(kOrange-3);
  mass->SetFillStyle(3002);
  mass->SetLineColor(kOrange-3);
  mass->SetLineWidth(3);
  mass->SetLineStyle(2);
  
  TF1* massSwap = new TF1(Form("fmassSwap_%.0f_%.0f",ptmin,ptmax),"[0]*(1-[2])*Gaus(x,[1],[3]*(1+[4]))/(sqrt(2*3.14159)*[3]*(1+[4]))");
  massSwap->SetParameters(f->GetParameter(0),f->GetParameter(1),f->GetParameter(7),f->GetParameter(8),f->GetParameter(11));
  massSwap->SetParError(0,f->GetParError(0));
  massSwap->SetParError(1,f->GetParError(1));
//.........这里部分代码省略.........

void draw_obj(std::string var="ptTau1", std::string cut="(ptTau2>45 && abs(etaTau1)<4.0 && abs(etaTau2)<4.0)", int chan_cat=1, int yTau=0, int yB=1, int nbins=10, double xmin=0, double xmax=500,std::string xtitle="leading #tau p_{T}", std::string ytitle="Events")
{
  SetStyle(); gStyle->SetLineStyleString(11,"20 10");
  TH1::SetDefaultSumw2(1);

  std::string dir = "/afs/cern.ch/work/j/jlawhorn/public/comb_ntuples/";
  double sigscale = 1000;
  double sigscale1 = 10; 
  std::stringstream scale; scale << sigscale;
  std::stringstream scale1; scale1 << sigscale1;

  //Cut definitions
  double luminosity = 3000;
  std::stringstream lumi; lumi << luminosity;

  std::string ttcut = "(tauCat1==1 && tauCat2==1 && bTag1>0 && bTag2>0 && ptTau1>45 && ptTau2>45)";
  std::string mtcut = "((tauCat1==1 && tauCat2==3)||(tauCat2==1 && tauCat1==3))*(bTag1>0 && bTag2>0 && ptTau1>30 && ptTau2>30)";
  std::string etcut = "((tauCat1==1 && tauCat2==2)||(tauCat2==1 && tauCat1==2))*(bTag1>0 && bTag2>0 && ptTau1>30 && ptTau2>30)";
  std::string emcut = "((tauCat1==3 && tauCat2==2)||(tauCat2==3 && tauCat1==2))*(bTag1>0 && bTag2>0 && ptTau1>30 && ptTau2>30)";

  std::string taucut = "(abs(etaTau1)<4.0 && abs(etaTau2)<4.0)";
  std::string bcut = "(ptB1>30 && ptB2>30 && abs(etaB1)<4.0 && abs(etaB2)<4.0)";

  std::string optcut = cut;

  if (chan_cat==0) { //tt
    optcut = cut+"*"+ttcut;
    if (yTau==1) optcut+="*"+taucut+"*(ptTau1>45 && ptTau2>45)";
  }
  else if (chan_cat==1) { //mt
    optcut = cut+"*"+mtcut;
    if (yTau==1) optcut+="*"+taucut+"*(ptTau1>30 && ptTau2>30)";
  }
  else if (chan_cat==2) { //et
    optcut = cut+"*"+etcut;
    if (yTau==1) optcut+="*"+taucut+"*(ptTau1>30 && ptTau2>30)";
  }
  else if (chan_cat==3) { //em
    optcut = cut+"*"+emcut;
    if (yTau==1) optcut+="*"+taucut+"*(ptTau1>30 && ptTau2>30)";
  }

  if (yB==1) optcut+="*"+bcut;

  //cout << optcut << endl;
  
  //signal region
  std::string hhcut = optcut+"*"+scale.str()+"*eventWeight*(1)*"+lumi.str();
  std::string ttbarcut = optcut+"*eventWeight*(1)*"+lumi.str();
  std::string hcut = optcut+"*eventWeight*(1)*"+lumi.str();
  std::string zjetcut = optcut+"*eventWeight*(eventType==4)*"+lumi.str();
  std::string ewkcut = optcut+"*eventWeight*(eventType!=4)*"+lumi.str();
  
  //--------------------------------------------------------------------------
  
  //Get the trees
  TTree *ttbartree = load(dir+"tt.root"); 
  TTree *htree = load(dir+"H.root"); 
  TTree *ewktree = load(dir+"EWK.root"); 
  TTree *sigtree = load(dir+"HHToTTBB_14TeV.root"); 

  //-------------------------------------------------------------------------
  
  //Get histograms
  TCanvas *canv0 = MakeCanvas("canv", "histograms", 600, 600);
  canv0->cd();
  std::string vardraw;
  TH1F *sig = new TH1F("HH","",nbins,xmin,xmax);
  vardraw = var+">>"+"HH";
  sigtree->Draw(vardraw.c_str(),hhcut.c_str());
  InitSignal(sig);
  sig->SetLineColor(kBlack);
  TH1F *ttbar = new TH1F("TTbar","",nbins,xmin,xmax);
  vardraw = var+">>"+"TTbar";
  ttbartree->Draw(vardraw.c_str(),ttbarcut.c_str());
  InitHist(ttbar, xtitle.c_str(), ytitle.c_str(), TColor::GetColor(155,152,204), 1001);
  TH1F *hbg = new TH1F("H","",nbins,xmin,xmax);
  vardraw = var+">>"+"H";
  htree->Draw(vardraw.c_str(),hcut.c_str());
  InitHist(hbg, xtitle.c_str(), ytitle.c_str(), TColor::GetColor(141,201,159), 1001);
  TH1F *zjet = new TH1F("Zjets","",nbins,xmin,xmax);
  vardraw = var+">>"+"Zjets";
  ewktree->Draw(vardraw.c_str(),zjetcut.c_str());
  InitHist(zjet, xtitle.c_str(), ytitle.c_str(),  TColor::GetColor(222,90,106), 1001);
  TH1F *ewk = new TH1F("EWK","",nbins,xmin,xmax);
  vardraw = var+">>"+"EWK";
  ewktree->Draw(vardraw.c_str(),ewkcut.c_str());
  InitHist(ewk, xtitle.c_str(), ytitle.c_str(),  TColor::GetColor(248,206,104), 1001);

  cout << sig->GetEntries() << endl;
  cout << ttbar->GetEntries() << endl;
  cout << hbg->GetEntries() << endl;
  cout << zjet->GetEntries() << endl;
  cout << ewk->GetEntries() << endl;
  
  delete canv0;
  //-----------------------------------------------------------------------
  //Draw the histograms
  TCanvas *canv = MakeCanvas("canv", "histograms", 600, 600);
  canv->cd();
//.........这里部分代码省略.........

//.........这里部分代码省略.........
    Dibosons->Add(WJets);
  }
  QCD->Add(Dibosons);
  TTJ->Add(QCD);
  ZTT->Add(TTJ);
  ZMM->Add(ZTT);
  if(log){
#ifdef MSSM
    ggH  ->Add(bbH);
#else
#ifndef DROP_SIGNAL
    qqH  ->Add(VH );
    ggH  ->Add(qqH);
#endif
#endif
  }
  else{
#ifdef MSSM
    bbH  ->Add(ZMM);
    ggH  ->Add(bbH);
#else
#ifndef DROP_SIGNAL
    VH   ->Add(ZMM);
    qqH  ->Add(VH );
    ggH  ->Add(qqH);
#endif
#endif
  }


  /*
    mass plot before and after fit
  */
  TCanvas* canv = MakeCanvas("canv", "histograms", 600, 600);
  canv->cd();
  if(log){ canv->SetLogy(1); }
#if defined MSSM
  if(!log){ data->GetXaxis()->SetRange(0, data->FindBin(345)); } else{ data->GetXaxis()->SetRange(0, data->FindBin(UPPER_EDGE)); };
#else
  data->GetXaxis()->SetRange(0, data->FindBin(345));
#endif
  data->SetNdivisions(505);
  data->SetMinimum(min);
#ifndef DROP_SIGNAL
  data->SetMaximum(max>0 ? max : std::max(std::max(maximum(data, log), maximum(ZMM, log)), maximum(ggH, log)));
#else
  data->SetMaximum(max>0 ? max : std::max(maximum(data, log), maximum(ZMM, log)));
#endif
  data->Draw("e");

  TH1F* errorBand = (TH1F*)ZMM->Clone("errorBand");
  errorBand->SetMarkerSize(0);
  errorBand->SetFillColor(1);
  errorBand->SetFillStyle(3013);
  errorBand->SetLineWidth(1);
  for(int idx=0; idx<errorBand->GetNbinsX(); ++idx){
    if(errorBand->GetBinContent(idx)>0){
      std::cout << "Uncertainties on summed background samples: " << errorBand->GetBinError(idx)/errorBand->GetBinContent(idx) << std::endl;
      break;
    }
  }
  if(log){
    ZMM->Draw("histsame");
    ZTT->Draw("histsame");
    TTJ->Draw("histsame");
    QCD->Draw("histsame");

void MultiHistoOverlap(TString namesandlabels, Int_t nOfFiles, const TString& outDir="./") {

    gROOT->Reset();
    gROOT->ProcessLine(".L tdrstyle.C");
    gROOT->ProcessLine("setTDRStyle()");

//  gSystem->Load("libRooFit");
//  using namespace RooFit;
// preamble
    TPaveText *cmsprel = new TPaveText(0.19, 0.95, 0.95, 0.99, "NDC");
    cmsprel->SetTextSize(0.03);
    cmsprel->SetTextFont(42);
    cmsprel->SetFillColor(0);
    cmsprel->SetBorderSize(0);
    cmsprel->SetMargin(0.01);
    cmsprel->SetTextAlign(12); // align left
    TString text = "CMS Preliminary 2011";
    cmsprel->AddText(0.0, 0.5,text);
    TString text2 = "#sqrt{s} = 7 TeV  |#eta_{#mu}|<2.4";
    cmsprel->AddText(0.8, 0.5, text2);


    TList* FileList  = new TList();
    TList* LabelList = new TList();
    TObjArray *nameandlabelpairs = namesandlabels.Tokenize(",");
    for (Int_t i = 0; i < nameandlabelpairs->GetEntries(); ++i) {
        TObjArray *aFileLegPair = TString(nameandlabelpairs->At(i)->GetName()).Tokenize("=");
        if(aFileLegPair->GetEntries() == 2) {
            FileList->Add( TFile::Open(aFileLegPair->At(0)->GetName())  );
            LabelList->Add( aFileLegPair->At(1) );
        } else {
            std::cout << "Please give file name and legend entry in the following form:\n" 		<< " filename1=legendentry1,filename2=legendentry2\n";
        }
    }


    Int_t NOfFiles =  FileList->GetSize();
    if ( NOfFiles!=nOfFiles ) {
        std::cout<<"&MSG-e: NOfFiles = "<<nOfFiles<<std::endl;
        return;
    }


    TString LegLabels[nOfFiles];
    for(Int_t j=0; j < nOfFiles; j++) {
        TObjString* legend = (TObjString*)LabelList->At(j);
        LegLabels[j] = legend->String();
        std::cout<<"LegLabels["<<j<<"]"<<LegLabels[j]<<std::endl;
    }

    TLegend *leg=0;

    TCanvas* c0 = new TCanvas("c0", "c0",50, 20, 800,600);
    TCanvas* c1 = new TCanvas("c1", "c1",50, 20, 800,600);
    TCanvas* c2 = new TCanvas("c2", "c2",50, 20, 800,600);
    TCanvas* c3 = new TCanvas("c3", "c3",50, 20, 800,600);
    TCanvas* c4 = new TCanvas("c4", "c4",50, 20, 800,600);
    TCanvas* c5 = new TCanvas("c5", "c5",50, 20, 1200,800);
    TCanvas* c6 = new TCanvas("c6", "c6",50, 20, 1200,800);

    TCanvas* c0s = new TCanvas("c0s", "c0s",50, 20, 800,600);
    TCanvas* c1s = new TCanvas("c1s", "c1s",50, 20, 800,600);
    TCanvas* c2s = new TCanvas("c2s", "c2s",50, 20, 800,600);
    TCanvas* c3s = new TCanvas("c3s", "c3s",50, 20, 800,600);

    TCanvas* cFit = new TCanvas("cFit", "cFit",50, 20, 1600, 800);


//----------------- CANVAS C0 --------------//
    c0->SetFillColor(0);
    c0->cd();

    leg = new TLegend(0.50,0.25,0.90,0.40);
    leg->SetBorderSize(1);
    leg->SetFillColor(0);
    leg->SetTextFont(42);

// Mass VS muon phi plus -------------------------------
    TH1D *histoMassVsPhiPlus[nOfFiles];
    for(Int_t j=0; j < nOfFiles; j++) {

        TFile *fin = (TFile*)FileList->At(j);
        if ( histoMassVsPhiPlus[j] = (TH1D*)fin->Get("MassVsPhiPlus/allHistos/meanHisto")) {
            histoMassVsPhiPlus[j]->SetLineStyle(linestylelist[j]);
            histoMassVsPhiPlus[j]->SetMarkerColor(colorlist[j]);
            histoMassVsPhiPlus[j]->SetLineColor(colorlist[j]);
            histoMassVsPhiPlus[j]->SetMarkerStyle(markerstylelist[j]);
            //     histoMassVsPhiPlus[j]->SetMarkerSize(0.75);
            if ( j == 0 ) {
                histoMassVsPhiPlus[j]->GetXaxis()->SetTitle("positive muon #phi (rad)");
                histoMassVsPhiPlus[j]->GetYaxis()->SetTitle("M_{#mu#mu} (GeV)");
                //       histoMassVsPhiPlus[j]->GetYaxis()->SetRangeUser(88.5,93.5);
                histoMassVsPhiPlus[j]->GetYaxis()->SetRangeUser(90.0,91.5);
                histoMassVsPhiPlus[j]->GetXaxis()->SetRangeUser(-3.14,3.14);
                histoMassVsPhiPlus[j]->Draw();
            } else {
                histoMassVsPhiPlus[j]->Draw("SAME");
            }
            leg->AddEntry(histoMassVsPhiPlus[j],LegLabels[j],"PL");
        }
//.........这里部分代码省略.........

//---------------------------------------------------------------
double *unbinnedFit(TString fiName, vector<double> xlim, double mtop, TLatex *channel_tex, TString outdir, double lumi)
//---------------------------------------------------------------
{
  using namespace RooFit;
  // reduce RooFit's verbosity on the INFO level
  RooMsgService::instance().getStream(1).removeTopic(Minimization);
  RooMsgService::instance().getStream(1).removeTopic(Plotting);
  RooMsgService::instance().getStream(1).removeTopic(ObjectHandling);
  RooMsgService::instance().getStream(1).removeTopic(Eval);
  RooMsgService::instance().getStream(1).removeTopic(Fitting);
  RooMsgService::instance().setSilentMode(true);
  RooMsgService::instance().setGlobalKillBelow(RooFit::FATAL);

  TFile *res = TFile::Open(fiName);

  RooRealVar mtl("mass", "M_{J/#psi+l}", 0., 250., "GeV");
  RooRealVar weight("weight", "weight", 0., 2.);
  //RooRealVar mean("mean", "mass", (xlim[1]+xlim[0])/2., (xlim[1]+3.*xlim[0])/4., (3.*xlim[1]+xlim[0])/4.);
  RooRealVar mean("mean", "mass", 70., 60., 80.);
  RooRealVar width("width", "width", 25., 15., 40.);

  TTree *tree = (TTree*)res->Get("MTriLept");
  RooDataSet *dataset = new RooDataSet("dataset", "dataset", RooArgSet(mtl, weight), Import(*tree), WeightVar(weight));
  RooGaussian pdf("gaus", "gaus", mtl, mean, width);
  pdf.fitTo(*dataset, Range(xlim[0], xlim[1]), SumW2Error(kTRUE), PrintLevel(-1), PrintEvalErrors(-1));

  TCanvas *cn = new TCanvas("cn", "cn", 800, 800);
  cn->cd();
  RooPlot *massframe = mtl.frame();
  if (mtop < 1e-6)
    dataset->plotOn(massframe, Binning(25), DataError(RooAbsData::SumW2));
  else
    dataset->plotOn(massframe, Binning(50), DataError(RooAbsData::SumW2));
  pdf.plotOn(massframe, Range(xlim[0], xlim[1]));
  massframe->Draw();
  TLegend *leg = new TLegend(0.58,0.82,0.93,0.92,NULL,"brNDC");
  if (mtop < 1e-6)
    leg->SetHeader(TString::Format("#tilde{M}_{J/#psi+l} = (%3.1f #pm %3.1f) GeV", mean.getVal(), mean.getError()));
  else
    leg->SetHeader(TString::Format("#splitline{M_{t}^{gen} = %3.1f GeV}{#tilde{M}_{J/#psi+l} = (%3.1f #pm %3.1f) GeV}", mtop+0.5, mean.getVal(), mean.getError()));
  leg_myStyle(leg);
  leg->Draw("same");
  channel_tex->Draw("same");
  if (mtop < 1e-6)
    cms_myStyle(lumi, true);
  else
    cms_myStyle(lumi, false);
  TString outFile = outdir;
  if (mtop < 1e-6)
    outFile += "GausianUnbinnedFit_Data";
  else
    outFile += TString::Format("GaussianUnbinnedFit_%d_5", (int)mtop);
  cn->SaveAs(outFile+".pdf");
  cn->SaveAs(outFile+".C");
  cn->SaveAs(outFile+".jpg");
  cn->SaveAs(outFile+".eps");

  res->Close();

  double *mean_err = new double[2];
  mean_err[0] = mean.getVal();
  mean_err[1] = mean.getError();

  return mean_err;
}

//---------------------------------------------------------------------------
void __fastcall TFMain_11011981::ShowUnits(bool showUnk)
{
    int oldItemIdx = lbUnits->ItemIndex;
    DWORD selAdr = 0;
    if (oldItemIdx != -1)
    {
        String item = lbUnits->Items->Strings[oldItemIdx];
        sscanf(item.c_str() + 1, "%lX", &selAdr);
    }
    int oldTopIdx = lbUnits->TopIndex;

    lbUnits->Clear();
    lbUnits->Items->BeginUpdate();

    if (UnitsNum)
    {
        switch (UnitSortField)
        {
        case 0:
            Units->Sort(SortUnitsByAdr);
            break;
        case 1:
            Units->Sort(SortUnitsByOrd);
            break;
        case 2:
            Units->Sort(SortUnitsByNam);
            break;
        }
    }

    int     stringLen;
    int     newItemIdx = -1;
    int     wid, maxwid = 0;
    TCanvas *canvas = lbUnits->Canvas;
    char    ci, cf, orderS[256];

    for (int i = 0; i < UnitsNum; i++)
    {
        PUnitRec recU = (UnitRec*)Units->Items[i];
        if (recU->fromAdr == selAdr) newItemIdx = i;
        ci = (!recU->trivialIni) ? 'I' : ' ';
        cf = (!recU->trivialFin) ? 'F' : ' ';
        stringLen = sprintf(StringBuf, " %08lX #%03d %c%c ", (int)recU->fromAdr, recU->iniOrder, ci, cf);
        if (recU->names->Count)
        {
            for (int u = 0; u < recU->names->Count; u++)
            {
                if (stringLen + recU->names->Strings[u].Length() >= 256)
                {
                    stringLen += sprintf(StringBuf + stringLen, "...");
                    break;
                }
                if (u)
                {
                    StringBuf[stringLen] = ';';
                    stringLen++;
                }
                stringLen += sprintf(StringBuf + stringLen, "%s", recU->names->Strings[u].c_str());
            }
        }
        else
            stringLen += sprintf(StringBuf + stringLen, "_Unit%d", recU->iniOrder);

        if (i != UnitsNum - 1)
        {
            //Trivial units
            if (recU->trivial)
                StringBuf[0] = TRIV_UNIT;
            else if (!recU->kb)
                StringBuf[0] = USER_UNIT;
        }
        //Last unit is user's
        else
        	StringBuf[0] = USER_UNIT;

        //Unit has undefined bytes
        if (showUnk && ContainsUnexplored(recU)) StringBuf[0] |= UNEXP_UNIT;
        String line = String(StringBuf, stringLen);
        lbUnits->Items->Add(line);
        wid = canvas->TextWidth(line); if (wid > maxwid) maxwid = wid;
    }
    if (newItemIdx == -1)
        lbUnits->TopIndex = oldTopIdx;
    else
    {
        if (newItemIdx != oldItemIdx)
        {
            lbUnits->ItemIndex = newItemIdx;
            int newTopIdx = newItemIdx - (oldItemIdx - oldTopIdx);
            if (newTopIdx < 0) newTopIdx = 0;
            lbUnits->TopIndex = newTopIdx;
        }
        else
        {
            lbUnits->ItemIndex = newItemIdx;
            lbUnits->TopIndex = oldTopIdx;
        }
    }
    lbUnits->ItemHeight = lbUnits->Canvas->TextHeight("T");
//.........这里部分代码省略.........

//---------------------------------------------------------------------------
void __fastcall TFMain_11011981::lbUnitItemsDrawItem(TWinControl *Control, int Index, TRect &Rect, TOwnerDrawState State)
{
    int         flags;
    TColor      _color;
    TListBox    *lb;
    TCanvas     *canvas;
    String      text, str;

    lb = (TListBox*)Control;
    canvas = lb->Canvas;
    SaveCanvas(canvas);

    if (Index < lb->Count)
    {
        flags = Control->DrawTextBiDiModeFlags(DT_SINGLELINE | DT_VCENTER | DT_NOPREFIX);
        if (!Control->UseRightToLeftAlignment())
            Rect.Left += 2;
        else
            Rect.Right -= 2;
        canvas->FillRect(Rect);

        text = lb->Items->Strings[Index];
        //lb->ItemHeight = canvas->TextHeight(text);
        str = text.SubString(2, text.Length() - 1);
        //Procs with Xrefs
        if (text[1] & 6)
        {
            //Xrefs from user units
            if (text[1] & 4)
            {
                if (!State.Contains(odSelected))
                    _color = TColor(0x00B000); //Green
                else
                    _color = TColor(0xBBBBBB); //LightGray
            }
            //No Xrefs from user units, only from KB units 
            else
            {
                if (!State.Contains(odSelected))
                    _color = TColor(0xC08000); //Blue
                else
                    _color = TColor(0xBBBBBB); //LightGray
            }
        }
        //Unresolved items
        else if (text[1] & 1)
        {
            if (!State.Contains(odSelected))
                _color = TColor(0x8080FF); //Red
            else
                _color = TColor(0xBBBBBB); //LightGray
        }
        //Other
        else
        {
            if (!State.Contains(odSelected))
                _color = TColor(0);        //Black
            else
                _color = TColor(0xBBBBBB); //LightGray
        }
        Rect.Right = Rect.Left;
        DrawOneItem(str, canvas, Rect, _color, flags);
    }
    RestoreCanvas(canvas);
}

//---------------------------------------------------------------------------
void __fastcall TFMain_11011981::lbUnitsDrawItem(
      TWinControl *Control, int Index, TRect &Rect, TOwnerDrawState State)
{
    char        *s, *pos;
    int         flags, len;
    TColor      _color;
    TListBox    *lb;
    TCanvas     *canvas;
    String      text, str1, str2;

    lb = (TListBox*)Control;
    canvas = lb->Canvas;
    SaveCanvas(canvas);

    if (Index < lb->Count)
    {
        flags = Control->DrawTextBiDiModeFlags(DT_SINGLELINE | DT_VCENTER | DT_NOPREFIX);
        if (!Control->UseRightToLeftAlignment())
            Rect.Left += 2;
        else
            Rect.Right -= 2;

        text = lb->Items->Strings[Index];
        //lb->ItemHeight = canvas->TextHeight(text);
        canvas->FillRect(Rect);

        s = text.c_str();
        //*XXXXXXXX #XXX XX NAME
        pos = strrchr(s, ' ');
        len = pos - s;
        str1 = text.SubString(2, len - 1);
        str2 = text.SubString(len + 1, text.Length() - len);

        if (!State.Contains(odSelected))
            _color = TColor(0);        //Black
        else
            _color = TColor(0xBBBBBB); //LightGray
        Rect.Right = Rect.Left;
        DrawOneItem(str1, canvas, Rect, _color, flags);

        //Unit name
        //Trivial unit - red
        if (text[1] & TRIV_UNIT)
        {
            if (!State.Contains(odSelected))
                _color = TColor(0x0000B0); //Red
            else
                _color = TColor(0xBBBBBB); //LightGray
        }
        else
        {
            //User unit - green
            if (text[1] & USER_UNIT)
            {
                if (!State.Contains(odSelected))
                {
                	if (text[1] & UNEXP_UNIT)
                    	_color = TColor(0xC0C0FF); //Light Red
                    else
                    	_color = TColor(0x00B000); //Green
                }
                else
                    _color = TColor(0xBBBBBB); //LightGray
            }
            //From knowledge base - blue
            else
            {
                if (!State.Contains(odSelected))
                    _color = TColor(0xC08000); //Blue
                else
                    _color = TColor(0xBBBBBB); //LightGray
            }
        }
        DrawOneItem(str2, canvas, Rect, _color, flags);
    }
    RestoreCanvas(canvas);
}

void DataRedux( double  db = -25., char * outfile = "output/sat_power.root", char * name_file = "SAT",  double  number_files = 5,   double xnbins = 1000,  double xbin_min = 0,  double xbin_max = 15,  double ynbins = 1000,  double ybin_min = -6,  double ybin_max = -8 )
//void DataRedux( double  db = -25., char * outfile = "output/sun_power.root", char * name_file = "SUN",  double  number_files = 5,  int xnbins = 1000,  double xbin_min = 0,  double xbin_max = 20.,  int ynbins = 50,  double ybin_min = -6.8,  double ybin_max = -7.4)
{

/**********************************************
* Global Variables
* ********************************************/

  const int numberfiles = number_files;
  const char namefile[1000] = name_file;
  const char output[1000] = outfile;
  
  char name[1000];
  char title[10000];
  
  TH2D* h[numberfiles];
  TH1D* h1D[numberfiles];
  TH1D* h1Dcut[numberfiles];
  
  ifstream f[numberfiles]; 
  
  TFile *output_root = new TFile(output,"RECREATE");


/***********************************************
 Opening ascii files and creating the root file 
************************************************/
 for(int i = 0; i < number_files; i++)
 {
   sprintf(name,"%s/%s%i.txt",namefile,name_file,i+1);
   f[i].open(name);
 }
 
  
/************************************************
  Filling histograms
*************************************************/

  for(int i = 0; i < numberfiles; i++)
  {
    double t, v, p, nlines=0, first;
    
    sprintf(name,"%s_%i", namefile, i+1);
    sprintf(title,"%s - Baseline %i", namefile, i+1);
    h[i] = new TH2D(name,title,xnbins, xbin_min, xbin_max, ynbins, ybin_min, ybin_max);

    sprintf(name,"%s_1D_%i", namefile, i+1);
    h1D[i] = new TH1D(name,title,xnbins, xbin_min, xbin_max);    
    
    sprintf(name,"%s_1D_cut_%i", namefile, i+1);
    h1Dcut[i] = new TH1D(name,title,xnbins, xbin_min, xbin_max);
    
    while (1) {
      f[i] >> t >> v ; 
      if (!f[i].good()) break;
      if(nlines==0) first = t; 
      p = VToP(v,db);
      h[i]->Fill(t-first,p);
      h1D[i]->Fill(t-first,p);
      Cut(h1Dcut[i],i,t-first,p);
      nlines++; 
    }
  }
  
 

/****************************************
 * Drawing and Saving 
 *****************************************/
  
  sprintf(name,"%s - P(t) vs t",namefile);
  TCanvas *c = new TCanvas("c",name,900,600);  
  
  TPad *pad1 = new TPad("pad1","",0.03,0.62,0.50,0.92,21);
  pad1->Draw();
  pad1->SetGridy();
  c->Clear();
  c->Divide(2,5);
 
  Style_Cosmetics();
  
  int canvas_number = 1;
  for(i=0; i<numberfiles; i++){
    c->SetGridy();
    // Plot Power vc t
    gStyle->SetOptLogy(1);
    c->cd(canvas_number);
      
      Style_Cosmetics();
    Power_Cosmetics(h1D[i], i);
    c->Update();

    // Plot FFT Magnitude
    gStyle->SetOptLogy(0);
    c->cd(canvas_number+1);
    TH1 *hm =0;
    FFT_Mag(h1Dcut[i],hm, i, namefile);
    c->Update();
    
    // Increment number camvas
//.........这里部分代码省略.........

void makePassFail_QCDMC()
{

	//jetbins
	//
	vector<pair<unsigned, unsigned> > jetBins;
	vector<pair<float, float> > htBins, mhtBins;
	pair<float, float> jetbin1(2,2);	
	pair<float, float> jetbin2(3,5);	
	pair<float, float> jetbin3(6,7);	
	pair<float, float> jetbin4(8,1000);	

	pair<float, float> htbin1(500,750);	
	pair<float, float> htbin2(750,1000);	
	pair<float, float> htbin3(1000,1250);	
	pair<float, float> htbin4(1250,1500);	
	pair<float, float> htbin5(1500,8000);	


	pair<float, float> mhtbin1(0,8000);	
	
	jetBins.push_back(jetbin1);
	jetBins.push_back(jetbin2);
	jetBins.push_back(jetbin3);
	jetBins.push_back(jetbin4);
	
	htBins.push_back(htbin1);
	htBins.push_back(htbin2);
	htBins.push_back(htbin3);
	htBins.push_back(htbin4);

	mhtBins.push_back(mhtbin1);

	vector<string> dphibins;
//	dphibins.push_back("0.15");
	dphibins.push_back("0.20");
//	dphibins.push_back("0.25");
//	dphibins.push_back("0.30");
//	dphibins.push_back("0.35");
//	dphibins.push_back("0.40");

	OpenFiles();
	TCanvas *c = new TCanvas("print");	
	c->Draw();
	c->Print("ratios.eps[");
	
/*	for (unsigned d = 0; d < dphibins.size(); ++d)
	{	
		for (unsigned jetbin = 0; jetbin < jetBins.size(); ++jetbin)
		{	
			for (unsigned htbin = 0; htbin < htBins.size(); ++htbin)
			{
				for (unsigned mhtbin = 0; mhtbin < mhtBins.size(); ++mhtbin)
				{
					stringstream htrange;
					htrange << htBins.at(htbin).first << "<HT<" << htBins.at(htbin).second << " GeV, "
						<< mhtBins.at(mhtbin).first << "<MHT<" << mhtBins.at(mhtbin).second << " GeV";

					stringstream folder;
					folder << "Hist/Njet" << jetBins.at(jetbin).first << "to" << jetBins.at(jetbin).second 
						<< "HT"   << htBins.at(htbin).first << "to" << htBins.at(htbin).second 
						<< "MHT"  << mhtBins.at(mhtbin).first << "to" << mhtBins.at(mhtbin).second;

					stringstream njetlabel, title, numeHistName, denoHistName, signalHistName, controlHistName;
					njetlabel << jetBins.at(jetbin).first << "-" << jetBins.at(jetbin).second << " Jets";
					title << "Njet " << njetlabel.str() << ", " << htrange.str() << ", #Delta #Phi _{min}<" << dphibins.at(d);
					numeHistName << folder.str() << "/smear_signal";
					//denoHistName << folder.str() << "/smeared_fail1";
					denoHistName << folder.str() << "/smeared_fail" << d;
					signalHistName << folder.str() << "/smear_signalFineBin";
					controlHistName << folder.str() << "/smeared_failFineBin" << d;
					//makePassFail_QCDMC(numeHistName.str(), denoHistName.str(), 
					//		title.str(), htrange.str(), signalHistName.str(), controlHistName.str()); 
					makePassFail_QCDMC(numeHistName.str(), denoHistName.str(), 
							title.str(), htrange, signalHistName.str(), controlHistName.str(), jetBins.at(jetbin), 50, 150); 
				}

			}
		}
	}
*/	


	for (unsigned d = 0; d < dphibins.size(); ++d)
	{	
		for (unsigned jetbin = 0; jetbin < jetBins.size(); ++jetbin)
		{	
			stringstream folder;
			folder << "Hist/Njet" << jetBins.at(jetbin).first << "to" << jetBins.at(jetbin).second; 
			string htrange("HT>500 GeV");
			cout << folder.str() << endl;

			stringstream njetlabel, title, numeHistName, denoHistName, signalHistName, controlHistName;
			njetlabel << jetBins.at(jetbin).first << "-" << jetBins.at(jetbin).second << " Jets";
			title << "Njet " << njetlabel.str() << ", " << htrange << ", #Delta #Phi _{min}<" << dphibins.at(d);
			numeHistName << folder.str() << "/smear_signal";
			denoHistName << folder.str() << "/smeared_fail" << d;
			signalHistName << folder.str() << "/smear_signalFineBin";
			//signalHistName << folder.str() << "/smear_signal";
			controlHistName << folder.str() << "/smeared_failFineBin" << d;
//.........这里部分代码省略.........

void fillTree(TTree*& tree, TGraph*& graph, double& limit, unsigned int itype, std::map<double, std::string>& tanb_values, bool upper_exclusion, unsigned int verbosity)
{
  double value=-99;
  double tanb_help=-99;
  unsigned int ibin=0;
  // fill graph with scanned points
  for(std::map<double, std::string>::const_iterator tanb = tanb_values.begin(); tanb!=tanb_values.end(); ++tanb){
    value = singlePointLimit(tanb->second, tanb->first, itype, verbosity);
    if( value>0 ){
      graph->SetPoint(ibin++, tanb->first, value); 
    }
    tanb_help=tanb->first;
  }
  // determine smooth curve on graph for interpolation
  TSpline3* spline = new TSpline3("spline", graph, "r", 3., 10.);
  // linear polarisation func
  TF1 *fnc = 0;
  // determine all crossing points with y==1 
  std::vector<CrossPoint> points = crossPoints(graph);

  int dist = 1;
  bool filled = false;
  unsigned int np = 0;
  unsigned int steps = 10e6; 
  if(points.size()>0) limit = graph->GetX()[upper_exclusion ? points.begin()->first : points.end()->first];

  for(std::vector<CrossPoint>::const_reverse_iterator point = points.rbegin(); point!=points.rend(); ++point, ++np){
  //for(std::vector<CrossPoint>::iterator point = points.begin(); point!=points.end(); ++point, ++np){
    //double min = (point->first-dist)>0 ? graph->GetX()[point->first-dist] : graph->GetX()[0]; 
    double min = (point->first)>0 ? graph->GetX()[point->first] : graph->GetX()[0]; 
    double max = (point->first+dist)<graph->GetN() ? graph->GetX()[point->first+dist] : graph->GetX()[graph->GetN()-1];

    //double y_min = (point->first-dist)>0 ? graph->GetY()[point->first-dist] : graph->GetY()[0]; 
    double y_min = (point->first)>0 ? graph->GetY()[point->first] : graph->GetY()[0];
    double y_max = (point->first+dist)<graph->GetN() ? graph->GetY()[point->first+dist] : graph->GetY()[graph->GetN()-1];

    vector<double> crossing;
    crossing.push_back((min-max-y_max*min+y_min*max)/(y_min-y_max));
    //double crossing;
    //crossing = (1.-y_min)/(y_max-y_min)*(max-min); 
				  
    double deltaM = -999.;
    double offset = min; double step_size = (max-min)/steps;
    for(unsigned int scan=0; scan<=steps; ++scan){
      if(deltaM<0 || fabs(spline->Eval(offset+scan*step_size)-1.)<deltaM){
	limit=offset+scan*step_size;
	deltaM=fabs(spline->Eval(offset+scan*step_size)-1.);
      }
    }
    std::cout << "****************************************************************" << std::endl;
    std::cout << "* [" << np+1 << "|" << point->second << "] asymptotic limit(";
    std::cout << limitType(itype) << ") :" << crossing[np] << " -- " << limit << " deltaM : " << deltaM;
//     if(((upper_exclusion && point->second) || (!upper_exclusion && !(point->second))) && !filled){
//       //std::cout << "limit is taken from linear interpolation at the moment" << std::endl;
//       //limit = crossing;
//       std::cout << "    [-->to file]"; filled=true; tree->Fill();
//     }
    if(np==0){
      fnc = new TF1("fnc", "[0]*x+[1]", min, max);
      fnc->SetParameter(0, (y_min-y_max)/(min-max));
      fnc->SetParameter(1, (y_max*min-y_min*max)/(min-max));
      std::cout << std::endl;
      std::cout << "limit is taken from linear interpolation at the moment" << std::endl;
      limit = crossing[np];
      std::cout << limit << std::endl;
      std::cout << "    [-->to file]"; filled=true; tree->Fill();
    }
    std::cout << endl;
    std::cout << "****************************************************************" << std::endl;
  }
  // catch cases where no crossing point was found
  if(!filled){
    if(value<1)
      {
	std::cout << "WARNING: no crossing found - all tanb values excluded: " << value << std::endl;
	if(itype == observed)     { limit=3.00; }
	if(itype == plus_2sigma)  { limit=5.00; }
	if(itype == plus_1sigma)  { limit=4.00; }
	if(itype == expected)     { limit=3.00; }
	if(itype == minus_1sigma) { limit=2.50; }
	if(itype == minus_2sigma) { limit=2.00; }
	tree->Fill();
      }
    else
      {
	std::cout << "WARNING: no crossing found - no tanb value excluded: " << value << " -- " << tanb_help << std::endl;
	if(itype == observed)     { limit=tanb_help*value; }
	if(itype == plus_2sigma)  { limit=tanb_help*value; }
	if(itype == plus_1sigma)  { limit=tanb_help*value; }
	if(itype == expected)     { limit=tanb_help*value; }
	if(itype == minus_1sigma) { limit=tanb_help*value; }
	if(itype == minus_2sigma) { limit=tanb_help*value; }
	tree->Fill();
      }
  }
  if( verbosity>0 ){
    std::string monitor = std::string("SCAN-")+limitType(itype);
    TCanvas* canv = new TCanvas(monitor.c_str(), monitor.c_str(), 600, 600);
    TH1F* frame = canv->DrawFrame(graph->GetX()[0]-0.1, 0., graph->GetX()[graph->GetN()-1]+0.1, 10.);
    canv->SetGridx(1); canv->SetGridy(1); canv->cd(); 
//.........这里部分代码省略.........

//.........这里部分代码省略.........
       else if(i==15) ngroup = 8;
     } else if( fTitle.find("2360") != string::npos ) {
       if(i==10) ngroup = 4;
       else if(i==11) ngroup = 4;
       else if(i==12) ngroup = 4;
       else if(i==13) ngroup = 4;
       else if(i==14) ngroup = 8;
       else if(i==15) ngroup = 8;
     }
     
     h_temp->Rebin(ngroup);
  
     TF1* f1= new TF1("f1", "gaus", h_temp->GetXaxis()->GetXmin(), h_temp->GetXaxis()->GetXmax());
     f1->SetParameters(h_temp->GetBinContent(h_temp->GetMaximumBin()),h_temp->GetMean(),h_temp->GetRMS());
     h_temp->Fit("f1","R", "",h_temp->GetXaxis()->GetXmin(),h_temp->GetXaxis()->GetXmax());
  
     h[1]->SetBinContent(i, f1->GetParameter(2));
     h[1]->SetBinError(i, f1->GetParError(2));
      
     delete f1;
     delete h_temp;
   }
   
   h[1]->SetTitle(fTitle.c_str());
   h[1]->GetXaxis()->SetTitle(fXAxisLabel.c_str());
   h[1]->GetYaxis()->SetTitle(fYAxisLabel.c_str());
   h[1]->GetXaxis()->SetRangeUser(fXmin,fXmax);
   h[1]->GetYaxis()->SetRangeUser(fYmin,fYmax);

   h[1]->SetTitleOffset(1.,"X");
   h[1]->SetTitleOffset(1.,"Y");
//    h[1]->GetXaxis()->SetTitleSize(0.04);
//    h[1]->GetYaxis()->SetTitleSize(0.04);
   
   TCanvas *c = new TCanvas("c","",1200,800);
//    TCanvas *c = new TCanvas("c");
   c->cd();

   TF1 *fit = new TF1("fit","sqrt(pow([0],2)+pow([1],2)*(x-[3])+pow([2]*(x-[3]),2))",fXmin,(fXmax+1));
   fit->SetParName(0, "A");
   fit->SetParameter(0, 0.5);
   fit->SetParLimits(0, 0.1, 2);
   fit->SetParName(1, "B");
   fit->SetParameter(1, 0.5);
   fit->SetParLimits(1, 0, 1.5);
   fit->SetParName(2, "C");
   fit->SetParameter(2, 0.02);
   fit->SetParLimits(2, 0, 0.2);
   fit->SetParName(3, "D");
   fit->SetParameter(3, 3);
   fit->SetParLimits(3, 2, 5);
//    fit->FixParameter(3, 0);

   if(doFit==1) h[1]->Fit("fit","RV");

   h[0]->SetLineWidth(2);
   h[0]->SetLineStyle(2);
   h[0]->SetLineColor(kRed);
   h[0]->SetFillColor(kRed);
//    h[0]->SetMarkerSize(.6);
   h[0]->SetMarkerStyle(26);
   h[0]->SetMarkerColor(kRed);
   h[1]->SetLineWidth(2);
//    h[1]->SetLineStyle(3);
   h[1]->SetLineColor(kBlack);
//    h[1]->SetMarkerSize(.8);
   h[1]->SetMarkerStyle(20);
   h[1]->SetMarkerColor(kBlack);
   h[1]->Draw("");
   h[0]->Draw("same");
   if(doFit==1) fit->Draw("same");
   
   TLegend *legend = new TLegend(.65,.25,.95,.4);
   legend->SetBorderSize(0);
   legend->SetFillColor(0);
   legend->SetFillStyle(0);
   legend->SetTextFont(62);
   legend->AddEntry(h[1],"Data","lp");
   legend->AddEntry(h[0],"Simulation","lp");
   legend->Draw();
   
   TLatex l;
   l.SetTextAlign(12);
   l.SetTextSize(0.05);
   l.SetTextFont(62);
   l.SetNDC();
   l.DrawLatex(0.14,0.84,"CMS 2009 Preliminary");
   l.DrawLatex(0.14,0.79,("#sqrt{s}="+fTitle).c_str());
   l.DrawLatex(0.12,0.5,"f(x)=#sqrt{A^{2}+B^{2}(x-D)+C^{2}(x-D)^{2}}");

   c->SetGridx();
   c->SetGridy();

   c->SaveAs(fName.c_str());
   
   delete legend;
   delete h[0];
   delete h[1];
   delete c;
}

void showFunctions1d(std::vector<TF1*>& functions, const std::vector<std::string>& legendEntries,
		     const TString& xAxisTitle, double yMin, double yMax, const TString& yAxisTitle, 
		     const std::string& outputFileName)
{
  assert(functions.size() == legendEntries.size());

  TCanvas* canvas = new TCanvas("canvas", "canvas", 800, 600);
  canvas->SetFillColor(10);
  canvas->SetBorderSize(2);
  
  canvas->SetLeftMargin(0.12);
  canvas->SetBottomMargin(0.12);

  assert(functions.size() > 0);
  TF1* refFunction = functions[0];

  TH1* dummyHistogram = new TH1F("dummyHistogram", "dummyHistogram", 10, refFunction->GetXmin(), refFunction->GetXmax());
  dummyHistogram->SetStats(false);
  dummyHistogram->SetTitle("");
  dummyHistogram->SetMinimum(-TMath::Log(yMax));
  dummyHistogram->SetMaximum(-TMath::Log(yMax) + 5.);
  
  TAxis* xAxis = dummyHistogram->GetXaxis();
  xAxis->SetTitle(xAxisTitle.Data());
  xAxis->SetTitleOffset(1.15);

  TAxis* yAxis = dummyHistogram->GetYaxis();
  yAxis->SetTitle(Form("-log(%s)", yAxisTitle.Data()));
  yAxis->SetTitleOffset(1.30);

  dummyHistogram->Draw("axis");

  int colors[] = { 1, 2, 3, 4, 5, 6, 7, 15 };
  int numFunctions = functions.size();
  if ( numFunctions > 8 ) {
    std::cerr << "<showFunctions1d>:" << std::endl;
    std::cerr << "Number of functions must not exceed 8 !!" << std::endl;
    assert(0);
  }
  
  for ( int iFunction = 0; iFunction < numFunctions; ++iFunction ) {
    TF1* function = functions[iFunction];
    function->SetLineColor(colors[iFunction]);
    function->SetLineWidth(2);
    function->Draw("same");
  }

  TLegend* legend = new TLegend(0.68, 0.89 - (0.03 + 0.040*numFunctions), 0.89, 0.89, "", "brNDC"); 
  legend->SetBorderSize(0);
  legend->SetFillColor(0);
  legend->SetTextSize(0.035);
  for ( int iFunction = 0; iFunction < numFunctions; ++iFunction ) {
    TF1* function = functions[iFunction];
    const std::string& legendEntry = legendEntries[iFunction];
    legend->AddEntry(function, legendEntry.data(), "l");
  }
  legend->Draw();

  canvas->Update();
  size_t idx = outputFileName.find_last_of('.');
  std::string outputFileName_plot = std::string(outputFileName, 0, idx);
  if ( idx != std::string::npos ) canvas->Print(std::string(outputFileName_plot).append(std::string(outputFileName, idx)).data());
  canvas->Print(std::string(outputFileName_plot).append(".png").data());
  canvas->Print(std::string(outputFileName_plot).append(".pdf").data());

  delete dummyHistogram;
  delete legend;
  delete canvas;
}