Files
JFCalibration/CuFluo_analysis_sc_singlethread_fits.cpp

828 lines
32 KiB
C++

// fit fluo spectrum per pixel for a given storage cell, only G0
// (can be called multiple times to handle multi-storage-cell data, each time with different Sc number as argument, i.e. supports manual multi-threading on the same data set)
// and save peak position and uncertainty
// Takes pedestal-corrected data (in root file) as input
// Call CuFluo_analysis_sc_singlethread_data before to analyse raw data files and create root file
#include "../sls_detector_calibration/jungfrauCommonHeader.h"
#include "../sls_detector_calibration/jungfrauCommonFunctions.h"
#include "../sls_detector_calibration/jungfrauFile.C"
#include "../sls_detector_calibration/jungfrauPixelMask.C"
#include "../sls_detector_calibration/jungfrauPedestal.C"
#include "../sls_detector_calibration/energyCalibration.h"
#include "../sls_detector_calibration/energyCalibration.cpp"
#include "TGraph.h"
#include "TGraphErrors.h"
#include "TF1.h"
#include "TFile.h"
#include "TPaveStats.h"
#include "TLegend.h"
#include "TPaveText.h"
#include <sys/stat.h>
#include <vector>
#include <iostream>
#include <string>
int main(int argc, char* argv[]) {
jungfrauStyle();
gStyle->SetOptFit(11);
if (argc != 8) {
cout << "Correct usage:" << endl;
cout << "arg 1: specify module number" << endl;
cout << "arg 2: specify HG0 or G0" << endl;
cout << "arg 3: specify data location folder" << endl;
cout << "arg 4: specify pede file prefix" << endl;
cout << "arg 5: specify data file prefix" << endl;
cout << "arg 6: specify location of analysis root files" << endl;
cout << "arg 7: specify storage cell number" << endl;
cout << " " << endl;
exit(1);
}
string module_str = argv[1];
string gain_str = argv[2];
bool isJF11=false;
if (gain_str == "HG0JF11") {
gain_str = "HG0";
isJF11=true;
}
string data_loc = argv[3];
string pede_file = argv[4];
string data_file = argv[5];
string anadata_loc = argv[6];
bool save_to_afs = false;
int createHistoFile = 1;
char histoname[256];
char savename[256];
char plotfolder[256];
char rootdatafolder[256];
int filen = 352;
int pedefilen = 5;
int this_storagecell = stoi(argv[7]); // needs compiler flag -std=c++11
// create necessary directories with permissions drwxrwxr-x
// data/Mxxx
sprintf(savename,"data/M%s", module_str.c_str());
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
// plots/Mxxx
sprintf(savename,"plots/M%s", module_str.c_str());
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
// plots/Mxxx/CuFluo
sprintf(savename,"plots/M%s/CuFluo", module_str.c_str());
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
// plots/Mxxx/CuFluo/G0
sprintf(savename,"plots/M%s/CuFluo/G0", module_str.c_str());
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
// plots/Mxxx/CuFluo/HG0
sprintf(savename,"plots/M%s/CuFluo/HG0", module_str.c_str());
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
sprintf(savename,"%s", anadata_loc.c_str());
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
if ( save_to_afs ) {
cout << "Saving to afs." << endl;
sprintf( plotfolder, "plots/M%s/CuFluo/%s", module_str.c_str(), gain_str.c_str() );
sprintf( rootdatafolder, "data/M%s", module_str.c_str() );
} else {
cout << "Not saving to afs." << endl;
// data/Mxxx
sprintf( savename,"%s/data", anadata_loc.c_str() );
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
// plots/Mxxx
sprintf( savename,"%s/plots", anadata_loc.c_str() );
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
// plots/Mxxx/CuFluo
sprintf( savename,"%s/plots/CuFluo", anadata_loc.c_str() );
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
// plots/Mxxx/CuFluo/G0
sprintf( savename,"%s/plots/CuFluo/G0", anadata_loc.c_str() );
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
// plots/Mxxx/CuFluo/HG0
sprintf( savename,"%s/plots/CuFluo/HG0", anadata_loc.c_str() );
mkdir(savename, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
sprintf( plotfolder, "%s/plots/CuFluo/%s", anadata_loc.c_str(), gain_str.c_str() );
sprintf( rootdatafolder, "%s/data", anadata_loc.c_str() );
}
jungfrauPixelMask *pixelMaskObject = new jungfrauPixelMask();
bool pixel_mask [NCH];
pixelMaskObject->initialisePixelMask(pixel_mask);
TCanvas *mapcanvas = new TCanvas("mapcanvas","",150,10,800,400);
mapcanvas->SetLeftMargin(0.1);
mapcanvas->SetRightMargin(0.13);
mapcanvas->SetTopMargin(0.08);
mapcanvas->SetBottomMargin(0.15);
TCanvas* c1 = new TCanvas("c1","");
/*****************************************************************************************************************/
/* NOTE: Since arrays, as they were used before, tend to cause memory leaks if they grow too large on the stack, */
/* it might be smarter to realize all array objects with std::vector instead. */
/*****************************************************************************************************************/
if (createHistoFile) {
jungfrauFile *thisfile = new jungfrauFile();
jungfrauPedestal* pedestalObject_SC = new jungfrauPedestal();
pedestalObject_SC->pedestalSetNFrames(100);
for (int pedefilei = 0; pedefilei < pedefilen; ++pedefilei) {
// open pede file
sprintf( savename, "%s/%s_%%6.6d.dat", data_loc.c_str(), pede_file.c_str() );
thisfile->open( (char*)savename, pedefilei );
//count events in file
int nevents = 0;
while ( thisfile->readNextFrame() ) {
nevents++;
}
thisfile->rewind();
cout << "read " << nevents << " events" << endl;
while ( thisfile->readNextFrame() ) {
// calculate pixel mask
pixelMaskObject->maskIfGainNot( 0, thisfile->getFrameDataHandle(), pixel_mask );
// caluclate pedestals
if( thisfile->currentSCnumber() == (uint64_t)this_storagecell ) {
pedestalObject_SC->addFrameToPedestalCalculation( thisfile->getFrameDataHandle() );
}
}
thisfile->close();
} // end of loops over files
cout << "after chip mask, n masked pixels is " << pixelMaskObject->getNMasked(pixel_mask) << endl;
} // end if i creatHistoFile
// We have saved the data into root-file by first calling CuFluo_analysis_sc_singlethread_data.cpp. Now, we read from that same root-file and perform the fitting.
c1->cd();
sprintf(savename,"%s/CuFluo_%s_sc%d_file0to%d.root", anadata_loc.c_str(), gain_str.c_str(), this_storagecell, filen-1); //uncomment for VH 210906
TFile* comb_file = new TFile((const char *)(savename),"READ");
pixelMaskObject->initialisePixelMask(pixel_mask);
int low_ADU_peak = 0;
int high_ADU_peak = 0;
if (gain_str == "HG0") {
low_ADU_peak = 700;
high_ADU_peak = 900;
if (isJF11) {
low_ADU_peak = 850;
high_ADU_peak = 1350;
}
} else if (gain_str == "G0") {
low_ADU_peak = 200; // was 250
high_ADU_peak = 400;
}
//declare
TH1F* fit_par3;
TH1F* fit_par4;
TH1F* fit_par5;
TH1F* fit_par6;
TH1F* fit_par7;
TH2F* fit_par3_2d;
TH2F* fit_par4_2d;
TH2F* fit_par5_2d;
TH2F* fit_par6_2d;
TH2F* fit_par7_2d;
TH1F* peak_fit_pos;
TH1F* peak_fit_poserr;
TH2F* peak_fit_pos_2d;
TH2F* peak_fit_poserr_2d;
TH1F* noise_fit_pos;
TH1F* noise_fit_poserr;
TH2F* noise_fit_pos_2d;
TH2F* noise_fit_poserr_2d;
TH1F* gain_fit;
TH1F* gain_fiterr;
TH2F* gain_fit_2d;
TH2F* gain_fiterr_2d;
TH1F* gain_fit_isEdge;
TH1F* gain_fit_isInnerEdge;
TH1F* gain_fit_isDouble;
TH1F* gain_fit_isNextToDouble;
TH1F* gain_fit_isQuad;
TH1F* gain_fit_isBulk;
TH2F* gain_ADUper1keV_2d;
TH2F* gainerr_ADUper1keV_2d;
//initialize
Char_t *fithistoname = new Char_t[50];
snprintf( fithistoname, 50, "fit_par3_sc%d", this_storagecell );
fit_par3 = new TH1F( fithistoname, "", 100, 0, 50 );
snprintf( fithistoname, 50, "fit_par4_sc%d", this_storagecell );
fit_par4 = new TH1F( fithistoname, "", 100, 0, 500 );
snprintf( fithistoname, 50, "fit_par5_sc%d", this_storagecell );
fit_par5 = new TH1F( fithistoname, "", 100, 0, 0.5 );
snprintf( fithistoname, 50, "fit_par6_sc%d", this_storagecell );
fit_par6 = new TH1F( fithistoname, "", 100, 1.05, 1.3 ); // was 1.05, 1.25
snprintf( fithistoname, 50, "fit_par7_sc%d", this_storagecell );
fit_par7 = new TH1F( fithistoname, "", 100, 0, 0.4 );
snprintf( fithistoname, 50, "fit_par3_2d_sc%d", this_storagecell );
fit_par3_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "fit_par4_2d_sc%d", this_storagecell );
fit_par4_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "fit_par5_2d_sc%d", this_storagecell );
fit_par5_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "fit_par6_2d_sc%d", this_storagecell );
fit_par6_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "fit_par7_2d_sc%d", this_storagecell );
fit_par7_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "peak_fit_pos_sc%d", this_storagecell );
peak_fit_pos = new TH1F( fithistoname, "", 100, low_ADU_peak, high_ADU_peak );
snprintf( fithistoname, 50, "peak_fit_poserr_sc%d", this_storagecell );
peak_fit_poserr = new TH1F( fithistoname, "", 100, 0, 2 );
snprintf( fithistoname, 50, "peak_fit_pos_2d_sc%d", this_storagecell );
peak_fit_pos_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "peak_fit_poserr_2d_sc%d", this_storagecell );
peak_fit_poserr_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "noise_fit_pos_sc%d", this_storagecell );
noise_fit_pos = new TH1F( fithistoname, "", 100, -10, 10 );
snprintf( fithistoname, 50, "noise_fit_poserr_sc%d", this_storagecell );
noise_fit_poserr = new TH1F( fithistoname, "", 100, 0, 0.15 ); // was 0, 0.1
snprintf( fithistoname, 50, "noise_fit_pos_2d_sc%d", this_storagecell );
noise_fit_pos_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "noise_fit_poserr_2d_sc%d", this_storagecell );
noise_fit_poserr_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "gain_fit_sc%d", this_storagecell );
gain_fit = new TH1F( fithistoname, "", 100, low_ADU_peak, high_ADU_peak );
snprintf( fithistoname, 50, "gain_fiterr_sc%d", this_storagecell );
gain_fiterr = new TH1F( fithistoname, "", 100, 0, 2 );
snprintf( fithistoname, 50, "gain_fit_2d_sc%d", this_storagecell );
gain_fit_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "gain_fiterr_2d_sc%d", this_storagecell );
gain_fiterr_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "gain_fit_isEdge_sc%d", this_storagecell );
gain_fit_isEdge = new TH1F( fithistoname, "", 100, low_ADU_peak, high_ADU_peak );
snprintf( fithistoname, 50, "gain_fit_isInnerEdge_sc%d", this_storagecell );
gain_fit_isInnerEdge = new TH1F( fithistoname, "", 100, low_ADU_peak, high_ADU_peak );
snprintf( fithistoname, 50, "gain_fit_isDouble_sc%d", this_storagecell );
gain_fit_isDouble = new TH1F( fithistoname, "", 100, low_ADU_peak, high_ADU_peak );
snprintf( fithistoname, 50, "gain_fit_isNextToDouble_sc%d", this_storagecell );
gain_fit_isNextToDouble = new TH1F( fithistoname, "", 100, low_ADU_peak, high_ADU_peak );
snprintf( fithistoname, 50, "gain_fit_isQuad_sc%d", this_storagecell );
gain_fit_isQuad = new TH1F( fithistoname, "", 100, low_ADU_peak, high_ADU_peak );
snprintf( fithistoname, 50, "gain_fit_isBulk_sc%d", this_storagecell );
gain_fit_isBulk = new TH1F( fithistoname, "", 100, low_ADU_peak, high_ADU_peak );
snprintf( fithistoname, 50, "gain_ADUper1keV_2d_sc%d", this_storagecell );
gain_ADUper1keV_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
snprintf( fithistoname, 50, "gainerr_ADUper1keV_2d_sc%d", this_storagecell );
gainerr_ADUper1keV_2d = new TH2F( fithistoname, "", NC, -0.5, NC-0.5, NR, -0.5, NR-0.5 );
int low_bin_noise = 101;
int high_bin_noise = 301;
int low_bin_peak = 0;
int high_bin_peak = 0;
if (gain_str == "HG0") {
low_bin_peak = 701;
high_bin_peak = 1200;
if (isJF11) {
low_bin_peak = 801;
high_bin_peak = 1451;
}
} else if (gain_str == "G0") {
low_bin_peak = 301;
high_bin_peak = 651;
}
for ( int j = 1; j < 9; ++j ) {
cout << "SC = " << this_storagecell << " slice " << j << endl;
sprintf( histoname, "adc2d_%d_sc%d", j, this_storagecell );
TH2I* adc2d_j = (TH2I*)( comb_file->Get( (const char *)(histoname) ) );
cout << adc2d_j->GetEntries() << endl;
adc2d_j->Draw("colz");
c1->Update();
for (int i=(65536*(j-1)); i<(65536*(j)); ++i) {
if (i%10000==0){cout << "another 10k" << endl;}
TH1D* proj = adc2d_j->ProjectionX("bin1",i-(65536*(j-1))+1,i-(65536*(j-1))+1);
if ( proj->Integral( low_bin_noise, high_bin_noise ) != 0 && proj->Integral( low_bin_peak, high_bin_peak ) != 0 ) {
// noise
TH1D *proj_noise = dynamic_cast<TH1D*>(proj->Rebin(4,"proj_noise"));
proj_noise->SetStats(kTRUE);
proj_noise->GetXaxis()->SetRangeUser(proj->GetBinLowEdge(low_bin_noise),proj->GetBinLowEdge(high_bin_noise+1));
proj_noise->Fit("gaus","Q");
TF1 *fit = proj_noise->GetFunction("gaus");
noise_fit_pos->Fill(fit->GetParameter(1));
noise_fit_pos_2d->Fill(i%NC,i/NC,fit->GetParameter(1));
noise_fit_poserr->Fill(fit->GetParError(1));
noise_fit_poserr_2d->Fill(i%NC,i/NC,fit->GetParError(1));
// peak
TH1D *proj_peak = dynamic_cast<TH1D*>(proj->Rebin(4,"proj_peak"));
proj_peak->SetStats(kTRUE);
proj_peak->GetXaxis()->SetRangeUser(proj->GetBinLowEdge(low_bin_peak),proj->GetBinLowEdge(high_bin_peak+1));
Double_t mypar[8];
mypar[0] = 0.0;
mypar[1] = 0.0;
mypar[2] = proj_peak->GetBinCenter(proj_peak->GetMaximumBin());
if (gain_str == "G0") {
mypar[3] = 16.;
} else if (gain_str == "HG0") {
mypar[3] = 29.;
}
mypar[4] = proj_peak->GetBinContent(proj_peak->GetMaximumBin());
if (gain_str == "G0") {
mypar[5] = 0.17;
} else if (gain_str == "HG0") {
mypar[5] = 0.14;
}
mypar[6] = 1.12;
if (gain_str == "G0") {
mypar[7] = 0.12;
} else if (gain_str == "HG0") {
mypar[7] = 0.14;
}
Double_t emypar[8];
energyCalibration *thiscalibration = new energyCalibration();
thiscalibration->setScanSign(1);
thiscalibration->setStartParametersKb(mypar);
thiscalibration->fixParameter(0,0.); // no background
thiscalibration->fixParameter(1,0.);
TF1* fittedfun = thiscalibration->fitSpectrumKb( proj_peak, mypar, emypar );
fit_par3->Fill( mypar[3] );
fit_par4->Fill( mypar[4] );
fit_par5->Fill( mypar[5] );
fit_par6->Fill( mypar[6] );
fit_par7->Fill( mypar[7] );
fit_par3_2d->Fill( i%NC, i/NC, mypar[3] );
fit_par4_2d->Fill( i%NC, i/NC, mypar[4] );
fit_par5_2d->Fill( i%NC, i/NC, mypar[5] );
fit_par6_2d->Fill( i%NC, i/NC, mypar[6] );
fit_par7_2d->Fill( i%NC, i/NC, mypar[7] );
peak_fit_pos->Fill( mypar[2] );
peak_fit_poserr->Fill( emypar[2] );
peak_fit_pos_2d->Fill( i%NC, i/NC, mypar[2] );
peak_fit_poserr_2d->Fill( i%NC, i/NC, emypar[2] );
if ( ( i >= 58000 && i < 58000+10 ) || // bulk
( i >= 10 && i < 10+10 ) || // edge
( i >= 1024+10 && i < 1024+10+10 ) || // inner edge
( i >= ( 256*1024 )+10 && i < ( 256*1024 )+10+10 ) || // double
( i >= ( 257*1024 )+10 && i < ( 257*1024 )+10+10 ) || // next to double
( i == ( 255*1024 )+255 ) // quad
) {
string pixel_type = "x";
if ( i >= 58000 && i < 58000+10 ) {
pixel_type = "b";
} else if ( i >= 10 && i < 10+10 ) {
pixel_type = "e";
} else if ( i >= 1024+10 && i < 1024+10+10 ) {
pixel_type = "ie";
} else if ( i >= ( 256*1024 )+10 && i < ( 256*1024 )+10+10 ) {
pixel_type = "d";
} else if ( i >= ( 257*1024 )+10 && i < ( 257*1024 )+10+10 ) {
pixel_type = "ntd";
} else if ( i == ( 255*1024 )+255 ) {
pixel_type = "q";
}
proj_noise->Draw();
c1->Update();
proj_noise->GetXaxis()->SetTitle( "Pedestal corrected ADC [ADU]" );
proj_noise->GetXaxis()->SetRangeUser( -100, 150 );
fit->SetParNames( "N_{#gamma}", "Peak pos", "Noise RMS" );
TPaveStats *st0 = (TPaveStats*)proj_noise->FindObject( "stats" );
st0->SetX1NDC(0.53);
st0->SetX2NDC(0.94);
st0->SetY1NDC(0.75);
st0->SetY2NDC(0.94);
st0->SetBorderSize(0);
st0->SetTextSize(0.04);
sprintf( savename, "%s/noise_%s_sc%d_%d_%s_M%s.png", plotfolder, pixel_type.c_str(), this_storagecell, i, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
TF1 *gaus_Ka = new TF1( "gaus_Ka", "gaus", proj->GetBinLowEdge(low_bin_peak), proj->GetBinLowEdge(high_bin_peak+1) );
gaus_Ka->SetParameters( mypar[4], mypar[2], mypar[3] );
gaus_Ka->SetLineColor( kBlue );
TF1 *erfc_Ka = new TF1( "erfc_Ka", "[0]/2.*(TMath::Erfc(([1]*(x-[2])/[3])/(TMath::Sqrt(2.))))", proj->GetBinLowEdge(low_bin_peak), proj->GetBinLowEdge(high_bin_peak+1) );
erfc_Ka->SetParameters( mypar[4]*mypar[5], 1, mypar[2], mypar[3] );
erfc_Ka->SetLineColor( kOrange );
TF1 *gaus_Kb = new TF1( "gaus_Kb", "gaus", proj->GetBinLowEdge(low_bin_peak), proj->GetBinLowEdge(high_bin_peak+1) );
gaus_Kb->SetParameters( mypar[4]*mypar[7], mypar[6]*mypar[2], mypar[3] );
gaus_Kb->SetLineColor( kGreen+2 );
TF1 *erfc_Kb = new TF1( "erfc_Kb", "[0]/2.*(TMath::Erfc(([1]*(x-[2])/[3])/(TMath::Sqrt(2.))))", proj->GetBinLowEdge(low_bin_peak), proj->GetBinLowEdge(high_bin_peak+1) );
erfc_Kb->SetParameters( mypar[4]*mypar[7]*mypar[5], 1, mypar[6]*mypar[2], mypar[3] );
erfc_Kb->SetLineColor( kOrange+7 );
proj_peak->Draw();
erfc_Kb->Draw("same");
erfc_Ka->Draw("same");
gaus_Kb->Draw("same");
gaus_Ka->Draw("same");
fittedfun->Draw("same");
c1->Update();
proj_peak->GetXaxis()->SetTitle( "Pedestal corrected ADC [ADU]" );
fittedfun->SetParNames( "Bkg height", "Bkg grad", "K_{#alpha} pos", "Noise RMS", "K_{#alpha} height", "CS", "K_{#beta}/K_{#alpha} pos", "K_{#beta} frac" );
TPaveStats *st = (TPaveStats*)proj_peak->FindObject( "stats" );
st->SetX1NDC(0.22);
st->SetX2NDC(0.62);
st->SetY1NDC(0.7);
st->SetY2NDC(0.94);
st->SetBorderSize(0);
st->SetTextSize(0.04);
sprintf( savename, "%s/peak_%s_sc%d_%d_%s_M%s.png", plotfolder, pixel_type.c_str(), this_storagecell, i, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
}
// gain
gain_fit->Fill( mypar[2] - fit->GetParameter(1) );
gain_fiterr->Fill( sqrt( pow(emypar[2],2) + pow(fit->GetParError(1),2) ) );
gain_fit_2d->Fill( i%NC, i/NC, mypar[2] - fit->GetParameter(1) );
gain_fiterr_2d->Fill( i%NC, i/NC, sqrt( pow(emypar[2],2) + pow(fit->GetParError(1),2) ) );
gain_ADUper1keV_2d->Fill( i%NC, i/NC, ( mypar[2] - fit->GetParameter(1) ) / 8.0 );
gainerr_ADUper1keV_2d->Fill( i%NC, i/NC, sqrt( pow( emypar[2],2 ) + pow( fit->GetParError(1),2) ) / 8.0 );
if (isEdge(i)) {
gain_fit_isEdge->Fill( mypar[2] - fit->GetParameter(1) );
}
if (isInnerEdge(i)) {
gain_fit_isInnerEdge->Fill( mypar[2] - fit->GetParameter(1) );
}
if (isDouble(i)) {
gain_fit_isDouble->Fill( mypar[2] - fit->GetParameter(1) );
}
if (isNextToDouble(i)) {
gain_fit_isNextToDouble->Fill( mypar[2] - fit->GetParameter(1) );
}
if (isQuad(i)) {
gain_fit_isQuad->Fill( mypar[2] - fit->GetParameter(1) );
}
if (isBulk(i)) {
gain_fit_isBulk->Fill( mypar[2] - fit->GetParameter(1) );
}
delete thiscalibration;
delete proj_noise;
delete proj_peak;
delete proj;
} else {
pixel_mask[i] = false;
} //end if can be fitted
} // end for loop pixel numbers (i)
} // end for loop sensor slices (j)
sprintf( savename, "%s/pixelmask_afterfit_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
pixelMaskObject->plotPixelMask( pixel_mask, savename );
TPaveText *pave = new TPaveText( 0.86, 0.95, 0.91, 0.98, "blNDC" );
pave->SetBorderSize(0);
pave->SetFillStyle(0);
pave->SetTextSize(0.06);
pave->SetTextAlign(32);
mapcanvas->cd();
fit_par3_2d->GetXaxis()->SetTitle("Column");
fit_par3_2d->GetYaxis()->SetTitle("Row");
fit_par3_2d->GetYaxis()->SetTitleOffset(0.7);
fit_par3_2d->Draw("colz");
fit_par3_2d->GetZaxis()->SetRangeUser(0,50);
sprintf( savename, "%s/fit_par3_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
fit_par3_2d->Delete();
fit_par4_2d->GetXaxis()->SetTitle("Column");
fit_par4_2d->GetYaxis()->SetTitle("Row");
fit_par4_2d->GetYaxis()->SetTitleOffset(0.7);
fit_par4_2d->Draw("colz");
fit_par4_2d->GetZaxis()->SetRangeUser(0,500);
sprintf( savename, "%s/fit_par4_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
fit_par4_2d->Delete();
fit_par5_2d->GetXaxis()->SetTitle("Column");
fit_par5_2d->GetYaxis()->SetTitle("Row");
fit_par5_2d->GetYaxis()->SetTitleOffset(0.7);
fit_par5_2d->Draw("colz");
fit_par5_2d->GetZaxis()->SetRangeUser(0,0.5);
sprintf( savename, "%s/fit_par5_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
fit_par5_2d->Delete();
fit_par6_2d->GetXaxis()->SetTitle("Column");
fit_par6_2d->GetYaxis()->SetTitle("Row");
fit_par6_2d->GetYaxis()->SetTitleOffset(0.7);
fit_par6_2d->Draw("colz");
fit_par6_2d->GetZaxis()->SetRangeUser(1.0,1.3); //was (1.0,1.25)
sprintf( savename, "%s/fit_par6_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
fit_par6_2d->Delete();
fit_par7_2d->GetXaxis()->SetTitle("Column");
fit_par7_2d->GetYaxis()->SetTitle("Row");
fit_par7_2d->GetYaxis()->SetTitleOffset(0.7);
fit_par7_2d->Draw("colz");
fit_par7_2d->GetZaxis()->SetRangeUser(0.,0.4);
sprintf( savename, "%s/fit_par7_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
peak_fit_pos_2d->GetXaxis()->SetTitle("Column");
peak_fit_pos_2d->GetYaxis()->SetTitle("Row");
peak_fit_pos_2d->GetYaxis()->SetTitleOffset(0.7);
peak_fit_pos_2d->Draw("colz");
peak_fit_pos_2d->GetZaxis()->SetRangeUser( low_ADU_peak, high_ADU_peak );
sprintf( savename, "%s/peak_fit_pos_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
peak_fit_pos_2d->Delete();
peak_fit_poserr_2d->GetXaxis()->SetTitle("Column");
peak_fit_poserr_2d->GetYaxis()->SetTitle("Row");
peak_fit_poserr_2d->GetYaxis()->SetTitleOffset(0.7);
peak_fit_poserr_2d->Draw("colz");
peak_fit_poserr_2d->GetZaxis()->SetRangeUser(0,2);
sprintf( savename, "%s/peak_fit_poserr_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
peak_fit_poserr_2d->Delete();
noise_fit_pos_2d->GetXaxis()->SetTitle("Column");
noise_fit_pos_2d->GetYaxis()->SetTitle("Row");
noise_fit_pos_2d->GetYaxis()->SetTitleOffset(0.7);
noise_fit_pos_2d->Draw("colz");
noise_fit_pos_2d->GetZaxis()->SetRangeUser(-5,5);
sprintf( savename, "%s/noise_fit_pos_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
noise_fit_pos_2d->Delete();
noise_fit_poserr_2d->GetXaxis()->SetTitle("Column");
noise_fit_poserr_2d->GetYaxis()->SetTitle("Row");
noise_fit_poserr_2d->GetYaxis()->SetTitleOffset(0.7);
noise_fit_poserr_2d->Draw("colz");
if (gain_str == "HG0") {
noise_fit_poserr_2d->GetZaxis()->SetRangeUser(0,0.1);
} else if (gain_str == "G0") {
noise_fit_poserr_2d->GetZaxis()->SetRangeUser(0,0.1); // was 0,0.05
}
sprintf( savename, "%s/noise_fit_poserr_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
noise_fit_poserr_2d->Delete();
gain_fit_2d->GetXaxis()->SetTitle("Column");
gain_fit_2d->GetYaxis()->SetTitle("Row");
gain_fit_2d->GetYaxis()->SetTitleOffset(0.7);
gain_fit_2d->Draw("colz");
sprintf( savename, "%s [ADU/8 keV]", gain_str.c_str() );
pave->AddText( (const char *)(savename) );
pave->Draw();
gain_fit_2d->GetZaxis()->SetRangeUser( low_ADU_peak, high_ADU_peak );
sprintf( savename, "%s/gain_fit_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
gain_fit_2d->Delete();
gain_fiterr_2d->GetXaxis()->SetTitle("Column");
gain_fiterr_2d->GetYaxis()->SetTitle("Row");
gain_fiterr_2d->GetYaxis()->SetTitleOffset(0.7);
gain_fiterr_2d->Draw("colz");
gain_fiterr_2d->GetZaxis()->SetRangeUser(0,2);
sprintf( savename, "%s/gain_fiterr_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
gain_fiterr_2d->Delete();
gain_ADUper1keV_2d->GetXaxis()->SetTitle("Column");
gain_ADUper1keV_2d->GetYaxis()->SetTitle("Row");
gain_ADUper1keV_2d->GetYaxis()->SetTitleOffset(0.7);
gain_ADUper1keV_2d->Draw("colz");
if (gain_str == "HG0") {
gain_ADUper1keV_2d->GetZaxis()->SetRangeUser(80,120);
} else if (gain_str == "G0") {
gain_ADUper1keV_2d->GetZaxis()->SetRangeUser(25,45); //was 35,50
}
sprintf( savename, "%s/gain_ADUper1keV_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
//gain_ADUper1keV_2d->Delete();
gainerr_ADUper1keV_2d->GetXaxis()->SetTitle("Column");
gainerr_ADUper1keV_2d->GetYaxis()->SetTitle("Row");
gainerr_ADUper1keV_2d->GetYaxis()->SetTitleOffset(0.7);
gainerr_ADUper1keV_2d->Draw("colz");
if (gain_str == "HG0") {
gainerr_ADUper1keV_2d->GetZaxis()->SetRangeUser(0,0.5);
} else if (gain_str == "G0") {
gainerr_ADUper1keV_2d->GetZaxis()->SetRangeUser(0,0.25);
}
sprintf( savename, "%s/gainerr_ADUper1keV_2d_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
mapcanvas->SaveAs( (const char *)(savename) );
//gainerr_ADUper1keV_2d[this_storagecell]->Delete();
c1->cd();
fit_par3->GetXaxis()->SetTitle("Fit par 3");
fit_par3->Draw();
sprintf( savename, "%s/fit_par3_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
fit_par3->Delete();
fit_par4->GetXaxis()->SetTitle("Fit par 4");
fit_par4->Draw();
sprintf( savename, "%s/fit_par4_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
fit_par4->Delete();
fit_par5->GetXaxis()->SetTitle("Fit par 5");
fit_par5->Draw();
sprintf( savename, "%s/fit_par5_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
fit_par5->Delete();
fit_par6->GetXaxis()->SetTitle("Fit par 6");
fit_par6->Draw();
sprintf( savename, "%s/fit_par6_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
fit_par6->Delete();
fit_par7->GetXaxis()->SetTitle("Fit par 7");
fit_par7->Draw();
sprintf( savename, "%s/fit_par7_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
peak_fit_pos->GetXaxis()->SetTitle("Peak position [ADU]");
peak_fit_pos->Draw();
sprintf( savename, "%s/peak_fit_pos_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
peak_fit_pos->Delete();
peak_fit_poserr->GetXaxis()->SetTitle("Peak position uncert [ADU]");
peak_fit_poserr->Draw();
sprintf( savename, "%s/peak_fit_poserr_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
peak_fit_poserr->Delete();
noise_fit_pos->GetXaxis()->SetTitle("Noise position [ADU]");
noise_fit_pos->Draw();
sprintf( savename, "%s/noise_fit_pos_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
noise_fit_pos->Delete();
noise_fit_poserr->GetXaxis()->SetTitle("Noise position uncert [ADU]");
noise_fit_poserr->Draw();
sprintf( savename, "%s/noise_fit_poserr_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
noise_fit_poserr->Delete();
sprintf( savename, "Gain %s [ADU / 8 keV]", gain_str.c_str() );
gain_fit->GetXaxis()->SetTitle( (const char *)(savename) );
gain_fit->Draw();
sprintf( savename, "%s/gain_fit_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
//gain_fit->Delete();
gain_fit->GetXaxis()->SetRangeUser( low_ADU_peak+30, high_ADU_peak );
gain_fit->Fit("gaus");
gain_fit->Draw();
c1->Update();
TPaveText *pave2 = new TPaveText( 0.6, 0.8, 0.94, 0.94, "blNDC" );
pave2->SetBorderSize(0);
pave2->SetFillStyle(0);
pave2->SetTextSize(0.04);
pave2->SetTextAlign(32);
TF1* gain_fit_gaus = gain_fit->GetFunction("gaus");
sprintf( savename, "Mean %0.2f #pm %0.2f", gain_fit_gaus->GetParameter(1), gain_fit_gaus->GetParError(1) );
pave2->AddText( (const char *)(savename) );
sprintf( savename, "Sigma %0.2f #pm %0.2f", gain_fit_gaus->GetParameter(2), gain_fit_gaus->GetParError(2) );
pave2->AddText( (const char *)(savename) );
pave2->Draw();
gain_fit->SetStats(kFALSE);
sprintf( savename, "%s/gain_fit_fit_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
gain_fit->Delete();
gain_fiterr->GetXaxis()->SetTitle("Gain uncert [ADU / 8 keV]");
gain_fiterr->Draw();
sprintf( savename, "%s/gain_fiterr_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
gain_fiterr->Delete();
gain_fit_isEdge->SetLineColor(kBlue);
gain_fit_isInnerEdge->SetLineColor(kCyan);
gain_fit_isDouble->SetLineColor(kGreen+2);
gain_fit_isNextToDouble->SetLineColor(kRed);
gain_fit_isQuad->SetLineColor(kOrange);
gain_fit_isEdge->Scale( 1./gain_fit_isEdge->GetEntries() );
gain_fit_isInnerEdge->Scale( 1./gain_fit_isInnerEdge->GetEntries() );
gain_fit_isDouble->Scale( 1./gain_fit_isDouble->GetEntries() );
gain_fit_isNextToDouble->Scale( 1./gain_fit_isNextToDouble->GetEntries() );
gain_fit_isQuad->Scale( 1./gain_fit_isQuad->GetEntries() );
gain_fit_isBulk->Scale( 1./gain_fit_isBulk->GetEntries() );
TLegend *leg = new TLegend( 0.62, 0.6, 0.93, 0.93 );
leg->AddEntry( gain_fit_isBulk, "Normal", "l" );
leg->AddEntry( gain_fit_isDouble, "Double", "l" );
leg->AddEntry( gain_fit_isNextToDouble, "Next to D", "l" );
leg->AddEntry( gain_fit_isEdge, "Edge", "l" );
leg->AddEntry( gain_fit_isInnerEdge, "Inner E", "l" );
sprintf( savename, "Gain %s [ADU / 8 keV]", gain_str.c_str() );
gain_fit_isDouble->GetXaxis()->SetTitle( (const char *)(savename) );
gain_fit_isDouble->GetYaxis()->SetTitle("Normalised");
gain_fit_isDouble->GetYaxis()->SetTitleOffset(1.3);
gain_fit_isDouble->SetMinimum(0.0);
gain_fit_isDouble->SetMaximum(0.16);
gain_fit_isDouble->Draw();
gain_fit_isEdge->Draw("same");
gain_fit_isInnerEdge->Draw("same");
gain_fit_isNextToDouble->Draw("same");
gain_fit_isBulk->Draw("same");
leg->Draw("same");
sprintf( savename, "%s/gain_fit_perType_sc%d_%s_M%s.png", plotfolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
c1->SaveAs( (const char *)(savename) );
gain_fit_isDouble->Delete();
gain_fit_isEdge->Delete();
gain_fit_isInnerEdge->Delete();
gain_fit_isNextToDouble->Delete();
gain_fit_isBulk->Delete();
sprintf( savename, "%s/CuFluo_gain_sc%d_%s_M%s.root", rootdatafolder, this_storagecell, gain_str.c_str(), module_str.c_str() );
TFile* saved_file = new TFile( (const char *)(savename), "RECREATE" );
gain_ADUper1keV_2d->Write();
gainerr_ADUper1keV_2d->Write();
saved_file->Close();
gain_ADUper1keV_2d->Delete();
gainerr_ADUper1keV_2d->Delete();
}