C++ RooDataSet::printMultiline方法代码示例

   void ws_profile_interval1( const char* wsfile = "ws-test1.root", const char* parName = "mu_susy_sig", double alpha = 0.10, double mu_susy_sig_val = 0., double xmax = -1. ) {




       TFile* wstf = new TFile( wsfile ) ;

       RooWorkspace* ws = dynamic_cast<RooWorkspace*>( wstf->Get("ws") );
       ws->Print() ;






 ////  ModelConfig* modelConfig = (ModelConfig*) ws->obj( "SbModel" ) ;

 ////  printf("\n\n\n  ===== SbModel ====================\n\n") ;
 ////  modelConfig->Print() ;







       RooDataSet* rds = (RooDataSet*) ws->obj( "ra2b_observed_rds" ) ;
       printf("\n\n\n  ===== RooDataSet ====================\n\n") ;

       rds->Print() ;
       rds->printMultiline(cout, 1, kTRUE, "") ;









       printf("\n\n\n  ===== Grabbing %s rrv ====================\n\n", parName ) ;
       RooRealVar* rrv_par = ws->var( parName ) ;
       if ( rrv_par == 0x0 ) {
          printf("\n\n\n *** can't find %s in workspace.  Quitting.\n\n\n", parName ) ; return ;
       } else {
          printf(" current value is : %8.3f\n", rrv_par->getVal() ) ; cout << flush ;
       }
       if ( xmax > 0 ) { rrv_par->setMax( xmax ) ; }

       printf("\n\n\n  ===== Grabbing mu_susy_sig rrv ====================\n\n") ;
       RooRealVar* rrv_mu_susy_sig = ws->var("mu_susy_sig") ;
       if ( rrv_mu_susy_sig == 0x0 ) {
          printf("\n\n\n *** can't find mu_susy_sig in workspace.  Quitting.\n\n\n") ; return ;
       }
       if ( strcmp( parName, "mu_susy_sig" ) != 0 ) {
          if ( mu_susy_sig_val >= 0. ) {
             printf(" current value is : %8.3f\n", rrv_mu_susy_sig->getVal() ) ; cout << flush ;
             printf(" fixing to %8.2f.\n", mu_susy_sig_val ) ;
             rrv_mu_susy_sig->setVal( mu_susy_sig_val ) ;
             rrv_mu_susy_sig->setConstant(kTRUE) ;
          } else {
             printf(" current value is : %8.3f\n", rrv_mu_susy_sig->getVal() ) ; cout << flush ;
             printf(" allowing mu_susy_sig to float.\n") ;
             rrv_mu_susy_sig->setConstant(kFALSE) ;
          }
       } else {
          printf("\n\n profile plot parameter is mu_susy_sig.\n") ;
          rrv_mu_susy_sig->setConstant(kFALSE) ;
       }

       printf("\n\n\n  ===== Grabbing likelihood pdf ====================\n\n") ;
       RooAbsPdf* likelihood = ws->pdf("likelihood") ;
       if ( likelihood == 0x0 ) {
          printf("\n\n\n *** can't find likelihood pdf in workspace.  Quitting.\n\n\n") ; return ;
       } else {
          printf("\n\n likelihood pdf: \n\n") ;
          likelihood->Print() ;
       }












       printf("\n\n\n  ===== Doing a fit ====================\n\n") ;

       likelihood->fitTo( *rds ) ;

       double mlValue = rrv_par->getVal() ;
       printf("  Maximum likelihood value of %s : %8.3f +/- %8.3f\n",
            parName, rrv_par->getVal(), rrv_par->getError() ) ;



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

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


       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") ;


    }

   void constrained_scan( const char* wsfile = "outputfiles/ws.root",
                          const char* new_poi_name="mu_bg_4b_msig_met1",
                          double constraintWidth=1.5,
                          int npoiPoints = 20,
                          double poiMinVal = 0.,
                          double poiMaxVal = 10.0,
                          double ymax = 9.,
                          int verbLevel=1  ) {

      TString outputdir("outputfiles") ;

      gStyle->SetOptStat(0) ;

      TFile* wstf = new TFile( wsfile ) ;
      RooWorkspace* ws = dynamic_cast<RooWorkspace*>( wstf->Get("ws") );
      ws->Print() ;

      RooDataSet* rds = (RooDataSet*) ws->obj( "hbb_observed_rds" ) ;
      cout << "\n\n\n  ===== RooDataSet ====================\n\n" << endl ;
      rds->Print() ;
      rds->printMultiline(cout, 1, kTRUE, "") ;

      RooRealVar* rv_sig_strength = ws->var("sig_strength") ;
      if ( rv_sig_strength == 0x0 ) { printf("\n\n *** can't find sig_strength in workspace.\n\n" ) ; return ; }

      RooAbsPdf* likelihood = ws->pdf("likelihood") ;
      if ( likelihood == 0x0 ) { printf("\n\n *** can't find likelihood in workspace.\n\n" ) ; return ; }
      printf("\n\n Likelihood:\n") ;
      likelihood -> Print() ;



      /////rv_sig_strength -> setConstant( kFALSE ) ;
      rv_sig_strength -> setVal(0.) ;
      rv_sig_strength -> setConstant( kTRUE ) ;

      likelihood->fitTo( *rds, Save(false), PrintLevel(0), Hesse(true), Strategy(1) ) ;
      //RooFitResult* fitResult = likelihood->fitTo( *rds, Save(true), PrintLevel(0), Hesse(true), Strategy(1) ) ;
      //double minNllSusyFloat = fitResult->minNll() ;
      //double susy_ss_atMinNll = rv_sig_strength -> getVal() ;

      RooMsgService::instance().getStream(1).removeTopic(Minimization) ;
      RooMsgService::instance().getStream(1).removeTopic(Fitting) ;



     //-- Construct the new POI parameter.
      RooAbsReal* new_poi_rar(0x0) ;

      new_poi_rar = ws->var( new_poi_name ) ;
      if ( new_poi_rar == 0x0 ) {
         printf("\n\n New POI %s is not a variable.  Trying function.\n\n", new_poi_name ) ;
         new_poi_rar = ws->function( new_poi_name ) ;
         if ( new_poi_rar == 0x0 ) {
            printf("\n\n New POI %s is not a function.  I quit.\n\n", new_poi_name ) ;
            return ;
         } else {
            printf("\n Found it.\n\n") ;
         }
      } else {
         printf("\n\n     New POI %s is a variable with current value %.1f.\n\n", new_poi_name, new_poi_rar->getVal() ) ;
      }

       double startPoiVal = new_poi_rar->getVal() ;


       RooAbsReal* nll = likelihood -> createNLL( *rds, Verbose(true) ) ;

       RooRealVar* rrv_poiValue = new RooRealVar( "poiValue", "poiValue", 0., -10000., 10000. ) ;

       RooRealVar* rrv_constraintWidth = new RooRealVar("constraintWidth","constraintWidth", 0.1, 0.1, 1000. ) ;
       rrv_constraintWidth -> setVal( constraintWidth ) ;
       rrv_constraintWidth -> setConstant(kTRUE) ;


       RooMinuit* rminuit( 0x0 ) ;


       RooMinuit* rminuit_uc = new RooMinuit( *nll  ) ;

       rminuit_uc->setPrintLevel(verbLevel-1) ;
       rminuit_uc->setNoWarn() ;

       rminuit_uc->migrad() ;
       rminuit_uc->hesse() ;

       RooFitResult* rfr_uc = rminuit_uc->fit("mr") ;

       double floatParInitVal[10000] ;
       char   floatParName[10000][100] ;
       int nFloatParInitVal(0) ;
       RooArgList ral_floats = rfr_uc->floatParsFinal() ;
       TIterator* floatParIter = ral_floats.createIterator() ;
       {
          RooRealVar* par ;
          while ( (par = (RooRealVar*) floatParIter->Next()) ) {
             sprintf( floatParName[nFloatParInitVal], "%s", par->GetName() ) ;
             floatParInitVal[nFloatParInitVal] = par->getVal() ;
             nFloatParInitVal++ ;
          }
//.........这里部分代码省略.........

   void ws_constrained_profile3D( const char* wsfile = "rootfiles/ws-data-unblind.root",
                                   const char* new_poi_name = "n_M234_H4_3b",
                                   int npoiPoints = 20,
                                   double poiMinVal = 0.,
                                   double poiMaxVal = 20.,
                                   double constraintWidth = 1.5,
                                   double ymax = 10.,
                                   int verbLevel=0 ) {


     gStyle->SetOptStat(0) ;

     //--- make output directory.

     char command[10000] ;
     sprintf( command, "basename %s", wsfile ) ;
     TString wsfilenopath = gSystem->GetFromPipe( command ) ;
     wsfilenopath.ReplaceAll(".root","") ;
     char outputdirstr[1000] ;
     sprintf( outputdirstr, "outputfiles/scans-%s", wsfilenopath.Data() ) ;
     TString outputdir( outputdirstr ) ;


     printf("\n\n Creating output directory: %s\n\n", outputdir.Data() ) ;
     sprintf(command, "mkdir -p %s", outputdir.Data() ) ;
     gSystem->Exec( command ) ;


     //--- Tell RooFit to shut up about anything less important than an ERROR.
      RooMsgService::instance().setGlobalKillBelow(RooFit::ERROR) ;



       if ( verbLevel > 0 ) { printf("\n\n Verbose level : %d\n\n", verbLevel) ; }


       TFile* wstf = new TFile( wsfile ) ;

       RooWorkspace* ws = dynamic_cast<RooWorkspace*>( wstf->Get("ws") );

       if ( verbLevel > 0 ) { ws->Print() ; }






       RooDataSet* rds = (RooDataSet*) ws->obj( "ra2b_observed_rds" ) ;

       if ( verbLevel > 0 ) {
          printf("\n\n\n  ===== RooDataSet ====================\n\n") ;
          rds->Print() ;
          rds->printMultiline(cout, 1, kTRUE, "") ;
       }





       ModelConfig* modelConfig = (ModelConfig*) ws->obj( "SbModel" ) ;
       RooAbsPdf* likelihood = modelConfig->GetPdf() ;

       RooRealVar* rrv_mu_susy_all0lep = ws->var("mu_susy_all0lep") ;
       if ( rrv_mu_susy_all0lep == 0x0 ) {
          printf("\n\n\n *** can't find mu_susy_all0lep in workspace.  Quitting.\n\n\n") ;
          return ;
       }





       //-- do BG only.
       rrv_mu_susy_all0lep->setVal(0.) ;
       rrv_mu_susy_all0lep->setConstant( kTRUE ) ;










       //-- do a prefit.

       printf("\n\n\n ====== Pre fit with unmodified nll var.\n\n") ;

       RooFitResult* dataFitResultSusyFixed = likelihood->fitTo(*rds, Save(true),Hesse(false),Minos(false),Strategy(1),PrintLevel(verbLevel));
       int dataSusyFixedFitCovQual = dataFitResultSusyFixed->covQual() ;
       if ( dataSusyFixedFitCovQual < 2 ) { printf("\n\n\n *** Failed fit!  Cov qual %d.  Quitting.\n\n", dataSusyFixedFitCovQual ) ; return ; }
       double dataFitSusyFixedNll = dataFitResultSusyFixed->minNll() ;

       if ( verbLevel > 0 ) {
          dataFitResultSusyFixed->Print("v") ;
       }

       printf("\n\n Nll value, from fit result : %.3f\n\n", dataFitSusyFixedNll ) ;

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

   void fitqual_plots( const char* wsfile = "outputfiles/ws.root", const char* plottitle="" ) {

      TText* tt_title = new TText() ;
      tt_title -> SetTextAlign(33) ;

      gStyle -> SetOptStat(0) ;
      gStyle -> SetLabelSize( 0.06, "y" ) ;
      gStyle -> SetLabelSize( 0.08, "x" ) ;
      gStyle -> SetLabelOffset( 0.010, "y" ) ;
      gStyle -> SetLabelOffset( 0.010, "x" ) ;
      gStyle -> SetTitleSize( 0.07, "y" ) ;
      gStyle -> SetTitleSize( 0.05, "x" ) ;
      gStyle -> SetTitleOffset( 1.50, "x" ) ;
      gStyle -> SetTitleH( 0.07 ) ;
      gStyle -> SetPadLeftMargin( 0.15 ) ;
      gStyle -> SetPadBottomMargin( 0.15 ) ;
      gStyle -> SetTitleX( 0.10 ) ;

      gDirectory->Delete("h*") ;

      TFile* wstf = new TFile( wsfile ) ;

      RooWorkspace* ws = dynamic_cast<RooWorkspace*>( wstf->Get("ws") );
      ws->Print() ;

      int bins_of_met = TMath::Nint( ws->var("bins_of_met")->getVal()  ) ;
      printf("\n\n Bins of MET : %d\n\n", bins_of_met ) ;

      int bins_of_nb = TMath::Nint( ws->var("bins_of_nb")->getVal()  ) ;
      printf("\n\n Bins of nb : %d\n\n", bins_of_nb ) ;

      int nb_lookup[10] ;
      if ( bins_of_nb == 2 ) {
         nb_lookup[0] = 2 ;
         nb_lookup[1] = 4 ;
      } else if ( bins_of_nb == 3 ) {
         nb_lookup[0] = 2 ;
         nb_lookup[1] = 3 ;
         nb_lookup[2] = 4 ;
      }

      TCanvas* cfq1 = (TCanvas*) gDirectory->FindObject("cfq1") ;
      if ( cfq1 == 0x0 ) {
         if ( bins_of_nb == 3 ) {
            cfq1 = new TCanvas("cfq1","hbb fit", 700, 1000 ) ;
         } else if ( bins_of_nb == 2 ) {
            cfq1 = new TCanvas("cfq1","hbb fit", 700, 750 ) ;
         } else {
            return ;
         }
      }

      RooRealVar* rv_sig_strength = ws->var("sig_strength") ;
      if ( rv_sig_strength == 0x0 ) { printf("\n\n *** can't find sig_strength in workspace.\n\n" ) ; return ; }

      ModelConfig* modelConfig = (ModelConfig*) ws->obj( "SbModel" ) ;

      RooDataSet* rds = (RooDataSet*) ws->obj( "hbb_observed_rds" ) ;

      rds->Print() ;
      rds->printMultiline(cout, 1, kTRUE, "") ;

      RooAbsPdf* likelihood = modelConfig->GetPdf() ;

      ///RooFitResult* fitResult = likelihood->fitTo( *rds, Save(true), PrintLevel(0) ) ;
      RooFitResult* fitResult = likelihood->fitTo( *rds, Save(true), PrintLevel(3) ) ;
      fitResult->Print() ;


      char hname[1000] ;
      char htitle[1000] ;
      char pname[1000] ;




     //-- unpack observables.

      int obs_N_msig[10][50] ; // first index is n btags, second is met bin.
      int obs_N_msb[10][50]  ; // first index is n btags, second is met bin.

      const RooArgSet* dsras = rds->get() ;
      TIterator* obsIter = dsras->createIterator() ;
      while ( RooRealVar* obs = (RooRealVar*) obsIter->Next() ) {
         for ( int nbi=0; nbi<bins_of_nb; nbi++ ) {
            for ( int mbi=0; mbi<bins_of_met; mbi++ ) {
               sprintf( pname, "N_%db_msig_met%d", nb_lookup[nbi], mbi+1 ) ;
               if ( strcmp( obs->GetName(), pname ) == 0 ) { obs_N_msig[nbi][mbi] = TMath::Nint( obs -> getVal() ) ; }
               sprintf( pname, "N_%db_msb_met%d", nb_lookup[nbi], mbi+1 ) ;
               if ( strcmp( obs->GetName(), pname ) == 0 ) { obs_N_msb[nbi][mbi] = TMath::Nint( obs -> getVal() ) ; }
            } // mbi.
         } // nbi.
      } // obs iterator.


      printf("\n\n") ;
      for ( int nbi=0; nbi<bins_of_nb; nbi++ ) {
         printf(" nb=%d :  ", nb_lookup[nbi] ) ;
         for ( int mbi=0; mbi<bins_of_met; mbi++ ) {
            printf("  sig=%3d, sb=%3d  |", obs_N_msig[nbi][mbi], obs_N_msb[nbi][mbi] ) ;
//.........这里部分代码省略.........

float ComputeTestStat(TString wsfile, double mu_susy_sig_val) {

  gROOT->Reset();

  TFile* wstf = new TFile( wsfile ) ;

  RooWorkspace* ws = dynamic_cast<RooWorkspace*>( wstf->Get("ws") );
  ws->Print() ;
  
  ModelConfig* modelConfig = (ModelConfig*) ws->obj( "SbModel" ) ;
  
  modelConfig->Print() ;

  RooDataSet* rds = (RooDataSet*) ws->obj( "ra2b_observed_rds" ) ;
  
  rds->Print() ;
  rds->printMultiline(cout, 1, kTRUE, "") ;
  
  RooAbsPdf* likelihood = modelConfig->GetPdf() ;
  
  RooRealVar* rrv_mu_susy_sig = ws->var("mu_susy_all0lep") ;
  if ( rrv_mu_susy_sig == 0x0 ) {
    printf("\n\n\n *** can't find mu_susy_all0lep in workspace.  Quitting.\n\n\n") ;
    return ;
  } else {
    printf(" current value is : %8.3f\n", rrv_mu_susy_sig->getVal() ) ; cout << flush ;
    rrv_mu_susy_sig->setConstant(kFALSE) ;
  }

  /*
  // check the impact of varying the qcd normalization:

  RooRealVar *rrv_qcd_0lepLDP_ratioH1 = ws->var("qcd_0lepLDP_ratio_H1");
  RooRealVar *rrv_qcd_0lepLDP_ratioH2 = ws->var("qcd_0lepLDP_ratio_H2");
  RooRealVar *rrv_qcd_0lepLDP_ratioH3 = ws->var("qcd_0lepLDP_ratio_H3");
  
  rrv_qcd_0lepLDP_ratioH1->setVal(0.3);
  rrv_qcd_0lepLDP_ratioH2->setVal(0.3);
  rrv_qcd_0lepLDP_ratioH3->setVal(0.3);
  
  rrv_qcd_0lepLDP_ratioH1->setConstant(kTRUE);
  rrv_qcd_0lepLDP_ratioH2->setConstant(kTRUE);
  rrv_qcd_0lepLDP_ratioH3->setConstant(kTRUE);
  */
  
  printf("\n\n\n  ===== Doing a fit with SUSY component floating ====================\n\n") ;

  RooFitResult* fitResult = likelihood->fitTo( *rds, Save(true), PrintLevel(0) ) ;
  double logLikelihoodSusyFloat = fitResult->minNll() ;
  
  double logLikelihoodSusyFixed(0.) ;
  double testStatVal(-1.) ;
  if ( mu_susy_sig_val >= 0. ) {
    printf("\n\n\n  ===== Doing a fit with SUSY fixed ====================\n\n") ;
    printf(" fixing mu_susy_sig to %8.2f.\n", mu_susy_sig_val ) ;
    rrv_mu_susy_sig->setVal( mu_susy_sig_val ) ;
    rrv_mu_susy_sig->setConstant(kTRUE) ;
    
    fitResult = likelihood->fitTo( *rds, Save(true), PrintLevel(0) ) ;
    logLikelihoodSusyFixed = fitResult->minNll() ;
    testStatVal = 2.*(logLikelihoodSusyFixed - logLikelihoodSusyFloat) ;
    printf("\n\n\n ======= test statistic : -2 * ln (L_fixed / ln L_max) = %8.3f\n\n\n", testStatVal ) ;
  }


  return testStatVal ;

}

//.........这里部分代码省略.........
   ///    sprintf( sel, "1BL" ) ;
   /// } else if ( strstr( wsfile, "2BL" ) != 0 ) {
   ///    sprintf( sel, "2BL" ) ;
   /// } else if ( strstr( wsfile, "3B" ) != 0 ) {
   ///    sprintf( sel, "3B" ) ;
   /// } else if ( strstr( wsfile, "1BT" ) != 0 ) {
   ///    sprintf( sel, "1BT" ) ;
   /// } else if ( strstr( wsfile, "2BT" ) != 0 ) {
   ///    sprintf( sel, "2BT" ) ;
   /// } else {
   ///    printf("\n\n\n *** can't figure out which selection this is.  I quit.\n\n" ) ;
   ///    return ;
   /// }
   /// printf("\n\n selection is %s\n\n", sel ) ;




       TFile* wstf = new TFile( wsfile ) ;

       RooWorkspace* ws = dynamic_cast<RooWorkspace*>( wstf->Get("ws") );

       ws->Print() ;






       RooDataSet* rds = (RooDataSet*) ws->obj( "ra2b_observed_rds" ) ;
       printf("\n\n\n  ===== RooDataSet ====================\n\n") ;

       rds->Print() ;
       rds->printMultiline(cout, 1, kTRUE, "") ;





       ModelConfig* modelConfig = (ModelConfig*) ws->obj( "SbModel" ) ;
       RooAbsPdf* likelihood = modelConfig->GetPdf() ;

       const RooArgSet* nuisanceParameters = modelConfig->GetNuisanceParameters() ;

       RooRealVar* rrv_mu_susy_sig = ws->var("mu_susy_sig") ;
       if ( rrv_mu_susy_sig == 0x0 ) {
          printf("\n\n\n *** can't find mu_susy_sig in workspace.  Quitting.\n\n\n") ;
          return ;
       }






 ////  printf("\n\n\n  ===== Doing a fit ====================\n\n") ;

 ////  RooFitResult* preFitResult = likelihood->fitTo( *rds, Save(true) ) ;
 ////  const RooArgList preFitFloatVals = preFitResult->floatParsFinal() ;
 ////  {
 ////    TIterator* parIter = preFitFloatVals.createIterator() ;
 ////    while ( RooRealVar* par = (RooRealVar*) parIter->Next() ) {
 ////       printf(" %20s : %8.2f\n", par->GetName(), par->getVal() ) ;
 ////    }
 ////  }