save and plot the number of fitted pi0s

Calibration/Pi0Calibration/Pi0Calibration/Pi0MassFitter.h

@@ -42,6 +42,7 @@ namespace Calibration
                 void   likelihoodfit( const std::string& outputdir, const std::string& outfilename, bool verbose );
 		std::pair<double,double> getMean() {return m_mean; }
 		std::pair<double,double> getSigma() {return m_sigma; }
+		std::pair<double,double> getNSignal() {return m_nsignal; }
                 int    getStatus() {return m_status; }
 
 	    private:
@@ -53,6 +54,7 @@ namespace Calibration
 
 		std::pair<double,double> m_mean;
 		std::pair<double,double> m_sigma;
+		std::pair<double,double> m_nsignal;
                 int    m_status;
 	};
     }//namespace Pi0Calibration

Calibration/Pi0Calibration/scripts/call_do_calibration.py

 #!/bin/python
 import os
 outputdir = os.path.join( "/home3/xu/worksapce/CALO/pi0calibration/data2016/TEST" )
-datatype = 'TH2D'
+datatype = 'MMAP'
 if datatype=='ROOT':
     tuplefile = "/gpfs/LAPP-DATA/lhcb/xu/Second_data_2016_pi0.root"
 if datatype=='MMAP':
@@ -10,12 +10,12 @@ if datatype=='TH2D':
     tuplefile = "/gpfs/LAPP-DATA/lhcb/xu/Second_Worker1.root"
 
 if '__main__' == __name__ :
-    cmd = "~/cmtuser/AlignmentDev_v11r1p1/run python do_calibration.py --filename %s --filetype %s "%( tuplefile, datatype )
+    cmd = "~/cmtuser/AlignmentDev_v11r2/run python do_calibration.py --filename %s --filetype %s "%( tuplefile, datatype )
     cmd+= "--outdir %s "%outputdir
     cmd+= "--indir %s  "%outputdir
     cmd+= "--year 2016 "
 
     import subprocess
-    for iteration in xrange(1,7+1):
+    for iteration in xrange(1,1+1):
         cmd+= "--iteration %s "%iteration
         subprocess.call( cmd , shell=True) 

Calibration/Pi0Calibration/src/Pi0CalibrationAlg.cpp

@@ -138,8 +138,10 @@ StatusCode Pi0CalibrationAlg::finalize(){
 
   std::map<unsigned int, std::pair<double,double> > meanValues;
   std::map<unsigned int, std::pair<double,double> > sigmaValues;
+  std::map<unsigned int, std::pair<double,double> > npi0Values;
   meanValues.clear();
   sigmaValues.clear();
+  npi0Values.clear();
   m_skippedCells.clear();
   
   if(msgLevel(MSG::INFO)) info() << " fit the pi0 mass distribution to get the initial lambda values and the resolutions in each region " << endmsg; 
@@ -169,6 +171,7 @@ StatusCode Pi0CalibrationAlg::finalize(){
 
     meanValues[cell.index()]   = fitter.getMean();
     sigmaValues[cell.index()]  = fitter.getSigma();
+    npi0Values[cell.index()]   = fitter.getNSignal();
     
     auto lold   = lMap.get_lambda( cell.index() );
     if( !fitter.getStatus() ){
@@ -193,8 +196,10 @@ StatusCode Pi0CalibrationAlg::finalize(){
   saveCells( m_skippedCells );
   saveValues(  meanValues, "mean" );
   saveValues( sigmaValues, "sigma" );
+  saveValues( npi0Values, "npi0" );
   meanValues.clear();
   sigmaValues.clear();
+  npi0Values.clear();
   return sc;
 }
 

Calibration/Pi0Calibration/src/Pi0CalibrationMonitor.cpp

@@ -139,8 +139,8 @@ StatusCode Pi0CalibrationMonitor::execute(){
 StatusCode Pi0CalibrationMonitor::finalize(){
   StatusCode sc = StatusCode::SUCCESS;
 
-  std::map<std::string, double> expects = { {"mean",m_pdgPi0}, {"sigma",10.0}, {"lambda",1.0} };
-  for( auto varname : {"mean", "sigma", "lambda"} ){
+  std::map<std::string, double> expects = { {"mean",m_pdgPi0}, {"sigma",10.0}, {"lambda",1.0}, {"npi0",5000} };
+  for( auto varname : {"mean", "sigma", "lambda", "npi0"} ){
     std::map<unsigned int, std::pair<double,double> > Values;
     readValues( varname, Values );
     monitorPlots(  varname,  expects[varname], Values );

Calibration/Pi0Calibration/src/Pi0MassFitter.cpp

@@ -60,8 +60,9 @@ Pi0MassFitter::Pi0MassFitter(TH1* hist, TH1* hist_bg, double mpi0, double sigma)
   m_mpi0    = mpi0;
   m_sigma0 = sigma;
 
-  m_mean  = std::make_pair(0.,0.);
-  m_sigma = std::make_pair(0.,0.);
+  m_mean    = std::make_pair(0.,0.);
+  m_sigma   = std::make_pair(0.,0.);
+  m_nsignal = std::make_pair(0.,0.);
   m_status= -1;
 }
 
@@ -218,8 +219,11 @@ void Pi0MassFitter::chi2fit( const std::string& outputdir, const std::string& ou
     }
   }
 
-  m_mean  = std::make_pair( model->GetParameter(1), model->GetParError(1) );
-  m_sigma = std::make_pair( model->GetParameter(2), model->GetParError(2) );
+  m_mean    = std::make_pair( model->GetParameter(1), model->GetParError(1) );
+  m_sigma   = std::make_pair( model->GetParameter(2), model->GetParError(2) );
+  double nS  =  model->GetParameter(0)*model->GetParameter(3);
+  double nSE =  std::sqrt( model->GetParError(0)/model->GetParameter(0)*model->GetParError(0)/model->GetParameter(0) + model->GetParError(3)/model->GetParameter(3)*model->GetParError(3)/model->GetParameter(3) );
+  m_nsignal = std::make_pair( nS, nSE*nS ); 
   if( model->GetParameter(2) <15.0 && model->GetParameter(2)>5.0 && model->GetParameter(3)>0.01 ) m_status = 0;
   else { verbose = true; subdir = "bad"; }
   if(verbose){