C++ RooFitResult::constPars方法代码示例

  void  accessParams(int n, double array[] ){
    /*
    RooArgList fl=fitresult->floatParsFinal();
    RooArgList cl=fitresult->constPars();
    fl.Print();
    cl.Print();
    int fn=fl.getSize();
    int cn=cl.getSize();
    */    //    mean,area
    RooRealVar* par;
    RooRealVar* pas;
    if (n==1){
                     par = (RooRealVar*) fitresult->floatParsFinal().find("mean1");
    if (par==NULL){  par = (RooRealVar*) fitresult->constPars().find("mean1"); }
                     pas = (RooRealVar*) fitresult->floatParsFinal().find("area1");
    if (pas==NULL){  pas = (RooRealVar*) fitresult->constPars().find("area1"); }
    //    printf( "   mean1 == %f +- %f\n", par->getVal() ,  par->getError()  );
    array[0]= par->getVal();     array[1]= par->getError();
    array[2]= pas->getVal();     array[3]= pas->getError();
    }
    if (npeaks<2) return;

    if (n==2){
                     par = (RooRealVar*) fitresult->floatParsFinal().find("mean2");
    if (par==NULL){  par = (RooRealVar*) fitresult->constPars().find("mean2"); }
                     pas = (RooRealVar*) fitresult->floatParsFinal().find("area2");
    if (pas==NULL){  pas = (RooRealVar*) fitresult->constPars().find("area2"); }
    //    printf( "   mean2 == %f +- %f\n", par->getVal() ,  par->getError()  );
    array[0]= par->getVal();     array[1]= par->getError();
     array[2]= pas->getVal();     array[3]= pas->getError();
   }
    if (npeaks<3) return;
    
    if (n==3){
                     par = (RooRealVar*) fitresult->floatParsFinal().find("mean3");
    if (par==NULL){  par = (RooRealVar*) fitresult->constPars().find("mean3"); }
                     pas = (RooRealVar*) fitresult->floatParsFinal().find("area3");
    if (pas==NULL){  pas = (RooRealVar*) fitresult->constPars().find("area3"); }
    //    printf( "   mean2 == %f +- %f\n", par->getVal() ,  par->getError()  );
    array[0]= par->getVal();     array[1]= par->getError();
    array[2]= pas->getVal();     array[3]= pas->getError();
    }
    if (npeaks<4) return;
    
    if (n==4){
                     par = (RooRealVar*) fitresult->floatParsFinal().find("mean4");
    if (par==NULL){  par = (RooRealVar*) fitresult->constPars().find("mean4"); }
                     pas = (RooRealVar*) fitresult->floatParsFinal().find("area4");
    if (pas==NULL){  pas = (RooRealVar*) fitresult->constPars().find("area4"); }
    //    printf( "   mean2 == %f +- %f\n", par->getVal() ,  par->getError()  );
    array[0]= par->getVal();     array[1]= par->getError();
     array[2]= pas->getVal();     array[3]= pas->getError();
   }
    if (npeaks<5) return;
    
    if (n==5){
                     par = (RooRealVar*) fitresult->floatParsFinal().find("mean5");
    if (par==NULL){  par = (RooRealVar*) fitresult->constPars().find("mean5"); }
                     pas = (RooRealVar*) fitresult->floatParsFinal().find("area5");
    if (pas==NULL){  pas = (RooRealVar*) fitresult->constPars().find("area5"); }
    //    printf( "   mean2 == %f +- %f\n", par->getVal() ,  par->getError()  );
    array[0]= par->getVal();     array[1]= par->getError();
     array[2]= pas->getVal();     array[3]= pas->getError();
   }
    if (npeaks<6) return;
    


    return;
  }//----------access

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


       RooGaussian* pdf_acc_ee = new RooGaussian( "pdf_acc_ee", "Gaussian pdf for Z to ee acceptance",
                       *rv_acc_ee, RooConst( acc_ee_mean ), RooConst( acc_ee_err ) ) ;




       RooGaussian* pdf_eff_mm = new RooGaussian( "pdf_eff_mm", "Gaussian pdf for Z to mumu efficiency",
                       *rv_eff_mm, RooConst( eff_mm_mean ), RooConst( eff_mm_err ) ) ;



       RooGaussian* pdf_eff_ee = new RooGaussian( "pdf_eff_ee", "Gaussian pdf for Z to ee efficiency",
                       *rv_eff_ee, RooConst( eff_ee_mean ), RooConst( eff_ee_err ) ) ;



       RooPoisson* pdf_Nsbee   = new RooPoisson( "pdf_Nsbee"  , "Nsb , Z to ee Poisson PDF", *rv_Nsbee  , *rv_n_sbee  ) ;
       RooPoisson* pdf_Nsbmm   = new RooPoisson( "pdf_Nsbmm"  , "Nsb , Z to mm Poisson PDF", *rv_Nsbmm  , *rv_n_sbmm  ) ;
       RooPoisson* pdf_Nsigee  = new RooPoisson( "pdf_Nsigee" , "Nsig, Z to ee Poisson PDF", *rv_Nsigee , *rv_n_sigee ) ;
       RooPoisson* pdf_Nsigmm  = new RooPoisson( "pdf_Nsigmm" , "Nsig, Z to mm Poisson PDF", *rv_Nsigmm , *rv_n_sigmm ) ;

       RooArgSet pdflist ;
       pdflist.add( *pdf_acc_mm ) ;
       pdflist.add( *pdf_acc_ee ) ;
       pdflist.add( *pdf_eff_mm ) ;
       pdflist.add( *pdf_eff_ee ) ;
       pdflist.add( *pdf_Nsbee  ) ;
       pdflist.add( *pdf_Nsbmm  ) ;
       pdflist.add( *pdf_Nsigee  ) ;
       pdflist.add( *pdf_Nsigmm  ) ;

       RooProdPdf* znnLikelihood = new RooProdPdf( "znnLikelihood", "Z to nunu likelihood", pdflist ) ;


       RooArgSet observedParametersList ;
       observedParametersList.add( *rv_Nsbee ) ;
       observedParametersList.add( *rv_Nsbmm ) ;
       observedParametersList.add( *rv_Nsigee ) ;
       observedParametersList.add( *rv_Nsigmm ) ;

       RooDataSet* dsObserved = new RooDataSet( "ztonn_rds", "Z to nunu dataset", observedParametersList ) ;
       dsObserved->add( observedParametersList ) ;
       //// RooDataSet* dsObserved = znnLikelihood->generate( observedParametersList, 1) ;

       RooWorkspace* znnWorkspace = new RooWorkspace("ztonn_ws") ;
       znnWorkspace->import( *znnLikelihood ) ;
       znnWorkspace->import( *dsObserved ) ;

       znnWorkspace->Print() ;

       dsObserved->printMultiline(cout, 1, kTRUE, "") ;
       RooFitResult* fitResult = znnLikelihood->fitTo( *dsObserved, Verbose(true), Save(true) ) ;

       printf("\n\n----- Constant parameters:\n") ;
       RooArgList constPars = fitResult->constPars() ;
       for ( int pi=0; pi<constPars.getSize(); pi++ ) {
          constPars[pi].Print() ;
       } // pi.

       printf("\n\n----- Floating parameters:\n") ;
       RooArgList floatPars = fitResult->floatParsFinal() ;
       for ( int pi=0; pi<floatPars.getSize(); pi++ ) {
          floatPars[pi].Print() ;
       } // pi.
       printf("\n\n") ;



       ProfileLikelihoodCalculator plc_znn_sb(  *dsObserved, *znnLikelihood, RooArgSet( *rv_mu_Znnsb ) ) ;
       ProfileLikelihoodCalculator plc_znn_sig( *dsObserved, *znnLikelihood, RooArgSet( *rv_mu_Znnsig ) ) ;

       plc_znn_sb.SetTestSize(0.32) ;
       plc_znn_sig.SetTestSize(0.32) ;

       ConfInterval* sb_znn_interval = plc_znn_sb.GetInterval() ;
       float sbZnnLow  = ((LikelihoodInterval*) sb_znn_interval)->LowerLimit(*rv_mu_Znnsb) ;
       float sbZnnHigh = ((LikelihoodInterval*) sb_znn_interval)->UpperLimit(*rv_mu_Znnsb) ;
       printf("\n\n znn SB interval %6.1f to %6.1f\n", sbZnnLow, sbZnnHigh ) ;

       ConfInterval* sig_znn_interval = plc_znn_sig.GetInterval() ;
       float sigZnnLow  = ((LikelihoodInterval*) sig_znn_interval)->LowerLimit(*rv_mu_Znnsig) ;
       float sigZnnHigh = ((LikelihoodInterval*) sig_znn_interval)->UpperLimit(*rv_mu_Znnsig) ;
       printf("\n\n znn SIG interval %6.1f to %6.1f\n", sigZnnLow, sigZnnHigh ) ;



       TCanvas* c_prof_sb = new TCanvas("c_prof_sb","SB Z to nunu profile") ;
       LikelihoodIntervalPlot plot_znn_sb((LikelihoodInterval*)sb_znn_interval) ;
       plot_znn_sb.Draw() ;
       c_prof_sb->SaveAs("znn_sb_profile.png") ;

       TCanvas* c_prof_sig = new TCanvas("c_prof_sig","sig Z to nunu profile") ;
       LikelihoodIntervalPlot plot_znn_sig((LikelihoodInterval*)sig_znn_interval) ;
       plot_znn_sig.Draw() ;
       c_prof_sig->SaveAs("znn_sig_profile.png") ;


    }

//.........这里部分代码省略.........
  lefth->Sumw2();
  righth->Sumw2();
  righth23->Scale(lefth23->Integral()/righth23->Integral());
  lefth23->Sumw2();
  righth23->Sumw2();

  RooDataHist *lefthist = new RooDataHist("lefthist","lefthist",RooArgSet(m2Kpi_d0mass,m2Kpi0_d0mass),lefth);
  RooDataHist *righthist = new RooDataHist("righthist","righthist",RooArgSet(m2Kpi_d0mass,m2Kpi0_d0mass),righth);
  RooDataHist *lefthist23 = new RooDataHist("lefthist23","lefthist23",RooArgSet(m2Kpi_d0mass,m2pipi0_d0mass),lefth23);
  RooDataHist *righthist23 = new RooDataHist("righthist23","righthist23",RooArgSet(m2Kpi_d0mass,m2pipi0_d0mass),righth23);

  RooHistPdf leftpdf("leftpdf","leftpdf",RooArgSet(m2Kpi_d0mass,m2Kpi0_d0mass),*lefthist,4);
  RooHistPdf rightpdf("rightpdf","rightpdf",RooArgSet(m2Kpi_d0mass,m2Kpi0_d0mass),*righthist,4);
  RooHistPdf leftpdf23("leftpdf23","leftpdf23",RooArgSet(m2Kpi_d0mass,m2pipi0_d0mass),*lefthist23,4);
  RooHistPdf rightpdf23("rightpdf23","rightpdf23",RooArgSet(m2Kpi_d0mass,m2pipi0_d0mass),*righthist23,4);

  RooRealVar fcomb("fcomb","fcomb",0.738);
  RooAddPdf combdalitz("combdalitz","combdalitz",RooArgList(leftpdf,rightpdf),RooArgList(fcomb));
  RooAddPdf combdalitz23("combdalitz23","combdalitz23",RooArgList(leftpdf23,rightpdf23),RooArgList(fcomb));
  ///////////////////////////////////////

  RooAddPdf totpdf("totpdf","totpdf",RooArgList(*D0pdf,combdalitz),RooArgList(Nsig,Nbkg));
  RooAddPdf totpdf23("totpdf23","totpdf23",RooArgList(*D0pdf23,combdalitz23),RooArgList(Nsig,Nbkg));

  if(doFit){
    // Start Minuit session on Chi2
    RooChi2Var chi2("chi2","chi2",totpdf,*dbdalitz);
    RooMinuit m2(chi2);
    m2.migrad();
    m2.hesse();
    RooFitResult* fitRes = m2.save();
    fitRes->Print("v");
    RooArgSet results(fitRes->floatParsFinal());
    RooArgSet conresults(fitRes->constPars());
    results.add(conresults);
    results.writeToFile("fit_isobar_RS.txt");

    //save the stupid result
    TFile f("fit_RSDalitz_result.root","RECREATE");
    fitRes->Write();
    f.Close();
  }

  if(doFract) {

    cout << "Calculating fit fractions" << endl;
    TFile f("fit_RSDalitz_result.root");
    RooFitResult* fitRes = (RooFitResult*)f.Get("chi2");

    //now calculate the fit fractions
    const Int_t nRes = D0pdf->getManager()->getnRes();

    //recalculate the normalization if necessary
    D0pdf->getManager()->calNorm();

    EvtComplex normarray[nRes][nRes];

    const Int_t myRes = 12;
    TH1F fitty[myRes];

    //read the integral value from the cache file.
    //In this way we don't need to compute the normalization everytime during MIGRAD
    char int_name[50];
    D0pdf->getManager()->getFileName(int_name);

    ifstream f1;