C++ RooRealVar::getErrorLo方法代码示例

TGraphAsymmErrors *plotEffPt(RooDataSet *a, int aa) {
    const RooArgSet *set = a->get();
    RooRealVar *xAx = (RooRealVar*)set->find("pt");
    RooRealVar *eff = (RooRealVar*)set->find("efficiency");

    const int nbins = xAx->getBinning().numBins();

    double tx[nbins], txhi[nbins], txlo[nbins];
    double ty[nbins], tyhi[nbins], tylo[nbins];

    for (int i=0; i<nbins; i++) {
        a->get(i);
        ty[i] = eff->getVal();
        tx[i] = xAx->getVal();
        txhi[i] = fabs(xAx->getErrorHi());
        txlo[i] = fabs(xAx->getErrorLo());
        tyhi[i] = fabs(eff->getErrorHi());
        tylo[i] = fabs(eff->getErrorLo());
    }

    cout<<"NBins : "<<nbins<<endl;

    const double *x = tx;
    const double *xhi = txhi;
    const double *xlo = txlo;
    const double *y = ty;
    const double *yhi = tyhi;
    const double *ylo = tylo;

    TGraphAsymmErrors *b = new TGraphAsymmErrors();
    if(aa == 1) {
        *b = TGraphAsymmErrors(nbins,x,y,xlo,xhi,ylo,yhi);
    }
    if(aa == 0) {
        *b = TGraphAsymmErrors(nbins,x,y,0,0,ylo,yhi);
    }
    b->SetMaximum(1.1);
    b->SetMinimum(0.0);
    b->SetMarkerStyle(20);
    b->SetMarkerColor(kRed+2);
    b->SetMarkerSize(1.0);
    b->SetTitle("");
    b->GetXaxis()->SetTitleSize(0.05);
    b->GetYaxis()->SetTitleSize(0.05);
    b->GetXaxis()->SetTitle("p_{T} [GeV/c]");
    b->GetYaxis()->SetTitle("Efficiency");
    b->GetXaxis()->CenterTitle();
    //b->Draw("apz");

    for (int i=0; i<nbins; i++) {
        cout << x[i] << " " << y[i] << " " << yhi[i] << " " << ylo[i] << endl;
    }

    return b;

}

pair<double,double> bkgEvPerGeV(RooWorkspace *work, int m_hyp, int cat, int spin=false){
  
  RooRealVar *mass = (RooRealVar*)work->var("CMS_hgg_mass");
  if (spin) mass = (RooRealVar*)work->var("mass");
  mass->setRange(100,180);
  RooAbsPdf *pdf = (RooAbsPdf*)work->pdf(Form("pdf_data_pol_model_8TeV_cat%d",cat));
  RooAbsData *data = (RooDataSet*)work->data(Form("data_mass_cat%d",cat));
  RooPlot *tempFrame = mass->frame();
  data->plotOn(tempFrame,Binning(80));
  pdf->plotOn(tempFrame);
  RooCurve *curve = (RooCurve*)tempFrame->getObject(tempFrame->numItems()-1);
  double nombkg = curve->Eval(double(m_hyp));
 
  RooRealVar *nlim = new RooRealVar(Form("nlim%d",cat),"",0.,0.,1.e5);
  //double lowedge = tempFrame->GetXaxis()->GetBinLowEdge(FindBin(double(m_hyp)));
  //double upedge  = tempFrame->GetXaxis()->GetBinUpEdge(FindBin(double(m_hyp)));
  //double center  = tempFrame->GetXaxis()->GetBinUpCenter(FindBin(double(m_hyp)));

  nlim->setVal(nombkg);
  mass->setRange("errRange",m_hyp-0.5,m_hyp+0.5);
  RooAbsPdf *epdf = 0;
  epdf = new RooExtendPdf("epdf","",*pdf,*nlim,"errRange");
		
  RooAbsReal *nll = epdf->createNLL(*data,Extended(),NumCPU(4));
  RooMinimizer minim(*nll);
  minim.setStrategy(0);
  minim.setPrintLevel(-1);
  minim.migrad();
  minim.minos(*nlim);
  
  double error = (nlim->getErrorLo(),nlim->getErrorHi())/2.;
  data->Print(); 
  return pair<double,double>(nombkg,error); 
}

TH2F *plotEff2D(RooDataSet *a, TString b) {
    const RooArgSet *set = a->get();
    RooRealVar *yAx = (RooRealVar*)set->find("pt");
    RooRealVar *xAx = (RooRealVar*)set->find("eta");
    RooRealVar *eff = (RooRealVar*)set->find("efficiency");

//    const int xnbins = xAx->getBinning().numBins();
//    const int ynbins = yAx->getBinning().numBins();

    //double xbins[] = {-2.4, -1.6, -0.8, 0.0, 0.8, 1.6, 2.4};
    //double ybins[] = {0, 2, 3, 5, 8, 10, 20};

    const double *xvbins = xAx->getBinning().array();
    const double *yvbins = yAx->getBinning().array();

    TH2F* h = new TH2F(b, "", xAx->getBinning().numBins(), xvbins, yAx->getBinning().numBins(), yvbins);

    gStyle->SetPaintTextFormat("5.2f");
    gStyle->SetPadRightMargin(0.12);
    gStyle->SetPalette(1);
    h->SetOption("colztexte");
    h->GetZaxis()->SetRangeUser(-0.001,1.001);
    h->SetStats(kFALSE);
    h->GetYaxis()->SetTitle("p_{T} [GeV/c]");
    h->GetXaxis()->SetTitle("#eta");
    h->GetXaxis()->CenterTitle();
    h->GetYaxis()->CenterTitle();
    h->GetXaxis()->SetTitleSize(0.05);
    h->GetYaxis()->SetTitleSize(0.05);
    h->GetYaxis()->SetTitleOffset(0.8);
    h->GetXaxis()->SetTitleOffset(0.9);
    for(int i=0; i<a->numEntries(); i++) {
        a->get(i);
        h->SetBinContent(h->FindBin(xAx->getVal(), yAx->getVal()), eff->getVal());
        h->SetBinError(h->FindBin(xAx->getVal(), yAx->getVal()), (eff->getErrorHi()-eff->getErrorLo())/2.);
    }

    return h;

}

//.........这里部分代码省略.........
  RooArgSet DataFitParam(ResDataCombinedFit->floatParsFinal());
  RooRealVar* DataEffFit     = (RooRealVar*)(&DataFitParam["DataEfficiency"]);
  RooRealVar* DataNumSigFit  = (RooRealVar*)(&DataFitParam["DataNumSgn"]);
  RooRealVar* DataNumBkgPFit = (RooRealVar*)(&DataFitParam["DataNumBkgP"]);
  RooRealVar* DataNumBkgFFit = (RooRealVar*)(&DataFitParam["DataNumBkgF"]);

  RooPlot* DataFrameP = mass->frame(Bins(24),Title("Data: passing sample"));
  DataCombData.plotOn(DataFrameP,Cut("DataPassing==DataPassing::pass"));
  w->pdf("DataModel")->plotOn(DataFrameP,Slice(*(w->cat("DataPassing")),"pass"), ProjWData(*(w->cat("DataPassing")),DataCombData), LineColor(kBlue),Range(xLow_,xHigh_));
  w->pdf("DataModel")->plotOn(DataFrameP,Slice(*(w->cat("DataPassing")),"pass"), ProjWData(*(w->cat("DataPassing")),DataCombData), Components("DataSignalPdfP"), LineColor(kRed),Range(xLow_,xHigh_));
  w->pdf("DataModel")->plotOn(DataFrameP,Slice(*(w->cat("DataPassing")),"pass"), ProjWData(*(w->cat("DataPassing")),DataCombData), Components("DataBackgroundPdfP"), LineColor(kGreen),LineStyle(kDashed),Range(xLow_,xHigh_));
  
  RooPlot* DataFrameF = mass->frame(Bins(24),Title("Data: failing sample"));
  DataCombData.plotOn(DataFrameF,Cut("DataPassing==DataPassing::fail"));
  w->pdf("DataModel")->plotOn(DataFrameF,Slice(*(w->cat("DataPassing")),"fail"), ProjWData(*(w->cat("DataPassing")),DataCombData), LineColor(kBlue),Range(xLow_,xHigh_));
  w->pdf("DataModel")->plotOn(DataFrameF,Slice(*(w->cat("DataPassing")),"fail"), ProjWData(*(w->cat("DataPassing")),DataCombData), Components("DataSignalPdfF"), LineColor(kRed),Range(xLow_,xHigh_));
  w->pdf("DataModel")->plotOn(DataFrameF,Slice(*(w->cat("DataPassing")),"fail"), ProjWData(*(w->cat("DataPassing")),DataCombData), Components("DataBackgroundPdfF"), LineColor(kGreen),LineStyle(kDashed),Range(xLow_,xHigh_));
  ///////////////////////////////////////////////////////////////

 
  if(makeSoupFit_) c->Divide(2,2);
  else c->Divide(2,1);
 
  c->cd(1);
  DataFrameP->Draw();
  c->cd(2);
  DataFrameF->Draw();

  if(makeSoupFit_){
    c->cd(3);
    McFrameP->Draw();
    c->cd(4);
    McFrameF->Draw();
  }
 
  c->Draw();
 
  test->cd(Form("bin%f",binCenter_));
 
  c->Write();
 
  c2->Divide(2,1);
  c2->cd(1);
  TemplateFrameP->Draw();
  c2->cd(2);
  TemplateFrameF->Draw();
  c2->Draw();
 
  test->cd(Form("bin%f",binCenter_));
  c2->Write();


  // MINOS errors, otherwise HESSE quadratic errors
  float McErrorLo = 0;
  float McErrorHi = 0;
  if(makeSoupFit_){
    McErrorLo = McEffFit->getErrorLo()<0 ? McEffFit->getErrorLo() : (-1)*McEffFit->getError();
    McErrorHi = McEffFit->getErrorHi()>0 ? McEffFit->getErrorHi() : McEffFit->getError();
  }
  float DataErrorLo = DataEffFit->getErrorLo()<0 ? DataEffFit->getErrorLo() : (-1)*DataEffFit->getError();
  float DataErrorHi = DataEffFit->getErrorHi()>0 ? DataEffFit->getErrorHi() : DataEffFit->getError();
  float BinomialError = TMath::Sqrt(SGNtruePass/SGNtrue*(1-SGNtruePass/SGNtrue)/SGNtrue);
 
  Double_t* truthMC = new Double_t[6];
  Double_t* tnpMC   = new Double_t[6];
  Double_t* tnpData = new Double_t[6];

  truthMC[0] = binCenter_;
  truthMC[1] = binWidth_;
  truthMC[2] = binWidth_;
  truthMC[3] = SGNtruePass/SGNtrue;
  truthMC[4] = BinomialError;
  truthMC[5] = BinomialError;
  if(makeSoupFit_){
    tnpMC[0] = binCenter_;
    tnpMC[1] = binWidth_;
    tnpMC[2] = binWidth_;
    tnpMC[3] = McEffFit->getVal();
    tnpMC[4] = (-1)*McErrorLo;
    tnpMC[5] = McErrorHi;
  }
  tnpData[0] = binCenter_;
  tnpData[1] = binWidth_;
  tnpData[2] = binWidth_;
  tnpData[3] = DataEffFit->getVal();
  tnpData[4] = (-1)*DataErrorLo;
  tnpData[5] = DataErrorHi;

  out.push_back(truthMC);
  out.push_back(tnpData);
  if(makeSoupFit_) out.push_back(tnpMC);

  test->Close();

  //delete c; delete c2;

  if(verbose_) cout << "returning from bin " << bin_ << endl;
  return out;

}

//.........这里部分代码省略.........
  TString outputtext;//for useful text

  //plot things  
  RooPlot *framex = mass->frame();
  framex->GetYaxis()->SetTitle("Events/(5 MeV)");
  data->plotOn(framex,Name("Hist"),MarkerColor(kBlack),LineColor(kBlack),DataError(RooAbsData::SumW2));
  totalPdf.plotOn(framex,Name("curvetot"),LineColor(kBlue));
  RooArgSet* totalPdfComponents = totalPdf.getComponents();
  TIterator* itertPC = totalPdfComponents->createIterator();
  RooAddPdf* vartPC = (RooAddPdf*) itertPC->Next();
  vartPC = (RooAddPdf*) itertPC->Next();//skip totalPdf
  int i=0;//color index
  TLegend *leg = new TLegend(0.2, 0.02, .4, .42);  
  leg->SetTextSize(0.06);
  leg->AddEntry(framex->findObject("curvetot"),"Total PDF","l");
  while(vartPC){//loop over compotents of totalPdf
    TString vartPCtitle = vartPC->GetTitle();
    TIterator* itercompPars;//forward declare so it persists outside the if statement
    RooRealVar* varcompPars;
    if(!(vartPCtitle.Contains(":")||vartPCtitle.Contains("@"))){//only for non-sub-shapes
      while(i==0||i==10||i==4||i==1||i==5||(i>=10&&i<=27))i++;//avoid white and blue and black and yellow and horribleness
      RooArgSet* compPars = vartPC->getParameters(data);//set of the parameters of the component the loop is on
      itercompPars = compPars->createIterator();
      varcompPars = (RooRealVar*) itercompPars->Next();
    
      while(varcompPars){//write and print mean, sig, etc. of sub-shapes
        TString vartitle = varcompPars->GetTitle();
        double varval = varcompPars->getVal();
        TString varvalstring = Form("%f",varval);
        double hi = varcompPars->getErrorHi();
        
        TString varerrorstring = "[exact]";
        if(hi!=-1){
          double lo = varcompPars->getErrorLo();
          double varerror = TMath::Max(fabs(lo),hi);
          varerrorstring = Form("%E",varerror);
        }
        
        outputtext = vartitle+" = "+varvalstring+" +/- "+varerrorstring;
        textfile<<outputtext<<endl;
        cout<<outputtext<<endl;
        
        varcompPars = (RooRealVar*) itercompPars->Next(); 
      }
      totalPdf.plotOn(framex,Name(vartPC->GetName()),LineStyle(kDashed),LineColor(i),Components(vartPC->GetName()));
      leg->AddEntry(framex->findObject(vartPC->GetName()),vartPCtitle,"l");
    
      i++;
    }
    vartPC = (RooAddPdf*) itertPC->Next();
    itercompPars->Reset();//make sure it's ready for the next vartPC
  }
  
  // Calculate chi2/ndf
  RooArgSet *floatpar = totalPdf.getParameters(data);
  int floatpars = (floatpar->selectByAttrib("Constant",kFALSE))->getSize();
  Double_t chi2 = framex->chiSquare("curvetot","Hist",floatpars);
  TString chi2string = Form("%f",chi2);
  //create text box to list important parameters on the plot
  // TPaveText* txt = new TPaveText(0.1,0.5,0.7,0.9,"NBNDC");
  // txt->SetTextSize(0.06);
  // txt->SetTextColor(kBlack);
  // txt->SetBorderSize(0);
  // txt->SetFillColor(0);
  // txt->SetFillStyle(0);
  outputtext = "#chi^{2}/N_{DoF} = "+chi2string;

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

  //return;



  // binned combined dataset
  RooDataHist DataCombData("DataCombData","combined data",mass,Index(category),Import("pass", *(DataDataSetP.createHistogram("histoDataP",mass)) ) ,Import("fail", *(DataDataSetF.createHistogram("histoDataF",mass))), Weight(0.5) ) ;
  std::cout << "data dataHist Comb " << DataCombData.sumEntries() << "  " << std::endl;
  std::cout << "+++++++++++++++++++++++++++++++++++++++++++++" << std::endl;
  // unbinned combined dataset
  RooDataSet DataCombDataUnBinned("DataCombDataUnBinned","combined data",mass,Index(category),Import("pass", DataDataSetP ) ,Import("fail",DataDataSetF), Weight(0.5) ) ;
  std::cout << "data dataset Comb " << DataCombDataUnBinned.numEntries() << "  " << std::endl;
  std::cout << "+++++++++++++++++++++++++++++++++++++++++++++" << std::endl;
  //return out;
  
  RooSimultaneous DataSimPdf("DataSimPdf","simultaneous pdf",category) ;
  DataSimPdf.addPdf(DataModelP,"pass") ;
  DataSimPdf.addPdf(DataModelF,"fail") ;

  //mass.setBins( 10000, "fft" );
  mass.setBins( nBins_ );
  RooFitResult* ResDataCombinedFit =  0;
  if(doBinned_)  ResDataCombinedFit = DataSimPdf.fitTo(DataCombData , Extended(1), Minos(1), Save(1), NumCPU(4), /*ExternalConstraints( RooArgSet(alfaSgn_CPdf,nSgn_CPdf) )*/  SumW2Error(1));
  else ResDataCombinedFit = DataSimPdf.fitTo(DataCombDataUnBinned , Extended(1), Minos(1), Save(1), NumCPU(4),  /*ExternalConstraints( RooArgSet(alfaSgn_CPdf,nSgn_CPdf) )*/ SumW2Error(1));


  RooArgSet DataFitParam(ResDataCombinedFit->floatParsFinal());
  RooRealVar* DataEffFit      = (RooRealVar*)(&DataFitParam["DataEfficiency"]);
  RooRealVar* DataNumSigFit   = (RooRealVar*)(&DataFitParam["DataNumSgn"]);

  RooPlot* DataFrameP = mass.frame(Bins(40),Title("CMS Preliminary 2011  #sqrt{s}=7 TeV   L=XXX pb^{-1}:  passing probe"));
  DataCombData.plotOn(DataFrameP,Cut("category==category::pass"));
  DataSimPdf.plotOn(DataFrameP,Slice(category,"pass"), ProjWData(category,DataCombData), LineColor(kBlue));
  DataSimPdf.plotOn(DataFrameP,Slice(category,"pass"), ProjWData(category,DataCombData), Components("sgnPdfP"), LineColor(kRed), LineStyle(kSolid));
  DataSimPdf.plotOn(DataFrameP,Slice(category,"pass"), ProjWData(category,DataCombData), Components("bkgPdfP"), LineColor(kMagenta), LineStyle(kDashed));
 

  RooPlot* DataFrameF = mass.frame(Bins(40),Title("CMS Preliminary 2011  #sqrt{s}=7 TeV   L=XXX pb^{-1}:  failing probe"));
  DataCombData.plotOn(DataFrameF,Cut("category==category::fail"));
  DataSimPdf.plotOn(DataFrameF,Slice(category,"fail"), ProjWData(category,DataCombData), LineColor(kBlue));
  DataSimPdf.plotOn(DataFrameF,Slice(category,"fail"), ProjWData(category,DataCombData), Components("sgnPdfF"), LineColor(kRed), LineStyle(kSolid));
  DataSimPdf.plotOn(DataFrameF,Slice(category,"fail"), ProjWData(category,DataCombData), Components("bkgPdfF"), LineColor(kMagenta), LineStyle(kDashed));


  TCanvas *cPass = new TCanvas("fitCanvasP","canvas",10,30,650,600);
  cPass->SetGrid(0,0);
  cPass->SetFillStyle(4000);
  cPass->SetFillColor(10);
  cPass->SetTicky();
  cPass->SetObjectStat(0);

  cPass->cd();
  DataFrameP->Draw();
  string fileNameP = "fitCanvasPassElecTnP_"+tnp_+"_"+category_;
  cPass->SaveAs(Form("%s_%.2f.png",fileNameP.c_str(), binCenter_));

  TCanvas *cFail = new TCanvas("fitCanvasF","canvas",10,30,650,600);
  cFail->SetGrid(0,0);
  cFail->SetFillStyle(4000);
  cFail->SetFillColor(10);
  cFail->SetTicky();
  cFail->SetObjectStat(0);

  cFail->cd();
  DataFrameF->Draw();
  string fileNameF = "fitCanvasFailElecTnP_"+tnp_+"_"+category_;
  cFail->SaveAs(Form("%s_%.2f.png",fileNameF.c_str(), binCenter_));

  ResDataCombinedFit->printArgs(std::cout);
  cout << endl;
  ResDataCombinedFit->printValue(std::cout);
  cout << endl;

  float DataErrorLo = DataEffFit->getErrorLo()<0 ? DataEffFit->getErrorLo() : (-1)*DataEffFit->getError();
  float DataErrorHi = DataEffFit->getErrorHi()>0 ? DataEffFit->getErrorHi() : DataEffFit->getError();

  cout << DataEffFit->getVal() << " +/- " << DataEffFit->getError() << "  ( " << DataErrorLo << ", " << DataErrorHi << ")" <<  endl;

  float* out1 = new float[6];
  float* out2 = new float[6];

  out1[0]=(binCenter_);
  out1[1]=(binWidth_);
  out1[2]=(binWidth_);
  out1[3]=(McTruthEff);
  out1[4]=(BinomialError);
  out1[5]=(BinomialError);

  out2[0]=(binCenter_);
  out2[1]=(binWidth_);
  out2[2]=(binWidth_);
  out2[3]=(DataEffFit->getVal());
  out2[4]=((-1)*DataErrorLo);
  out2[5]=(DataErrorHi);

  out.push_back(out1);
  out.push_back(out2);

  return out;
}

//.........这里部分代码省略.........
  modelFail->plotOn(mframeFail,Components("bkgexpFail_shifted"),LineStyle(kDashed),LineColor(kGreen+2));
  mframeFail->GetYaxis()->SetTitle("");
  mframeFail->GetYaxis()->SetTitleOffset(1.2);
  mframeFail->GetXaxis()->SetTitle("m_{ee#gamma} [GeV/c^{2}]");
  mframeFail->GetXaxis()->SetTitleOffset(1.05);
  mframeFail->SetTitle("");
  mframeFail->Draw();
  cv->SaveAs(Form("ShiftedModel_%d.gif",Option));


  //*************************************************************************************
  //Do Toys
  //*************************************************************************************
  for(uint t=0; t < NToys; ++t) {

    RooDataSet *pseudoData_pass    = modelPass_default->generate(*m_default, randomnumber->Poisson(npass));
    RooDataSet *pseudoData_fail  = 0;
    pseudoData_fail    = modelFail->generate(*m_default, randomnumber->Poisson(nfail));
    RooDataSet *pseudoDataCombined = new RooDataSet("pseudoDataCombined","pseudoDataCombined",RooArgList(*m_default),
                                                    RooFit::Index(sample),
                                                    RooFit::Import("Pass",*pseudoData_pass),
                                                    RooFit::Import("Fail",*pseudoData_fail));

    pseudoDataCombined->write(Form("toy%d.txt",t));

    RooFitResult *fitResult=0;
    fitResult = totalPdf->fitTo(*pseudoDataCombined,
                                RooFit::Extended(),
                                RooFit::Strategy(2),
                                //RooFit::Minos(RooArgSet(eff)),
                                RooFit::Save());

    cout << "\n\n";
    cout << "Eff Fit: " << eff->getVal() << " -" << fabs(eff->getErrorLo()) << " +" << eff->getErrorHi() << endl;

    //Fill Tree
    varEff = eff->getVal();
    varEffErrL = fabs(eff->getErrorLo());
    varEffErrH = eff->getErrorHi();
    outTree->Fill();


//   //*************************************************************************************
//   //Plot Toys
//   //*************************************************************************************
//   TCanvas *cv = new TCanvas("cv","cv",800,600);
//   char pname[50];
//   char binlabelx[100];
//   char binlabely[100];
//   char yield[50];
//   char effstr[100];
//   char nsigstr[100];
//   char nbkgstr[100];
//   char chi2str[100];

//   //
//   // Plot passing probes
//   //

//   RooPlot *mframeFail_default = m.frame(Bins(Int_t(130-60)/2));
//   modelFail_default->plotOn(mframeFail_default);
//   modelFail_default->plotOn(mframeFail_default,Components("ebackgroundFail"),LineStyle(kDashed),LineColor(kRed));
//   modelFail_default->plotOn(mframeFail_default,Components("bkgexpFail"),LineStyle(kDashed),LineColor(kGreen+2));
//   mframeFail_default->Draw();
//   cv->SaveAs("DefaultModel.gif");

//.........这里部分代码省略.........
	// this is the new extended PDF whith the normalization restricted to the bin-area
        RooAbsPdf *extLimPdf = NULL;
	if( iscombcat ) {
	  extLimPdf = new RooSimultaneous("epdf","",finalcat);
	  // loop over the cats and generate temporary extended PDFs
	  for (int jcat=0; jcat<(catnames.size()-1); ++jcat) {
            RooRealVar *rvar = dynamic_cast<RooRealVar*>(normu_vec.at(jcat));
            if (rvar) rvar->setVal(fitcurves.at(jcat)->interpolate(center));
            RooExtendPdf *ecpdf = new RooExtendPdf(TString::Format("ecpdf%s",catnames.at(jcat).Data()),"",*pdfShape_vec.at(jcat),*normu_vec.at(jcat),"errRange");
            static_cast<RooSimultaneous*>(extLimPdf)->addPdf(*ecpdf,catnames.at(jcat));
          }
	} else
	  extLimPdf = new RooExtendPdf("extLimPdf","",*thisPdf,*nlim,"errRange");

        RooAbsReal *nll = extLimPdf->createNLL(*datanorm,Extended(),NumCPU(1));
        RooMinimizer minim(*nll);
        minim.setStrategy(0);
        double clone = 1.0 - 2.0*RooStats::SignificanceToPValue(1.0);
        double cltwo = 1.0 - 2.0*RooStats::SignificanceToPValue(2.0);
	
        if (iscombcat) minim.setStrategy(2);
        
        minim.migrad();
	
        if (!iscombcat) { 
          minim.minos(*nlim);
        }
        else {
          minim.hesse();
          nlim->removeAsymError();
        }

	if( verbose ) 
	  printf("errlo = %5f, errhi = %5f\n",nlim->getErrorLo(),nlim->getErrorHi());
        
        onesigma->SetPoint(i-1,center,nombkg);
        onesigma->SetPointError(i-1,0.,0.,-nlim->getErrorLo(),nlim->getErrorHi());
        
	// to get the 2-sigma bands...
        minim.setErrorLevel(0.5*pow(ROOT::Math::normal_quantile(1-0.5*(1-cltwo),1.0), 2)); // the 0.5 is because qmu is -2*NLL
                          // eventually if cl = 0.95 this is the usual 1.92!      
        
        if (!iscombcat) { 
          minim.migrad();
          minim.minos(*nlim);
        }
        else {
          nlim->setError(2.0*nlim->getError());
          nlim->removeAsymError();          
        }
	
        twosigma->SetPoint(i-1,center,nombkg);
        twosigma->SetPointError(i-1,0.,0.,-nlim->getErrorLo(),nlim->getErrorHi());      
        
        // for memory clean-up
        delete nll;
        delete extLimPdf;
      }
      
      hmass->setRange("errRange",massmin,massmax);

      if( verbose )
	onesigma->Print("V");
      
      // plot[icat] the error bands
      twosigma->SetLineColor(kGreen);

//.........这里部分代码省略.........
  pdfFail.plotOn(frame2,RooFit::ProjWData(*dataFail), 
  RooFit::Components(*bkgFailPdf),RooFit::LineColor(kRed));
  pdfFail.plotOn(frame2,RooFit::ProjWData(*dataFail));
  frame2->Draw("e0");

  plotlabel = new TPaveText(0.23,0.87,0.43,0.92,"NDC");
   plotlabel->SetTextColor(kBlack);
   plotlabel->SetFillColor(kWhite);
   plotlabel->SetBorderSize(0);
   plotlabel->SetTextAlign(12);
   plotlabel->SetTextSize(0.03);
   plotlabel->AddText("CMS Preliminary 2010");
 plotlabel2 = new TPaveText(0.23,0.82,0.43,0.87,"NDC");
   plotlabel2->SetTextColor(kBlack);
   plotlabel2->SetFillColor(kWhite);
   plotlabel2->SetBorderSize(0);
   plotlabel2->SetTextAlign(12);
   plotlabel2->SetTextSize(0.03);
   plotlabel2->AddText("#sqrt{s} = 7 TeV");
  plotlabel3 = new TPaveText(0.23,0.75,0.43,0.80,"NDC");
   plotlabel3->SetTextColor(kBlack);
   plotlabel3->SetFillColor(kWhite);
   plotlabel3->SetBorderSize(0);
   plotlabel3->SetTextAlign(12);
   plotlabel3->SetTextSize(0.03);
   sprintf(temp, "%.4f", LUMINOSITY);
   plotlabel3->AddText((string("#int#font[12]{L}dt = ") + 
                        temp + string(" pb^{ -1}")).c_str());
   plotlabel4 = new TPaveText(0.6,0.82,0.8,0.87,"NDC");
   plotlabel4->SetTextColor(kBlack);
   plotlabel4->SetFillColor(kWhite);
   plotlabel4->SetBorderSize(0);
   plotlabel4->SetTextAlign(12);
   plotlabel4->SetTextSize(0.03);
   err = ErrorInProduct(nsig, nErr, 1.0-e, eErr, corr);
   sprintf(temp, "Signal = %.2f #pm %.2f", NumSignalFail->getVal(), err);
   plotlabel4->AddText(temp);
   plotlabel5 = new TPaveText(0.6,0.77,0.8,0.82,"NDC");
   plotlabel5->SetTextColor(kBlack);
   plotlabel5->SetFillColor(kWhite);
   plotlabel5->SetBorderSize(0);
   plotlabel5->SetTextAlign(12);
   plotlabel5->SetTextSize(0.03);
   sprintf(temp, "Bkg = %.2f #pm %.2f", ParNumBkgFail->getVal(), ParNumBkgFail->getError());
   plotlabel5->AddText(temp);
   plotlabel6 = new TPaveText(0.6,0.87,0.8,0.92,"NDC");
   plotlabel6->SetTextColor(kBlack);
   plotlabel6->SetFillColor(kWhite);
   plotlabel6->SetBorderSize(0);
   plotlabel6->SetTextAlign(12);
   plotlabel6->SetTextSize(0.03);
   plotlabel6->AddText("Failing probes");
   plotlabel7 = new TPaveText(0.6,0.72,0.8,0.77,"NDC");
   plotlabel7->SetTextColor(kBlack);
   plotlabel7->SetFillColor(kWhite);
   plotlabel7->SetBorderSize(0);
   plotlabel7->SetTextAlign(12);
   plotlabel7->SetTextSize(0.03);
   sprintf(temp, "Eff = %.3f #pm %.3f", ParEfficiency->getVal(), ParEfficiency->getErrorHi(), ParEfficiency->getErrorLo());
   plotlabel7->AddText(temp);
   plotlabel8 = new TPaveText(0.6,0.72,0.8,0.66,"NDC");
   plotlabel8->SetTextColor(kBlack);
   plotlabel8->SetFillColor(kWhite);
   plotlabel8->SetBorderSize(0);
   plotlabel8->SetTextAlign(12);
   plotlabel8->SetTextSize(0.03);
   sprintf(temp, "#chi^{2}/DOF = %.3f", frame2->chiSquare());
   plotlabel8->AddText(temp);

//   plotlabel->Draw();
//   plotlabel2->Draw();
//   plotlabel3->Draw();
  plotlabel4->Draw();
  plotlabel5->Draw();
  plotlabel6->Draw();
  plotlabel7->Draw();
  plotlabel8->Draw();

  c2->SaveAs( cname + TString(".gif"));
//   c2->SaveAs( cname + TString(".eps"));
//   c2->SaveAs( cname + TString(".root"));
  canvasFile->WriteTObject(c2, c2->GetName(), "WriteDelete");

  canvasFile->Close();

  
  effTextFile.width(40);
  effTextFile << label;
  effTextFile.width(20);
  effTextFile  << setiosflags(ios::fixed) << setprecision(4) << left << ParEfficiency->getVal() ;
  effTextFile.width(20);
  effTextFile  << left << ParEfficiency->getErrorHi();
  effTextFile.width(20);
  effTextFile  << left << ParEfficiency->getErrorLo();
  effTextFile.width(14);
  effTextFile  << setiosflags(ios::fixed) << setprecision(2) << left << nSignalPass ;
  effTextFile.width(14);
  effTextFile << left << nSignalFail << endl;

}