Belle II Software  release-05-01-25
SpaceResolutionCalibration Class Reference

Class for Space resolution calibration. More...

#include <SpaceResolutionCalibration.h>

Collaboration diagram for SpaceResolutionCalibration:

Public Member Functions

 SpaceResolutionCalibration ()
 Constructor.
 
virtual ~SpaceResolutionCalibration ()
 Destructor.
 
virtual void setDebug (bool debug=false)
 Debug or not.
 
virtual void setUseDB (bool useDB=false)
 Use database or text mode.
 
virtual void setMinimumNDF (double minndf)
 minimum NDF required for track
 
virtual void setMinimumPval (double minPval)
 Minimum Pval required.
 
virtual void setBinWidth (double bw)
 Bin width of each slide.
 
virtual void BField (bool bfield)
 Work with B field or not;.
 
virtual void inputFileNames (std::string inputname)
 Input root file names, results of collector module.
 
virtual void setStoreHisto (bool storeHist=false)
 Store histograms durring the calibration or not.
 
virtual void ProfileFileNames (std::string profileFileName)
 File name describe theta/alpha bin, if don't want to use default from input sigma.
 
virtual void useProfileFromInputSigma (bool useProfileFromInputSigma)
 use sigma bin profile form input sigma or new one from input file
 
virtual void setSigmaFileName (std::string name)
 Output sigma file name, for text mode.
 
void execute ()
 execute all, make the interface the same as CAF
 

Protected Member Functions

virtual bool calibrate ()
 Run algo on data. More...
 
virtual void createHisto ()
 create histogram
 
virtual void readProfile ()
 read sigma bining (alpha, theta bining)
 
virtual void readSigma ()
 read sigma from previous calibration, (input sigma)
 
virtual void readSigmaFromDB ()
 read sigma from DB More...
 
virtual void readSigmaFromText ()
 read sigma from text file
 
virtual void storeHisto ()
 store histogram
 
virtual void write ()
 save calibration, in text file or db
 

Private Member Functions

double getUpperBoundaryForFit (TGraphErrors *graph)
 search max point at boundary region
 

Private Attributes

double m_ndfmin = 5
 Minimum NDF

 
double m_Pvalmin = 0.
 Minimum Prob(chi2) of track.
 
double m_binWidth = 0.05
 width of each bin, unit cm
 
bool m_debug = false
 Debug or not.
 
bool m_draw = false
 print out histogram in pdf file or not
 
bool m_storeHisto = false
 Store histogram or not.
 
bool m_useDB = false
 use db or text mode
 
bool m_useProfileFromInputSigma = true
 Use binning from old sigma or new one form input.
 
bool m_BField = true
 Work with BField, fit range and initial parameters is different incase B and noB.
 
double sigma_old [56][2][18][7][8]
 old sigma prameters.
 
double sigma_new [56][2][18][7][8]
 new sigma prameters.
 
TF1 * ffit [56][2][18][7]
 fitting function
 
TGraphErrors * gfit [56][2][18][7]
 sigma*sigma graph for fit
 
TGraphErrors * gr [56][2][18][7]
 sigma graph.
 
TH2F * hist_b [56][2][Max_nalpha][Max_ntheta]
 2D histogram of biased residual
 
TH2F * hist_u [56][2][Max_nalpha][Max_ntheta]
 2D histogram of unbiased residual
 
TH1F * hu_m [56][2][Max_nalpha][Max_ntheta]
 mean histogram biased residual
 
TH1F * hu_s [56][2][Max_nalpha][Max_ntheta]
 sigma histogram of biased residual
 
TH1F * hb_m [56][2][Max_nalpha][Max_ntheta]
 mean histogram of unbiased residual
 
TH1F * hb_s [56][2][Max_nalpha][Max_ntheta]
 sigma histogram of ubiased residual
 
int m_fitflag [56][2][Max_nalpha][Max_ntheta] = {{{{0}}}}
 Fit flag; 1:OK ; 0:error.
 
std::string m_outputSigmaFileName = "sigma_new.dat"
 Output sigma file name.
 
std::string m_inputRootFileNames = "rootfile/output*"
 Input root file names.
 
std::string m_ProfileFileName = "sigma_profile"
 Profile file name.
 
DBObjPtr< CDCSpaceResols > * m_sResolFromDB
 Database for sigma.
 
std::string m_sigmafile = "cdc/data/sigma.dat"
 Sigma file name, for text mode.
 
int m_firstExperiment
 First experiment.
 
int m_firstRun
 First run.
 
int m_lastExperiment
 Last experiment.
 
int m_lastRun
 Last run.
 
int m_nalpha
 number of alpha bins
 
int m_ntheta
 number of theta bins
 
double l_alpha [18]
 Lower boundays of alpha bins.
 
double u_alpha [18]
 Upper boundays of alpha bins.
 
double ialpha [18]
 represented alphas of alpha bins.
 
double l_theta [7]
 Lower boundays of theta bins.
 
double u_theta [7]
 Upper boundays of theta bins.
 
double itheta [7]
 represented alphas of theta bins.
 
int nalpha_old
 number of alpha bins from input
 
int ntheta_old
 number of theta bins from input
 
double l_alpha_old [18]
 Lower boundays of alpha bins from input.
 
double u_alpha_old [18]
 Upper boundays of alpha bins from input.
 
double ialpha_old [18]
 represented alphas of alpha bins from input.
 
double l_theta_old [7]
 Lower boundays of theta bins from input.
 
double u_theta_old [7]
 Upper boundays of theta bins from input.
 
double itheta_old [7]
 represented alphas of theta bins from input.
 
unsigned short m_sigmaParamMode_old
 sigma mode from input.
 
unsigned short m_sigmaParamMode = 1
 sigma mode for this calibration.
 

Static Private Attributes

static const int Max_nalpha = 18
 Maximum alpha bin.
 
static const int Max_ntheta = 7
 maximum theta bin

 
static const unsigned short Max_np = 40
 Maximum number of point =1/binwidth.
 

Detailed Description

Class for Space resolution calibration.

Definition at line 18 of file SpaceResolutionCalibration.h.

Member Function Documentation

◆ calibrate()

bool calibrate ( )
protectedvirtual

Run algo on data.

Upper limit of fitting.

Definition at line 278 of file SpaceResolutionCalibration.cc.

279 {
280  gROOT->SetBatch(1);
281  gErrorIgnoreLevel = 3001;
282  createHisto();
283  B2INFO("Start to calibrate");
284  gSystem->Exec("mkdir -p Sigma_Fit_Err"); //create a folder to store error histo
285 
286  TF1* func = new TF1("func", "[0]/(x*x + [1])+[2]* x+[3]+[4]*exp([5]*(x-[6])*(x-[6]))", 0, 1.);
287  TH1F* hprob = new TH1F("h1", "", 20, 0, 1);
288  double upFit;
289  double intp6;
290 
291  for (int i = 0; i < 56; ++i) {
292  for (int lr = 0; lr < 2; ++lr) {
293  for (int al = 0; al < m_nalpha; ++al) {
294  for (int th = 0; th < m_ntheta; ++th) {
295  if (m_fitflag[i][lr][al][th] != -1) { /*if graph exist, do fitting*/
296 
297 
298  upFit = getUpperBoundaryForFit(gfit[i][lr][al][th]);
299  intp6 = upFit + 0.2;
300  B2DEBUG(199, "xmax for fitting: " << upFit);
301  // func->SetLineColor(1 + lr + al * 2 + th * 3);
302  func->SetParameters(5E-6, 0.007, 1E-4, 1E-5, 0.00008, -30, intp6);
303  func->SetParLimits(0, 1E-7, 1E-4);
304  func->SetParLimits(1, 0.00001, 0.02);
305  func->SetParLimits(2, 1E-6, 0.0005);
306  func->SetParLimits(3, 1E-8, 0.0005);
307  func->SetParLimits(4, 0., 0.001);
308  func->SetParLimits(5, -40, 0.);
309  func->SetParLimits(6, intp6 - 0.5, intp6 + 0.3);
310  B2DEBUG(199, "FITTING for layer: " << i << "lr: " << lr << " ial" << al << " ith:" << th);
311  B2DEBUG(199, "Fit flag before fit:" << m_fitflag[i][lr][al][th]);
312  // if(!(gfit[i][lr][al][th]->isValid())) continue;
313  // m_fitflag[i][lr][al][th] = 0;
314  for (int j = 0; j < 10; j++) {
315 
316  B2DEBUG(199, "loop: " << j);
317  B2DEBUG(199, "Int p6: " << intp6);
318  B2DEBUG(199, "Number of Point: " << gfit[i][lr][al][th]->GetN());
319  Int_t stat = gfit[i][lr][al][th]->Fit("func", "MQE", "", 0.05, upFit);
320  B2DEBUG(199, "stat of fit" << stat);
321  std::string Fit_status = gMinuit->fCstatu.Data();
322  B2DEBUG(199, "FIT STATUS: " << Fit_status);
323  // stat=gfit[i]->Fit(Form("ffit[%d]",i),"M "+Q,"",0.0,cellsize(i)+0.05+j*0.005);
324  if (Fit_status == "OK" ||
325  Fit_status == "SUCCESSFUL" ||
326  Fit_status == "CALL LIMIT" ||
327  Fit_status == "PROBLEMS") {
328  if (fabs(func->Eval(0.3)) > 0.00035 || func->Eval(0.3) < 0) {
329  func->SetParameters(5E-6, 0.007, 1E-4, 1E-7, 0.0007, -30, intp6 + 0.05 * j);
330  // func->SetParameters(defaultparsmall);
331  m_fitflag[i][lr][al][th] = 0;
332  } else {
333  B2DEBUG(199, "Prob of fit: " << func->GetProb());
334  m_fitflag[i][lr][al][th] = 1;
335  break;
336  }
337  } else {
338  m_fitflag[i][lr][al][th] = 0;
339  func->SetParameters(5E-6, 0.007, 1E-4, 1E-7, 0.0007, -30, intp6 + 0.05 * j);
340  // upFit += 0.025;
341  if (j == 9) {
342  TCanvas* c1 = new TCanvas("c1", "", 600, 600);
343  gfit[i][lr][al][th]->Draw();
344  c1->SaveAs(Form("Sigma_Fit_Err/%d_%d_al%d_th%d_%s.png", i, lr, al, th, Fit_status.c_str()));
345  B2WARNING("Fit error: " << i << " " << lr << " " << al << " " << th);
346  }
347  }
348  }
349  if (m_fitflag[i][lr][al][th] == 1) {
350  B2DEBUG(199, "ProbFit: Lay_lr_al_th: " << i << " " << lr << " " << al << " " << th << func->GetProb());
351  hprob->Fill(func->GetProb());
352  func->GetParameters(sigma_new[i][lr][al][th]);
353  }
354  }
355  }
356  }
357  }
358  }
359  storeHisto();
360  write();
361 
362  return true;
363 }

◆ readSigmaFromDB()

void readSigmaFromDB ( )
protectedvirtual

read sigma from DB

< angle bin info.

Definition at line 564 of file SpaceResolutionCalibration.cc.


The documentation for this class was generated from the following files:
Belle2::CDC::SpaceResolutionCalibration::createHisto
virtual void createHisto()
create histogram
Definition: SpaceResolutionCalibration.cc:31
Belle2::CDC::SpaceResolutionCalibration::write
virtual void write()
save calibration, in text file or db
Definition: SpaceResolutionCalibration.cc:397
Belle2::CDC::SpaceResolutionCalibration::sigma_new
double sigma_new[56][2][18][7][8]
new sigma prameters.
Definition: SpaceResolutionCalibration.h:97
Belle2::CDC::SpaceResolutionCalibration::m_fitflag
int m_fitflag[56][2][Max_nalpha][Max_ntheta]
Fit flag; 1:OK ; 0:error.
Definition: SpaceResolutionCalibration.h:107
Belle2::CDC::SpaceResolutionCalibration::storeHisto
virtual void storeHisto()
store histogram
Definition: SpaceResolutionCalibration.cc:364
Belle2::CDC::SpaceResolutionCalibration::gfit
TGraphErrors * gfit[56][2][18][7]
sigma*sigma graph for fit
Definition: SpaceResolutionCalibration.h:99
Belle2::CDC::SpaceResolutionCalibration::m_ntheta
int m_ntheta
number of theta bins
Definition: SpaceResolutionCalibration.h:120
Belle2::CDC::SpaceResolutionCalibration::m_nalpha
int m_nalpha
number of alpha bins
Definition: SpaceResolutionCalibration.h:119
Belle2::CDC::SpaceResolutionCalibration::getUpperBoundaryForFit
double getUpperBoundaryForFit(TGraphErrors *graph)
search max point at boundary region
Definition: SpaceResolutionCalibration.h:144