From 094e38ef126ea9fd699757a0c71fa3b41fb10f85 Mon Sep 17 00:00:00 2001 From: Benjamin Linnik Date: Wed, 31 Aug 2016 21:45:27 +0200 Subject: [PATCH] Anylyzer: Before F0 into database --- MABS_run_analyzer/ChargeSpektrum.c | 47 ++++--- MABS_run_analyzer/ChargeSpektrumFunctions.c | 143 +++++++++++++++++++- MABS_run_analyzer/HistogramType.c | 29 ++-- MABS_run_analyzer/HistogramType.h | 10 +- MABS_run_analyzer/MAPS.c | 12 +- MABS_run_analyzer/Run.c | 128 +++++++++++++----- MABS_run_analyzer/help.h | 68 ++++++++++ 7 files changed, 363 insertions(+), 74 deletions(-) diff --git a/MABS_run_analyzer/ChargeSpektrum.c b/MABS_run_analyzer/ChargeSpektrum.c index 8527083..3133e50 100644 --- a/MABS_run_analyzer/ChargeSpektrum.c +++ b/MABS_run_analyzer/ChargeSpektrum.c @@ -96,30 +96,32 @@ void ChargeSpektrum(TString runnumber = "") // compareHistogramClassVector.push_back(runs[runi]->histogramthreshold); // compareHistogramVector.push_back(runs[runi]->histogramthreshold->Seed); -// Uncomment below to do analysis without RTS pixel - compareHistogramClassVector2.push_back(runs[runi]->histogramwoRTSthreshold); -// compareHistogramVector2.push_back(runs[runi]->histogramwoRTSthreshold->Seed); - compareHistogramClassVector3.push_back(runs[runi]->histogramwoRTSthreshold->calibrated); - compareHistogramVector4.push_back(runs[runi]->histogramwoRTSthreshold->calibrated->Seed); - -// RTS analysis - runs[runi]->compareHistogramClassVector.push_back(runs[runi]->histogramthreshold); - runs[runi]->compareHistogramClassVector.push_back(runs[runi]->histogramwoRTSthreshold); - runs[runi]->compareHistogramClassVector.push_back(runs[runi]->histogramwoRTSAggresivethreshold); - plotAllRuns(&runs[runi]->compareHistogramClassVector); - runs[runi]->compareHistogramClassVector.clear(); - - // RTS Hit analysis - runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogramwoRTSthreshold->pixeltimefired, "", 0); -// runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogramwoRTSthreshold->pixeltimefiredsorted, "", 0); -// compareHistogramVector3.push_back(runs[runi]->histogramwoRTSthreshold->pixeltimefiredsorted); - compareHistogramVector6.push_back(runs[runi]->histogramwoRTSthreshold->calibrated->pixeltimefiredsorted); +// // Uncomment below to do analysis without RTS pixel +// compareHistogramClassVector2.push_back(runs[runi]->histogramwoRTSthreshold); +// // compareHistogramVector2.push_back(runs[runi]->histogramwoRTSthreshold->Seed); +// compareHistogramClassVector3.push_back(runs[runi]->histogramwoRTSthreshold->calibrated); +// compareHistogramVector4.push_back(runs[runi]->histogramwoRTSthreshold->calibrated->Seed); +// +// // RTS analysis +// runs[runi]->compareHistogramClassVector.push_back(runs[runi]->histogramthreshold); +// runs[runi]->compareHistogramClassVector.push_back(runs[runi]->histogramwoRTSthreshold); +// // runs[runi]->compareHistogramClassVector.push_back(runs[runi]->histogramwoRTSAggresivethreshold); +// plotAllRuns(&runs[runi]->compareHistogramClassVector); +// runs[runi]->compareHistogramClassVector.clear(); +// +// // RTS Hit analysis +// runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogramwoRTSthreshold->pixeltimefired, "", 0); +// // runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogramwoRTSthreshold->pixeltimefiredsorted, "", 0); +// // compareHistogramVector3.push_back(runs[runi]->histogramwoRTSthreshold->pixeltimefiredsorted); +// compareHistogramVector6.push_back(runs[runi]->histogramwoRTSthreshold->calibrated->pixeltimefiredsorted); // Leakage current -// runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogramthreshold->LeakageCurrentInPixel, "", 0); + runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogramthreshold->calibrated->LeakageCurrentInPixel, "", 0); // runs[runi]->plot1DHistogram(runs[runi]->histogram, runs[runi]->histogramthreshold->LeakageCurrentInPixelSorted, "", 0); -// compareHistogramVector5.push_back(runs[runi]->histogramthreshold->LeakageCurrentInPixelSorted); - compareHistogramVector7.push_back(runs[runi]->histogramthreshold->calibrated->LeakageCurrentInPixelSorted); + compareHistogramVector5.push_back(runs[runi]->histogramthreshold->calibrated->LeakageCurrentInPixelSorted); + compareHistogramClassVector4.push_back(runs[runi]->histogramwoRTSthreshold->calibrated); + compareHistogramClassVector5.push_back(runs[runi]->histogramwoRTSthreshold); +// compareHistogramVector7.push_back(runs[runi]->histogramthreshold->calibrated->LeakageCurrentInPixelSorted); @@ -200,6 +202,7 @@ void ChargeSpektrum(TString runnumber = "") // runs[3]->setLabel("HR18, P13, 5.0 V"); printSummaryTable(&compareHistogramClassVector); printSummaryTable(&compareHistogramClassVector2); +printSummaryTable(&compareHistogramClassVector4); CompareHistograms(&compareHistogramVector); CompareHistograms(&compareHistogramVector2); CompareHistograms(&compareHistogramVector3); @@ -213,6 +216,8 @@ plotAllRuns(&compareHistogramClassVector2); plotAllRuns(&compareHistogramClassVector3); plotAllRuns(&compareHistogramClassVector4); plotAllRuns(&compareHistogramClassVector5); +CompareLeageCurrent(&compareHistogramClassVector4); +CompareLeageCurrent(&compareHistogramClassVector5); writeObservableToFile(); // plotAllRuns("seed threshold calibrated"); // setCustomPath("Excel/"); diff --git a/MABS_run_analyzer/ChargeSpektrumFunctions.c b/MABS_run_analyzer/ChargeSpektrumFunctions.c index 58b295c..a3c1078 100644 --- a/MABS_run_analyzer/ChargeSpektrumFunctions.c +++ b/MABS_run_analyzer/ChargeSpektrumFunctions.c @@ -263,7 +263,9 @@ Bool_t CompareHistograms(vector* ptCompareHistogramVector, TString titles leg1->SetTextSize(0.025); leg1->SetFillStyle(1001); leg1->SetTextFont(132); - leg1->SetFillColor(0); leg1->SetBorderSize(0); + leg1->SetFillColor(0); + leg1->SetBorderSize(0); + TString legendEntry; TPaveText *owntitle = new TPaveText(0.15,0.9,0.930401,0.995,"nbNDC"); @@ -307,8 +309,8 @@ Bool_t CompareHistograms(vector* ptCompareHistogramVector, TString titles leg1->Draw("SAME"); curhistogramclone->SetAxisRange(0,lastbin1*1.1,"X"); - curhistogramclone->GetYaxis()->SetRangeUser(5,heighestval1*4); - gPad->SetLogy(1); +// curhistogramclone->GetYaxis()->SetRangeUser(5,heighestval1*4); +// gPad->SetLogy(1); owntitle->Clear(); owntitle->AddText(trimRunnumberAtBegin(curhistogramclone->GetName())); owntitle->Draw("SAME"); @@ -331,6 +333,141 @@ Bool_t CompareHistograms(vector* ptCompareHistogramVector, TString titles return 1; } +Bool_t CompareLeageCurrent(vector* ptCompareHistogramClassVector) +{ + if (ptCompareHistogramClassVector->size() > 0) + { + gROOT->SetStyle("RadHard_NoTitle"); + if (testifMixingCalibration(ptCompareHistogramClassVector)) return 1; + + // legend entries + Float_t height = ptCompareHistogramClassVector->size() * 0.055; + TLegend* leg1 = new TLegend(0.50,0.89-height,0.95,0.89);//(0.6,0.7,0.89,0.89); + leg1->SetTextSize(0.035); + leg1->SetFillStyle(0); + leg1->SetTextFont(132); + leg1->SetFillColor(0); + leg1->SetBorderSize(0); + TString legendEntry; + + TPaveText *owntitle = new TPaveText(0.15,0.9,0.930401,0.995,"nbNDC"); + owntitle->SetFillStyle(0); + owntitle->SetBorderSize(0); + owntitle->SetTextAlign(22); + owntitle->SetTextSize(0.05); + owntitle->SetTextFont(42); + + TString canvastitle = "SummaryLeakeCurrent"+ptCompareHistogramClassVector->at(0)->histogramdescription; + for (UInt_t histogrami=0; histogrami < ptCompareHistogramClassVector->size(); histogrami++) + { + HistogramType* curhistogramclassp = ptCompareHistogramClassVector->at(histogrami); + canvastitle+= Form("_%d",curhistogramclassp->labbook->runnumber); + } + TTimeStamp* time = new TTimeStamp(); + TString canvasname = Form("%d",time->GetNanoSec()); + TCanvas* canvas = new TCanvas(canvasname, canvastitle, 1200, 800); + + int numberofhistogramclasses = ptCompareHistogramClassVector->size(); + double lowestxvalue = 1e50; + double heighestxvalue = 0; + + +// std::vector sortedXValues; +// for (vector::iterator curHistogramClass = ptCompareHistogramClassVector->begin(); curHistogramClass != ptCompareHistogramClassVector->end(); curHistogramClass++) { +// sortedXValues.push_back((*curHistogramClass)->labbook->radDoseNonIon); +// } +// std::sort(sortedXValues.begin(),sortedXValues.end()); +// double xbins[numberofhistogramclasses+1]; +// for(int i=0; i < numberofhistogramclasses; i++){ +// xbins[i] = sortedXValues.at(i) - 0.5; +// } +// xbins[numberofhistogramclasses]=sortedXValues.at(numberofhistogramclasses-1); + + for (vector::iterator curHistogramClass = ptCompareHistogramClassVector->begin(); curHistogramClass != ptCompareHistogramClassVector->end(); curHistogramClass++) { + lowestxvalue = lowestxvalue>(*curHistogramClass)->labbook->radDoseNonIon?(*curHistogramClass)->labbook->radDoseNonIon:lowestxvalue; + heighestxvalue = heighestxvalue<(*curHistogramClass)->labbook->radDoseNonIon?(*curHistogramClass)->labbook->radDoseNonIon:heighestxvalue; +// cout << heighestxvalue << " " << lowestxvalue << " " << (*curHistogramClass)->labbook->radDoseNonIon << endl; + } +// cout << "----" << endl; + //double stepsizex = (heighestxvalue-lowestxvalue)/(numberofhistogramclasses); + double stepsizex = 1; + Int_t numberbins = (heighestxvalue-lowestxvalue)/stepsizex; +// if (stepsizex == 0) +// stepsizex = 1; +// cout << "stepsizex: " << stepsizex << endl; + double xbins[numberbins+2]; + for(int i=0; i <= numberbins+1; i++){ + xbins[i] = lowestxvalue+stepsizex*i-stepsizex/2; + } + + double heighestyvalue = 0; + for (vector::iterator curHistogramClass = ptCompareHistogramClassVector->begin(); curHistogramClass != ptCompareHistogramClassVector->end(); curHistogramClass++) { + Int_t lastbinabovezero = (*curHistogramClass)->LeakageCurrentInPixelSorted->FindLastBinAbove(0); + double curheightyvalue = (*curHistogramClass)->LeakageCurrentInPixelSorted->GetBinContent(lastbinabovezero); + heighestyvalue = heighestyvalue::iterator curHistogramClass = ptCompareHistogramClassVector->begin(); curHistogramClass != ptCompareHistogramClassVector->end(); curHistogramClass++) { + for (Int_t pixeli=0; (*curHistogramClass)->LeakageCurrentInPixelSorted->GetBinContent(pixeli) != 0; pixeli++) { +// cout << pixeli << " " << (*curHistogramClass)->LeakageCurrentInPixelSorted->GetBinContent(pixeli) << " "; + hcandle->Fill((*curHistogramClass)->labbook->radDoseNonIon,(*curHistogramClass)->LeakageCurrentInPixelSorted->GetBinContent(pixeli)); +// hcandle->Fill((*curHistogramClass)->labbook->radDoseNonIon,1.5); + } + legendEntry = Form("%f", (*curHistogramClass)->labbook->radDoseNonIon); + // leg1->AddEntry(curhistogramclone, legendEntry, "l"); + } +// hcandle->GetXaxis()->SetTitle("Non ionizing radiation damage [10^{13} n_{eq}/cm^{2}]]"); + hcandle->GetXaxis()->SetTitle("Radiation damage [10^{13} n/cm^{2}]]"); + if (ptCompareHistogramClassVector->at(0)->iscalibrated) { + hcandle->GetYaxis()->SetTitle("Charge [fA]"); + } else { + hcandle->GetYaxis()->SetTitle("Charge [ADU]"); + } + hcandle->GetXaxis()->CenterTitle(); + hcandle->GetYaxis()->CenterTitle(); + hcandle->GetYaxis()->SetNdivisions(1010); + hcandle->SetBarWidth(1); + hcandle->SetLineWidth(2); + hcandle->Draw("CANDLEX(012311)"); + gPad->SetLogy(1); + gPad->SetGridy(1); + + + owntitle->Clear(); + owntitle->AddText("Leakage current comparison"); + owntitle->Draw("SAME"); + +// leg1->Draw("SAME"); + + canvas->Update(); + MSaveBigPNG(canvas,savepathresults + "/" + canvastitle + ".png"); + + TImageDump *img = new TImageDump(savepathresults + "/" + canvastitle + ".png"); + canvas->Paint(); + img->Close(); + + TFile *f = new TFile(savepathresults + "/" + canvastitle + ".root","RECREATE"); + f->cd(); + f->Append(canvas); + //f->Append(img); + f->Write(); + + // gROOT->SetStyle("RadHard_AutoTitle"); + return 0; + } + return 1; +} + Bool_t plotAllRuns() { return plotAllRuns(&compareHistogramClassVector); diff --git a/MABS_run_analyzer/HistogramType.c b/MABS_run_analyzer/HistogramType.c index 7ca45d0..e9b4f04 100644 --- a/MABS_run_analyzer/HistogramType.c +++ b/MABS_run_analyzer/HistogramType.c @@ -50,14 +50,14 @@ void HistogramType::initHistograms(Int_t gotcolor, Int_t gotstyle) { initHistogramCustom(pixeltimefired, "Pixel fired, used for " + histogramdescription, color, style, 0, cursensorinfo->columns*cursensorinfo->rows-1, cursensorinfo->columns*cursensorinfo->rows, "Pixel index", "Counts"); initHistogramCustom(pixeltimefiredsorted, "Fire counter" + histogramdescription, color, style, 0, labbook->frames_foundDB/100-1, labbook->frames_foundDB/100, "times fired", "Counts"); initHistogramCustom(LeakageCurrentInPixel, "Leakage current per pixel" + histogramdescription, color, style, 0, cursensorinfo->columns*cursensorinfo->rows-1, cursensorinfo->columns*cursensorinfo->rows, "Pixel index", "Average CDS"); - initHistogram(LeakageCurrentInPixelSorted, "Leakage current" + histogramdescription, color, style); + initHistogramCustom(LeakageCurrentInPixelSorted, "Leakage current per pixel sorted" + histogramdescription, color, style, 0, cursensorinfo->columns*cursensorinfo->rows-1, cursensorinfo->columns*cursensorinfo->rows, "Pixel index", "Average CDS"); +// initHistogram(LeakageCurrentInPixelSorted, "Leakage current" + histogramdescription, color, style); initHistogramCustom(pixelHadGoodVeto, "Number of times pixel had good veto" + histogramdescription, color, style, 0, cursensorinfo->columns*cursensorinfo->rows-1, cursensorinfo->columns*cursensorinfo->rows, "Pixel index", "Counter good veto"); - SeedPerc->GetXaxis()->SetRangeUser(0,50); Noise->SetBins(cursystempar->nbinsnoise, 0, cursystempar->maxbinnoise); NoiseEnd->SetBins(cursystempar->nbinsnoise, 0, cursystempar->maxbinnoise); - LeakageCurrentInPixelSorted->SetBins(cursystempar->nbinsnoise*10, 0, cursystempar->maxbinnoise*30); +// LeakageCurrentInPixelSorted->SetBins(cursystempar->nbinsnoise*10, 0, cursystempar->maxbinnoise*30); } void HistogramType::initHistogram(TH1F* &histogrampointer, TString prefix, Int_t color, Int_t style) { @@ -117,8 +117,8 @@ Bool_t HistogramType::calibrateHistograms( Float_t gotgain ) { if (pixeltimefiredsorted != 0) calibrated->pixeltimefiredsorted = (TH1F*)pixeltimefiredsorted->Clone(); if (SeedPerc != 0) calibrated->SeedPerc = SeedPerc; if (histAvgCluster != 0) calibrate2DHistogramCounts(calibrated->histAvgCluster, histAvgCluster); - if (LeakageCurrentInPixel != 0)calibrated->LeakageCurrentInPixel = (TH1F*)LeakageCurrentInPixel->Clone(); - if (LeakageCurrentInPixelSorted != 0) calibrateHistogram(calibrated->LeakageCurrentInPixelSorted, LeakageCurrentInPixelSorted); + if (LeakageCurrentInPixel != 0) calibrateHistogramYAxis(calibrated->LeakageCurrentInPixel, LeakageCurrentInPixel); + if (LeakageCurrentInPixelSorted != 0) calibrateHistogramYAxis(calibrated->LeakageCurrentInPixelSorted, LeakageCurrentInPixelSorted); calibrated->posSeed = posSeed * gain; calibrated->posSum = posSum * gain; @@ -140,9 +140,9 @@ Bool_t HistogramType::calibrateHistograms( Float_t gotgain ) { calibrated->RTSpixel = RTSpixel; calibrated->percentageofRTSpixel = percentageofRTSpixel; calibrated->avgLeakageCurrentInChip = avgLeakageCurrentInChip * gain; - calibrated->medianLeakageCurrent = medianLeakageCurrent * gain; - calibrated->medianLeakageCurrentMinus = medianLeakageCurrentMinus * gain; - calibrated->medianLeakageCurrentPlus = medianLeakageCurrentPlus * gain; + calibrated->medianLeakageCurrent = medianLeakageCurrent * gain / 3.7 / pow10(-3) * (1.6*pow10(-18)) * pow10(15); + calibrated->medianLeakageCurrentMinus = medianLeakageCurrentMinus * gain / 3.7 / pow10(-3) * (1.6*pow10(-18)) * pow10(15); + calibrated->medianLeakageCurrentPlus = medianLeakageCurrentPlus * gain / 3.7 / pow10(-3) * (1.6*pow10(-18)) * pow10(15); calibrated->iscalibrated = true; @@ -164,6 +164,19 @@ void HistogramType::calibrateHistogram(TH1F* &histogrampointernew, TH1F* &histog histogrampointernew->GetYaxis()->SetTitle(Form("Entries [1/%.1f e]",histogrampointernew->GetBinWidth(1))); } + +void HistogramType::calibrateHistogramYAxis(TH1F* &histogrampointernew, TH1F* &histogrampointerold) { + histogrampointernew = (TH1F*)histogrampointerold->Clone(); + histogrampointernew->SetName(Form("%s Calibr.", histogrampointerold->GetName())); + histogrampointernew->SetTitle(Form("%s Calibr.", histogrampointerold->GetTitle())); + histogrampointernew->GetYaxis()->SetTitle("Q_coll [e]"); + int nbins = histogrampointernew->GetXaxis()->GetNbins(); + for(int x=0; x <= nbins; x++){ + histogrampointernew->SetBinContent(x,histogrampointerold->GetBinContent(x)*gain); + } + histogrampointernew->GetYaxis()->SetTitle(Form("Entries [1/%.1f e]",histogrampointernew->GetBinWidth(1))); +} + void HistogramType::calibrate2DHistogramCounts(TH2F* &histogrampointernew, TH2F* &histogrampointerold ) { histogrampointernew = (TH2F*)histogrampointerold->Clone(); histogrampointernew->SetName(Form("%s Calibr.", histogrampointerold->GetName())); diff --git a/MABS_run_analyzer/HistogramType.h b/MABS_run_analyzer/HistogramType.h index 6291ce7..ff71db1 100644 --- a/MABS_run_analyzer/HistogramType.h +++ b/MABS_run_analyzer/HistogramType.h @@ -68,9 +68,13 @@ private: // PRIVATE METHODS APPLYABLE TO HISTOGRAMS //***************** /** - * @brief rescales one specific histogram from ADU to electrons */ + * @brief rescales one specific histogram from ADU to electrons along the x axis */ void calibrateHistogram(TH1F* &histogrampointernew, TH1F* &histogrampointerold); + /** + * @brief rescales one specific histogram from ADU to electrons along the y axis */ + void calibrateHistogramYAxis(TH1F* &histogrampointernew, TH1F* &histogrampointerold); + /** * @brief normalizes one specific histogram with #frames_found along the y axis, so it can be compared graphically */ void normalizeHistogram(TH1F* &histogrampointernew, TH1F* &histogrampointerold); @@ -162,9 +166,9 @@ public: /// Average Noise value Float_t medianLeakageCurrent = 0; - /// Positive Noise sigma + /// Positive Noise sigma/2 Float_t medianLeakageCurrentPlus = 0; - /// Negative Noise sigma + /// Negative Noise sigma/2 Float_t medianLeakageCurrentMinus = 0; /// Integral value, after integrating from #noisethresholdborder to maxbin. diff --git a/MABS_run_analyzer/MAPS.c b/MABS_run_analyzer/MAPS.c index 6ebc066..0ff60f3 100644 --- a/MABS_run_analyzer/MAPS.c +++ b/MABS_run_analyzer/MAPS.c @@ -69,6 +69,7 @@ Bool_t MAPS::initNewRootFile() { fNoiseTree->Branch("pixelnumber", &fNoiseInfo.pixelnumber , "pixelnumber/i" , 32000); fNoiseTree->Branch("noise" , &fNoiseInfo.fNoise[0] , "noise[pixelnumber]/F" , 32000); fNoiseTree->Branch("pedestal" , &fNoiseInfo.fPedestals[0] , "pedestal[pixelnumber]/F" , 32000); + fNoiseTree->Branch("F0" , &fNoiseInfo.fF0[0] , "F0[pixelnumber]/i" , 32000); initHistograms(); cout<<"-----------------------"<SetBranchAddress("frame" , &fFrameInfo.frame); fNoiseTree->SetBranchAddress("noise" , &fNoiseInfo.fNoise[0]); fNoiseTree->SetBranchAddress("pedestal" , &fNoiseInfo.fPedestals[0]); + fNoiseTree->SetBranchAddress("F0" , &fNoiseInfo.fF0[0]); initHistograms(); cout<<"-----------------------"< Pedestal suspiciously high!"< Noise suspiciously high!"<1000 && fFrameNumber< 1200 && RefillNoiseBranch) + if (fFrameNumber%(Frames) == 0 && RefillNoiseBranch && fFrameNumber >= 500) // Don't use first 500 frames { if(fSave) { diff --git a/MABS_run_analyzer/Run.c b/MABS_run_analyzer/Run.c index 43d0e05..e469c01 100644 --- a/MABS_run_analyzer/Run.c +++ b/MABS_run_analyzer/Run.c @@ -174,14 +174,14 @@ Run::Run(Int_t runnumber, Int_t loopi, TString customtitle) HistogramClassVector.push_back(histogramwoRTSthreshold); - histogramwoRTSAggresive = new HistogramType(" more RTS cleaned", &cursystemparam, &cursensorinfo, humanreadablestr, &labbook, rootcolors[plotStyle], rootlinestyle[plotStyle] ); - histogramwoRTSAggresive->maskRTSpixel = true; - histogramwoRTSAggresive->RTSthreshold = 1.5; - HistogramClassVector.push_back(histogramwoRTSAggresive); - histogramwoRTSAggresivethreshold = new HistogramType(" Threshold, more RTS cleaned", &cursystemparam, &cursensorinfo, humanreadablestr, &labbook, rootcolors[plotStyle], rootlinestyle[plotStyle] ); - histogramwoRTSAggresivethreshold->maskRTSpixel = true; - histogramwoRTSAggresivethreshold->RTSthreshold = 1.0; - HistogramClassVector.push_back(histogramwoRTSAggresivethreshold); +// histogramwoRTSAggresive = new HistogramType(" more RTS cleaned", &cursystemparam, &cursensorinfo, humanreadablestr, &labbook, rootcolors[plotStyle], rootlinestyle[plotStyle] ); +// histogramwoRTSAggresive->maskRTSpixel = true; +// histogramwoRTSAggresive->RTSthreshold = 1.5; +// HistogramClassVector.push_back(histogramwoRTSAggresive); +// histogramwoRTSAggresivethreshold = new HistogramType(" Threshold, more RTS cleaned", &cursystemparam, &cursensorinfo, humanreadablestr, &labbook, rootcolors[plotStyle], rootlinestyle[plotStyle] ); +// histogramwoRTSAggresivethreshold->maskRTSpixel = true; +// histogramwoRTSAggresivethreshold->RTSthreshold = 1.0; +// HistogramClassVector.push_back(histogramwoRTSAggresivethreshold); // histogram with pixel, which have a good veto spectrum @@ -396,7 +396,7 @@ Bool_t Run::analyzeRun(Bool_t force) processed->InitialDynNoise(); int start = 0; int nframes = processed->GetNumberFrames(); -// for(int i=0; i<5000;i++) // TODO remove 100000 run 342272 +// for(int i=0; i<200;i++) // TODO remove 100000 run 342272 for(int i=0; ipixeltimefiredsorted); - Double_t const probabilities[] = {1-erf(HistogramTypepointer->RTSthreshold/sqrt(2)), 0.5, erf(HistogramTypepointer->RTSthreshold/sqrt(2))}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17}; - HistogramTypepointer->pixeltimefiredsorted->GetQuantiles( 3, leakagequantiles, probabilities); - HistogramTypepointer->medianLeakageCurrent = leakagequantiles[1]; - HistogramTypepointer->medianLeakageCurrentPlus = leakagequantiles[2] - leakagequantiles[1]; - HistogramTypepointer->medianLeakageCurrentMinus = leakagequantiles[1] - leakagequantiles[0]; +// Double_t const probabilities[] = {1-erf(HistogramTypepointer->RTSthreshold/sqrt(2)), 0.5, erf(HistogramTypepointer->RTSthreshold/sqrt(2))}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17}; +// HistogramTypepointer->pixeltimefiredsorted->GetQuantiles( 3, leakagequantiles, probabilities); +// HistogramTypepointer->medianLeakageCurrent = leakagequantiles[1]; +// HistogramTypepointer->medianLeakageCurrentPlus = leakagequantiles[2] - leakagequantiles[1]; +// HistogramTypepointer->medianLeakageCurrentMinus = leakagequantiles[1] - leakagequantiles[0]; // cout << "Median is at: " << HistogramTypepointer->medianLeakageCurrent << endl; // cout << "Upper limit is at: " << leakagequantiles[2] << endl; // cout << "Lower limit is at: " << leakagequantiles[0] << endl; @@ -925,7 +925,7 @@ void Run::constructUpdateString(string *sqlupdatequery, const string databaseval // cout << colorred << databasevaluename << " : " << value << endlr; if (!std::isinf(value)) { - if (value>0) + if (abs(value)>0.0) { if ((*sqlupdatequery).length() > 0) *sqlupdatequery+= ", "; @@ -945,18 +945,18 @@ void Run::updateDatabase() { constructUpdateString(&sqlupdatequery, "VetoIntegral", histogramthreshold->normalized->integralVeto, 10); constructUpdateString(&sqlupdatequery, "SumIntegral", histogramthreshold->normalized->integralSum, 10); if (histogramthreshold->normalized->calibrated != 0) - if (histogramthreshold->normalized->calibrated->avgNoise < 100) +// if (histogramthreshold->normalized->calibrated->avgNoise < 100) constructUpdateString(&sqlupdatequery, "Avg.Noise", histogramthreshold->normalized->calibrated->avgNoise); if (histogramthreshold->normalized->calibrated != 0) - if (histogramthreshold->normalized->calibrated->avgNoisePlus < 100) +// if (histogramthreshold->normalized->calibrated->avgNoisePlus < 100) constructUpdateString(&sqlupdatequery, "Avg.Noise+", histogram->normalized->calibrated->avgNoisePlus, 2); if (histogramthreshold->normalized->calibrated != 0) - if (histogramthreshold->normalized->calibrated->avgNoiseMinus < 100) +// if (histogramthreshold->normalized->calibrated->avgNoiseMinus < 100) constructUpdateString(&sqlupdatequery, "Avg.Noise-", histogramthreshold->normalized->calibrated->avgNoiseMinus, 2); constructUpdateString(&sqlupdatequery, "CCE_1", histogramthreshold->normalized->CCE_in_Perc_1); constructUpdateString(&sqlupdatequery, "CCE_25", histogramthreshold->normalized->CCE_in_Perc_25); constructUpdateString(&sqlupdatequery, "StoN", histogramthreshold->normalized->StoN, 3); - if (histogramthreshold->normalized->avgNoise < 100) +// if (histogramthreshold->normalized->avgNoise < 100) constructUpdateString(&sqlupdatequery, "Avg.NoiseADC", histogramthreshold->normalized->avgNoise); constructUpdateString(&sqlupdatequery, "Frames_found", frames_found,100000000); if (histogramwoRTSthreshold != 0) { @@ -965,8 +965,11 @@ void Run::updateDatabase() { constructUpdateString(&sqlupdatequery, "LeakageCurMedADC", histogramwoRTSthreshold->medianLeakageCurrent); constructUpdateString(&sqlupdatequery, "LeakageCurAvgADC", histogramwoRTSthreshold->avgLeakageCurrentInChip); if (histogramwoRTSthreshold->calibrated != 0) { - constructUpdateString(&sqlupdatequery, "LeakageCurMedE", histogramwoRTSthreshold->calibrated->medianLeakageCurrent); - constructUpdateString(&sqlupdatequery, "LeakageCurAvgE", histogramwoRTSthreshold->calibrated->avgLeakageCurrentInChip); + constructUpdateString(&sqlupdatequery, "LeakageCurAvgE", histogramwoRTSthreshold->calibrated->avgLeakageCurrentInChip); + constructUpdateString(&sqlupdatequery, "LeakageCurfA", histogramwoRTSthreshold->calibrated->medianLeakageCurrent); + constructUpdateString(&sqlupdatequery, "LeakageCurfA+", histogramwoRTSthreshold->calibrated->medianLeakageCurrentPlus); + constructUpdateString(&sqlupdatequery, "LeakageCurfA-", histogramwoRTSthreshold->calibrated->medianLeakageCurrentMinus); + constructUpdateString(&sqlupdatequery, "CalibrationPeak", Fe55run.posVeto); } } if (histogramthreshold->normalized->calibrated != 0) @@ -1001,23 +1004,26 @@ string Run::to_str_w_prec(const Float_t a_value, int precision = 3) Bool_t Run::binLeakageCurrent(HistogramType* oneHistogramClass) { - Double_t const probabilities[] = {0.3415/2, 0.5, 1-0.3415/2}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17}; + Double_t const probabilities[] = {0.158655254, 0.5, 1-0.158655254}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17}; +// Double_t const probabilities[] = {0.25, 0.5, 0.75}; // sigma/2 from gaus to the left and to the right //{0.17, 0.5, 1-0.17}; Double_t pedestals [cursensorinfo.columns*cursensorinfo.rows]; - for (Int_t pixeli=0; pixelifNoiseTree->GetEntries(); framei++) { // loop over all frames processed->fNoiseTree->GetEntry(framei); - for (Int_t pixeli=0; pixelifNoiseInfo.fPedestals[pixeli]; + for (Int_t pixeli=0; pixelifNoiseInfo.fPedestals[pixeli]; } } Bool_t RTSpixel =false; u_int numberofconsideredpixel=0; u_int poscounter = 0; - for (u_int pixeli=0; (int)pixeliRTSpixel.size(); RTSpixeli++) { + oneHistogramClass->LeakageCurrentInPixel->SetBinContent(pixeli, 0); + for (u_int RTSpixeli=poscounter; RTSpixeli < oneHistogramClass->RTSpixel.size(); RTSpixeli++) { // loop over all RTS pixel if (pixeli == oneHistogramClass->RTSpixel[RTSpixeli]) { poscounter = RTSpixeli; RTSpixel = true; @@ -1026,7 +1032,6 @@ Bool_t Run::binLeakageCurrent(HistogramType* oneHistogramClass) pedestals[pixeli] /= processed->fNoiseTree->GetEntries(); if (!RTSpixel) { oneHistogramClass->LeakageCurrentInPixel->SetBinContent(pixeli, pedestals[pixeli]); - oneHistogramClass->LeakageCurrentInPixelSorted->Fill(pedestals[pixeli]); oneHistogramClass->avgLeakageCurrentInChip += pedestals[pixeli]; numberofconsideredpixel++; } @@ -1034,20 +1039,69 @@ Bool_t Run::binLeakageCurrent(HistogramType* oneHistogramClass) { // cout << colorcyan << "Excluded: " << pixeli << endlr; } - } -// for (u_int pixeli=0; (int)pixeli<111 ; pixeli++) { -// cout << pedestals[pixeli] << " "; + } +// for (u_int pixeli=0; (int)pixelifNoiseTree->GetEntries() << endl; - oneHistogramClass->avgLeakageCurrentInChip /= numberofconsideredpixel; - oneHistogramClass->LeakageCurrentInPixelSorted->GetQuantiles( 3, leakagequantiles, probabilities); - oneHistogramClass->medianLeakageCurrent = leakagequantiles[1]; - oneHistogramClass->medianLeakageCurrentPlus = leakagequantiles[2] - leakagequantiles[1]; - oneHistogramClass->medianLeakageCurrentMinus = leakagequantiles[1] - leakagequantiles[0]; +// TThread *th1, *th2, *th3, *th4; +// typedef struct { +// char name[64]; +// TRandom *r; +// Long64_t n; +// Int_t color; +// } args_t; +// args_t args1, args2, args3, args4; +// th1 = new TThread("th1", sortThread, (void*) &args1); +// th1->Run(); +// + std::vector sortedHistogram; + for (Int_t pixeli=0; pixeli < oneHistogramClass->LeakageCurrentInPixel->GetNbinsX(); pixeli++) { + if (oneHistogramClass->LeakageCurrentInPixel->GetBinContent(pixeli) != 0) // != 0 TODO + sortedHistogram.push_back(oneHistogramClass->LeakageCurrentInPixel->GetBinContent(pixeli)); +// cout << oneHistogramClass->LeakageCurrentInPixel->GetBinContent(pixeli) << " " << endl; + } + std::sort(sortedHistogram.begin(),sortedHistogram.end()); +// oneHistogramClass->LeakageCurrentInPixelSorted->SetBins(sortedHistogram.size(), 0, sortedHistogram.size()); + for (UInt_t pixeli=0; pixeli < sortedHistogram.size(); pixeli++) { + oneHistogramClass->LeakageCurrentInPixelSorted->SetBinContent(pixeli, sortedHistogram.at(pixeli)); + } + + oneHistogramClass->avgLeakageCurrentInChip /= numberofconsideredpixel; + cout << "avgLeakageCurrentInChip: " << oneHistogramClass->avgLeakageCurrentInChip << endl; + +// oneHistogramClass->LeakageCurrentInPixelSorted->GetQuantiles( 3, leakagequantiles, probabilities); + Int_t lastbinabovezero = oneHistogramClass->LeakageCurrentInPixelSorted->FindLastBinAbove(0); +// cout << "lastbinabovezero: " << lastbinabovezero << endl; + oneHistogramClass->medianLeakageCurrent = oneHistogramClass->LeakageCurrentInPixelSorted->GetBinContent(lastbinabovezero*probabilities[1]); + oneHistogramClass->medianLeakageCurrentPlus = oneHistogramClass->LeakageCurrentInPixelSorted->GetBinContent(lastbinabovezero*probabilities[2]); + oneHistogramClass->medianLeakageCurrentPlus = oneHistogramClass->medianLeakageCurrentPlus - oneHistogramClass->medianLeakageCurrent; + oneHistogramClass->medianLeakageCurrentMinus = oneHistogramClass->LeakageCurrentInPixelSorted->GetBinContent(lastbinabovezero*probabilities[0]); + oneHistogramClass->medianLeakageCurrentMinus = oneHistogramClass->medianLeakageCurrentMinus - oneHistogramClass->medianLeakageCurrent; + cout << "medianLeakageCurrent: " << oneHistogramClass->medianLeakageCurrent << endl; + cout << "medianLeakageCurrentPlus: " << oneHistogramClass->medianLeakageCurrentPlus << endl; + cout << "medianLeakageCurrentMinus: " << oneHistogramClass->medianLeakageCurrentMinus << endl; +// +// +// oneHistogramClass->medianLeakageCurrent = leakagequantiles[1]; +// oneHistogramClass->medianLeakageCurrentPlus = leakagequantiles[2]; +// oneHistogramClass->medianLeakageCurrentMinus = leakagequantiles[0]; +// cout << "medianLeakageCurrent: " << oneHistogramClass->medianLeakageCurrent << endl; +// cout << "medianLeakageCurrentPlus: " << oneHistogramClass->medianLeakageCurrentPlus << endl; +// cout << "medianLeakageCurrentMinus: " << oneHistogramClass->medianLeakageCurrentMinus << endl; // cout << colorcyan << "oneHistogramClass->avgLeakageCurrentInChip: " << oneHistogramClass->avgLeakageCurrentInChip << endlr; // cout << colorcyan << "oneHistogramClass->medianLeakageCurrent: " << oneHistogramClass->medianLeakageCurrent << endlr; + TCanvas *c1 = new TCanvas("c1","c1",600,400); + TH2F *hcandle = new TH2F("hcandle","Option CANDLE6 example ",40,-4,4,40,-20,20); + Float_t px, py; + for (Int_t pixeli=0; pixeli < oneHistogramClass->LeakageCurrentInPixel->GetNbinsX(); pixeli++) { + hcandle->Fill(1,oneHistogramClass->LeakageCurrentInPixel->GetBinContent(pixeli)); + } + hcandle->SetMarkerSize(0.5); + hcandle->SetBarWidth(1.0); + hcandle->Draw("CANDLEX5"); +// gPad->BuildLegend(0.6,0.7,0.7,0.8); return 0; } diff --git a/MABS_run_analyzer/help.h b/MABS_run_analyzer/help.h index e6a1ae5..bf71c36 100644 --- a/MABS_run_analyzer/help.h +++ b/MABS_run_analyzer/help.h @@ -10,6 +10,8 @@ #include #include #include +#include // std::sort + #include "TStopwatch.h" #include @@ -18,6 +20,8 @@ #include #include +#include "TThread.h" + #define MAXHITS 1000000 #define MAXPIXELS 100000 //#################################################################### @@ -146,6 +150,8 @@ struct noiseInfo{ Float_t* fNoise; //[MAXPIXELS]; /// Array with pedestial information of given pixel Float_t* fPedestals; //[MAXPIXELS]; + // Array with F0 information of given pixel + Int_t* fF0; }; @@ -513,5 +519,67 @@ TH1F * ShowOverflow(TH1F *h) return htmp; } +bool compareFloatvalues (Float_t i,Float_t j) { return (iSetMarkerStyle(20); +// gr->SetMarkerSize(0.7); +// gr->SetMarkerColor(args->color); +// gr->SetLineColor(args->color); +// gr->SetLineWidth(2); +// +// Int_t k = 0; +// Double_t diffpi; +// Long64_t npi = 0; +// Double_t pi = TMath::Pi(); +// const Int_t NR = 20000; +// const Int_t NR2 = NR/2; +// Double_t rn[NR]; +// Long64_t i = 0; +// while (i<=args->n) { +// i += NR2; +// args->r->RndmArray(NR,rn); +// for (Int_t j=0;jSetPoint(k,i,diffpi); +// if (k ==0) { +// TThread::Lock(); +// gr->Draw("lp"); +// TThread::UnLock(); +// } +// c1->Modified(); +// c1->Update(); +// k++; +// } +// // yield execution to another thread that is ready to run +// // MSDN : +// // Sleep(0) causes the thread to relinquish the remainder of its +// // time slice to any other thread of equal priority that is ready to run. +// // If there are no other threads of equal priority ready to run, +// // the function returns immediately, and the thread continues execution. +// gSystem->Sleep(0); +// } +// timer.Stop(); +// Double_t cpu = timer.CpuTime(); +// cputot += cpu; +// Double_t nanos = 1.e9*cpu/Double_t(args->n); +// legend->AddEntry(gr,Form("%-14s: %6.1f ns/call",args->name,nanos),"lp"); +// c1->Modified(); +// c1->Update(); +// TThread::Printf("RANDOM = %s : RT=%7.3f s, Cpu=%7.3f s\n",args->name,timer.RealTime(),cpu); + return 0; +} + //#################################################################### #endif \ No newline at end of file -- 2.43.0