FEDRA emulsion software from the OPERA Collaboration
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
EdbShowPVRQuality Class Reference

#include <EdbShowPVRQuality.h>

Inheritance diagram for EdbShowPVRQuality:
Collaboration diagram for EdbShowPVRQuality:

Public Member Functions

void CheckEdbPVRec ()
 
void CheckFilledXYSize ()
 
Bool_t CheckSegmentQualityInPattern_ConstBTDens (EdbPVRec *ali, Int_t PatternAtNr, EdbSegP *seg)
 
Bool_t CheckSegmentQualityInPattern_ConstQual (EdbPVRec *ali, Int_t PatternAtNr, EdbSegP *seg)
 
 ClassDef (EdbShowPVRQuality, 1)
 
void CreateEdbPVRec ()
 
 EdbShowPVRQuality ()
 
 EdbShowPVRQuality (EdbPVRec *ali)
 
 EdbShowPVRQuality (EdbPVRec *ali, Float_t BTDensityTargetLevel)
 
void Execute_ConstantBTDensity ()
 
void Execute_ConstantQuality ()
 
Int_t FindFirstBinAbove (TH1 *hist, Double_t threshold, Int_t axis)
 
Int_t FindLastBinAbove (TH1 *hist, Double_t threshold, Int_t axis)
 
Float_t * GetagreementChi2Cut ()
 
Float_t GetagreementChi2Cut (Int_t patNR)
 
Float_t GetagreementChi2CutMeanChi2 ()
 
Float_t GetagreementChi2CutMeanW ()
 
Float_t GetagreementChi2CutRMSChi2 ()
 
Float_t GetagreementChi2CutRMSW ()
 
Float_t GetBTDensityLevel ()
 
Bool_t GetBTDensityLevelCalcMethodMC ()
 
Int_t GetBTDensityLevelCalcMethodMCConfirmation ()
 
Int_t GetCutMethod ()
 
Bool_t GetCutMethodIsDone (Int_t type)
 
Float_t * GetCutp0 ()
 
Float_t GetCutp0 (Int_t patNR)
 
Float_t * GetCutp1 ()
 
Float_t GetCutp1 (Int_t patNR)
 
EdbPVRecGetEdbPVRec ()
 
EdbPVRecGetEdbPVRec (Bool_t NeedModified)
 
EdbPVRecGetEdbPVRec (Int_t EdbPVRecType)
 
EdbPVRecGetEdbPVRec_modified ()
 
EdbPVRecGetEdbPVRec_orig ()
 
TH2F * GetHistChi2W ()
 
TH2F * GetHistYX ()
 
TObjArray * GetTracksFromLinkedTracksRootFile ()
 
void Help ()
 
void Print ()
 
void PrintCutType ()
 
void PrintCutType0 ()
 
void PrintCutType1 ()
 
EdbPVRecRemove_DoubleBT (EdbPVRec *aliSource)
 
EdbPVRecRemove_Passing (EdbPVRec *aliSource)
 
EdbPVRecRemove_Segment (EdbSegP *seg, EdbPVRec *aliSource=NULL, Int_t Option=0)
 
EdbPVRecRemove_SegmentArray (TObjArray *segarray, EdbPVRec *aliSource=NULL, Int_t Option=0)
 
EdbPVRecRemove_Track (EdbTrackP *track, EdbPVRec *aliSource=NULL, Int_t Option=0)
 
EdbPVRecRemove_TrackArray (TObjArray *trackArray, EdbPVRec *aliSource=NULL, Int_t Option=0)
 
void SetBTDensityLevel (Float_t BTDensityLevel)
 
void SetBTDensityLevelCalcMethodMC (Bool_t BTDensityLevelCalcMethodMC)
 
void SetBTDensityLevelCalcMethodMCConfirmation (Int_t BTDensityLevelCalcMethodMCConfirmationNumber)
 
void SetCutMethod (Int_t CutMethod)
 
void SetEdbPVRec (EdbPVRec *Ali_orig)
 
void SetHistGeometry_MC ()
 
void SetHistGeometry_OPERA ()
 
void SetHistGeometry_OPERAandMC ()
 
TObjArray * TrackArrayToSegmentArray (TObjArray *trackArray)
 
TObjArray * TrackToSegmentArray (EdbTrackP *track)
 
virtual ~EdbShowPVRQuality ()
 

Protected Member Functions

void Init ()
 
void Set0 ()
 

Private Attributes

Float_t eAgreementChi2CutMeanChi2
 
Float_t eAgreementChi2CutMeanW
 
Float_t eAgreementChi2CutRMSChi2
 
Float_t eAgreementChi2CutRMSW
 
Float_t eAgreementChi2WDistCut [114]
 
Int_t eAli_maxNpatterns
 
EdbPVReceAli_modified
 
EdbPVReceAli_orig
 
Float_t eBTDensityLevel
 
Bool_t eBTDensityLevelCalcMethodMC
 
Int_t eBTDensityLevelCalcMethodMCConfirmationNumber
 
Float_t eCutDistChi2 [114]
 
Float_t eCutDistW [114]
 
Int_t eCutMethod
 
Bool_t eCutMethodIsDone [2]
 
Float_t eCutp0 [114]
 
Float_t eCutp1 [114]
 
TH2F * eHistChi2W
 
Int_t eHistGeometry
 
TH2F * eHistYX
 
Bool_t eIsSource
 
Bool_t eIsTarget
 
Bool_t eNeedModified
 
Float_t ePatternBTDensity_modified [114]
 
Float_t ePatternBTDensity_orig [114]
 
TProfile * eProfileBTdens_vs_PID_source
 
Float_t eProfileBTdens_vs_PID_source_meanX
 
Float_t eProfileBTdens_vs_PID_source_meanY
 
Float_t eProfileBTdens_vs_PID_source_rmsX
 
Float_t eProfileBTdens_vs_PID_source_rmsY
 
TProfile * eProfileBTdens_vs_PID_target
 
Float_t eProfileBTdens_vs_PID_target_meanX
 
Float_t eProfileBTdens_vs_PID_target_meanY
 
Float_t eProfileBTdens_vs_PID_target_rmsX
 
Float_t eProfileBTdens_vs_PID_target_rmsY
 
Float_t maxX
 
Float_t maxY
 
Float_t minX
 
Float_t minY
 
Int_t NbinsX
 
Int_t NbinsY
 

Constructor & Destructor Documentation

◆ EdbShowPVRQuality() [1/3]

EdbShowPVRQuality::EdbShowPVRQuality ( )

◆ EdbShowPVRQuality() [2/3]

EdbShowPVRQuality::EdbShowPVRQuality ( EdbPVRec ali)
29{
30 // Default Constructor with a EdbPVRec object.
31 // (the EdbPVRec object corresponds to the collection of plates scanned with the collection
32 // of basetracks (and additionally linked tracks)).
33 // This constructor automatically checks the original data for background level and
34 // creates a new EdbPVRec object that contains only segments that fulfill the quality
35 // cut in accordance with the desired (predefined) background level.
36 // If general basetrack density is lower than 20BT/mm2 then no cleaning is done.
37 cout << "EdbShowPVRQuality::EdbShowPVRQuality(EdbPVRec* ali) Constructor (does automatically all in one...)"<<endl;
38
39 Init();
40 Set0();
41
42 if (NULL == ali || ali->Npatterns()<1) {
43 cout << "WARNING EdbShowPVRQuality(EdbPVRec* ali) ali is no good EdbPVRec object! Check!" << endl;
44 return;
45 }
46
47 // Set ali as eAli_orig
49 eIsSource=kTRUE;
50
54 Print();
55}
Int_t Npatterns() const
Definition: EdbPattern.h:366
void Execute_ConstantBTDensity()
Definition: EdbShowPVRQuality.cxx:541
void CreateEdbPVRec()
Definition: EdbShowPVRQuality.cxx:1024
void Print()
Definition: EdbShowPVRQuality.cxx:408
void CheckEdbPVRec()
Definition: EdbShowPVRQuality.cxx:215
void Set0()
Definition: EdbShowPVRQuality.cxx:102
EdbPVRec * eAli_orig
Definition: EdbShowPVRQuality.h:37
void Init()
Definition: EdbShowPVRQuality.cxx:164
Bool_t eIsSource
Definition: EdbShowPVRQuality.h:40
EdbPVRec * ali
Definition: test_oracle.C:9
#define NULL
Definition: nidaqmx.h:84

◆ EdbShowPVRQuality() [3/3]

EdbShowPVRQuality::EdbShowPVRQuality ( EdbPVRec ali,
Float_t  BTDensityTargetLevel 
)
60{
61 // Default Constructor with a EdbPVRec object and the desired basetrack density target level.
62 // Does same as constructor EdbShowPVRQuality::EdbShowPVRQuality(EdbPVRec* ali) but now with adjustable
63 // background level.
64 // If general basetrack density is lower than BTDensityTargetLevel BT/mm2 then no cleaning is done.
65 cout << "EdbShowPVRQuality::EdbShowPVRQuality(EdbPVRec* ali) Constructor (does automatically all in one...)"<<endl;
66
67 Init();
68 Set0();
69
70 if (NULL == ali || ali->Npatterns()<1) {
71 cout << "WARNING EdbShowPVRQuality(EdbPVRec* ali) ali is no good EdbPVRec object! Check!" << endl;
72 return;
73 }
74
75 // Set ali as eAli_orig
77 eIsSource=kTRUE;
78
79 // Set BTDensityTargetLevel
80 SetBTDensityLevel(BTDensityTargetLevel);
81 cout << "EdbShowPVRQuality::EdbShowPVRQuality(EdbPVRec* ali) GetBTDensityLevel() " << GetBTDensityLevel() << endl;
82
86 Print();
87}
void SetBTDensityLevel(Float_t BTDensityLevel)
Definition: EdbShowPVRQuality.h:96
Float_t GetBTDensityLevel()
Definition: EdbShowPVRQuality.h:166

◆ ~EdbShowPVRQuality()

EdbShowPVRQuality::~EdbShowPVRQuality ( )
virtual
92{
93 // Default Destructor
94 cout << "EdbShowPVRQuality::~EdbShowPVRQuality()"<<endl;
95 delete eHistChi2W;
96 delete eHistYX;
98}
TProfile * eProfileBTdens_vs_PID_source
Definition: EdbShowPVRQuality.h:63
TH2F * eHistYX
Definition: EdbShowPVRQuality.h:57
TH2F * eHistChi2W
Definition: EdbShowPVRQuality.h:56

Member Function Documentation

◆ CheckEdbPVRec()

void EdbShowPVRQuality::CheckEdbPVRec ( )
216{
217 // Main function to check if the EdbPVRec object of the scanned data is of low/high background.
218 // Following steps are carried out:
219 // Get plate, count number of basetracks in the unit area (1x1cm^2).
220 // Fill (draw if desired (like in EDA display)) histogram with the entries of the unit area.
221 // Get mean of the histogram, compare this value with the reference value.
222 // The histogram covers all the area of one emulsion. (for the record: the old ORFEO MC
223 // simulation gives not the same position as data does. The area of the histogramm was largely
224 // increased to cover both cases).
225 // It gives a good estimation of the density. Spikes in some plates, or in some zones are not
226 // checked for, this is on the todo list, but maybe not so important.
227
228 Log(2,"EdbShowPVRQuality::CheckEdbPVRec","EdbShowPVRQuality::CheckEdbPVRec");
229
230 if (!eIsSource) {
231 cout << "EdbShowPVRQuality::CheckEdbPVRec eIsSource= " << eIsSource << ". This means no source set. Return!" << endl;
232 return;
233 }
234 // Check the patterns of the EdbPVRec:
236// cout << "EdbShowPVRQuality::CheckEdbPVRec eAli_orig->Npatterns()= " << eAli_maxNpatterns << endl;
237 if (eAli_maxNpatterns>57) cout << " This tells us not yet if we do have one/two brick reconstruction done. A possibility could also be that the dataset was read with microtracks. Further investigation is needed! (On todo list)." << endl;
238 if (eAli_maxNpatterns>114) {
239 cout << "ERROR! EdbShowPVRQuality::CheckEdbPVRec eAli_orig->Npatterns()= " << eAli_maxNpatterns << " is greater than possible basetrack data two bricks. This class does (not yet) work with this large number of patterns. Set maximum patterns to 114!!!." << endl;
241 }
242
243 int Npat = eAli_maxNpatterns;
244 TH1F* histPatternBTDensity = new TH1F("histPatternBTDensity","histPatternBTDensity",200,0,200);
245
246 // Loop over the patterns of the EdbPVRec:
247 for (Int_t i=0; i<Npat; i++) {
248
249 if (i>56) {
250 cout << "WARNING EdbShowPVRQuality::CheckEdbPVRec() Your EdbPVRec object has more than 57 patterns! " << endl;
251 cout << "WARNING EdbShowPVRQuality::CheckEdbPVRec() Check it! " << endl;
252 }
253 if (gEDBDEBUGLEVEL>2) cout << "CheckEdbPVRec Doing Pattern " << i << endl;
254
255
256 eHistYX->Reset(); // important to clean the histogram
257 eHistChi2W->Reset(); // important to clean the histogram
258
260 Int_t npat=pat->N();
261
262 EdbSegP* seg=0;
263 // Loop over the segments of the pattern
264 for (Int_t j=0; j<npat; j++) {
265 seg=(EdbSegP*)pat->At(j);
266 // Very important:
267 // For the data case, we assume the following:
268 // Data (MCEvt<0) will be taken for BTdensity calculation
269 // Sim (MCEvt>0) will NOT be taken for BTdensity calculation
270 // We take it ONLY and ONLY into account if it is especially wished
271 // by the user!
272 // Therefore (s)he needs to know how many Gauge Coupling Parameters
273 // in the Standard Model exist (at least)...
274 Bool_t result=kTRUE;
275 if (seg->MCEvt()>0) {
277 result = kTRUE;
278 // cout << "result = kTRUE !! " << endl;
279 }
280 else {
281 result = kFALSE;
282 }
283 }
284
285 if (gEDBDEBUGLEVEL>4) cout << "Doing segment " << j << " result for bool query is: " << result << endl;
286
287 // Main decision for segment to be kept or not (seg is of MC or data type).
288 if ( kFALSE == result ) continue;
289
290 // For the check, fill the histograms in any case:
291 eHistYX->Fill(seg->Y(),seg->X());
292 eHistChi2W->Fill(seg->W(),seg->Chi2());
293 }
294
295 if (gEDBDEBUGLEVEL>2) cout << "I have filled the eHistYX Histogram. Entries = " << eHistYX->GetEntries() << endl;
296
297 // Important to reset histogram before it is filled.
298 histPatternBTDensity->Reset();
299
300 // Search for empty bins, because they can spoil the overall calulation
301 // of the mean value.
302 Int_t nbins=eHistYX->GetNbinsX()*eHistYX->GetNbinsY();
303 Int_t nemptybinsXY=0;
304 Int_t bincontentXY=0;
305 for (int k=1; k<nbins-1; k++) {
306 if (eHistYX->GetBinContent(k)==0) {
307 ++nemptybinsXY;
308 continue;
309 }
310 bincontentXY=eHistYX->GetBinContent(k);
311 histPatternBTDensity->Fill(bincontentXY);
312 eProfileBTdens_vs_PID_source->Fill(i,bincontentXY);
313 }
314
315 // failsafe warning in case that there are many bins with zero content.
316 // for now we print a error message:
318
319 // Save the density in the variable.
320 ePatternBTDensity_orig[i]=histPatternBTDensity->GetMean();
321
322 }
323
328
329 // No assignment for the eProfileBTdens_vs_PID_target histogram yet.
330 // This will be done in one of the two Execute_ functions.
331
332
333 // This will be commented when using in batch mode...
334 // For now its there for clarity reasons.
335 TCanvas* c1 = new TCanvas();
336 c1->Divide(2,2);
337 c1->cd(1);
338 eHistYX->DrawCopy("colz");
339 c1->cd(2);
340 eHistChi2W->DrawCopy("colz");
341 c1->cd(3);
342 histPatternBTDensity->DrawCopy("");
343 c1->cd(4);
344 eProfileBTdens_vs_PID_source->Draw("profileZ");
345 eProfileBTdens_vs_PID_source->GetXaxis()->SetRangeUser(0,eAli_maxNpatterns+2);
346 c1->cd();
347
348 eHistYX->Reset();
349 eHistChi2W->Reset();
350
351 Log(2,"EdbShowPVRQuality::CheckEdbPVRec","EdbShowPVRQuality::CheckEdbPVRec...done");
352 return;
353}
bool Log(int level, const char *location, const char *fmt,...)
Definition: EdbLog.cxx:75
Int_t npat
Definition: Xi2HatStartScript.C:33
Definition: EdbPattern.h:273
EdbPattern * GetPattern(int id) const
Definition: EdbPattern.cxx:1721
Definition: EdbSegP.h:21
Float_t X() const
Definition: EdbSegP.h:173
Float_t Chi2() const
Definition: EdbSegP.h:157
Float_t Y() const
Definition: EdbSegP.h:174
Float_t W() const
Definition: EdbSegP.h:151
Int_t MCEvt() const
Definition: EdbSegP.h:145
Int_t N() const
Definition: EdbPattern.h:86
EdbPoint * At(int i) const
Definition: EdbPattern.h:87
Float_t eProfileBTdens_vs_PID_source_meanY
Definition: EdbShowPVRQuality.h:64
Float_t ePatternBTDensity_orig[114]
Definition: EdbShowPVRQuality.h:52
Bool_t eBTDensityLevelCalcMethodMC
Definition: EdbShowPVRQuality.h:49
Int_t eBTDensityLevelCalcMethodMCConfirmationNumber
Definition: EdbShowPVRQuality.h:50
Int_t eAli_maxNpatterns
Definition: EdbShowPVRQuality.h:43
Float_t eProfileBTdens_vs_PID_source_rmsY
Definition: EdbShowPVRQuality.h:65
void CheckFilledXYSize()
Definition: EdbShowPVRQuality.cxx:1126
Float_t eProfileBTdens_vs_PID_source_rmsX
Definition: EdbShowPVRQuality.h:65
Float_t eProfileBTdens_vs_PID_source_meanX
Definition: EdbShowPVRQuality.h:64
gEDBDEBUGLEVEL
Definition: energy.C:7
TCanvas * c1
Definition: energy.C:13
new TCanvas()

◆ CheckFilledXYSize()

void EdbShowPVRQuality::CheckFilledXYSize ( )
1127{
1128 // Check the bins filled in the actual pattern.
1129 // The histogram of the XY distribution is analysed.
1130 // A warning is given if more than 10 percent of
1131 // the bins are empty.
1132 // In this case one should look closer at the specific
1133 // plate distribution.
1134
1135 if (gEDBDEBUGLEVEL>3) cout << "----- void EdbShowPVRQuality::CheckFilledXYSize() return maximum/minimum entries of eHistYX -----" << endl;
1136 Int_t nbx,nby=0;
1137 nbx=eHistYX->GetNbinsX();
1138 nby=eHistYX->GetNbinsY();
1139
1140 Int_t n1x= FindFirstBinAbove(eHistYX,0,1);
1141 Int_t n1y= FindFirstBinAbove(eHistYX,0,2);
1142 Int_t n2x= FindLastBinAbove(eHistYX,0,1);
1143 Int_t n2y= FindLastBinAbove(eHistYX,0,2);
1144 Int_t width_x=0;
1145 Int_t width_y=0;
1146 width_x=TMath::Abs(eHistYX->GetXaxis()->GetBinCenter(n1x)-eHistYX->GetXaxis()->GetBinCenter(n2x));
1147 width_y=TMath::Abs(eHistYX->GetYaxis()->GetBinCenter(n1y)-eHistYX->GetYaxis()->GetBinCenter(n2y));
1148
1149 // Now check the number of empty bins between! the filled area
1150 // within (FindFirstBinAbove,FindLastBinAbove)
1151 // This function is NOT optimized for speed :-)
1152 Int_t NonEmptyBins=0;
1153 Int_t nBins=0;
1154 for (Int_t i=n1x; i<n2x; ++i) {
1155 for (Int_t j=n1y; j<n2y; ++j) {
1156 ++nBins;
1157 if (eHistYX->GetBinContent(i,j)==0) continue;
1158 if (eHistYX->GetBinContent(i,j)==0) continue;
1159 ++NonEmptyBins;
1160 }
1161 }
1162 Float_t FractionOfEmptyBins=1-(Float_t(NonEmptyBins)/Float_t(nBins));
1163 if (FractionOfEmptyBins>0.1) cout << "WARNING: void EdbShowPVRQuality::CheckFilledXYSize() FractionOfEmptyBins = " << FractionOfEmptyBins << endl;
1164
1165 if (gEDBDEBUGLEVEL>3) {
1166 cout << "---- NonEmptyBins bins in covered area: = " << NonEmptyBins << endl;
1167 cout << "---- nBins totale bins in covered area: = " << nBins << endl;
1168 cout << "---- Extrem Center Endpoints of eHistYX --- " << endl;
1169 cout << "---- nbxMin= " << n1x << endl;
1170 cout << "---- nbxMax= " << n1y << endl;
1171 cout << "---- nbyMin= " << n2x << endl;
1172 cout << "---- nbyMax= " << n2y << endl;
1173 cout << "---- nbxMin= " << eHistYX->GetXaxis()->GetBinCenter(n1x) << endl;
1174 cout << "---- nbxMax= " << eHistYX->GetYaxis()->GetBinCenter(n1y) << endl;
1175 cout << "---- nbxMin= " << eHistYX->GetXaxis()->GetBinCenter(n2x) << endl;
1176 cout << "---- nbxMax= " << eHistYX->GetYaxis()->GetBinCenter(n2y) << endl;
1177 cout << "----- void EdbShowPVRQuality::CheckFilledXYSize() ... done. " << endl;
1178 }
1179 return;
1180}
Int_t FindLastBinAbove(TH1 *hist, Double_t threshold, Int_t axis)
Definition: EdbShowPVRQuality.cxx:1599
Int_t FindFirstBinAbove(TH1 *hist, Double_t threshold, Int_t axis)
Definition: EdbShowPVRQuality.cxx:1579

◆ CheckSegmentQualityInPattern_ConstBTDens()

Bool_t EdbShowPVRQuality::CheckSegmentQualityInPattern_ConstBTDens ( EdbPVRec ali,
Int_t  PatternAtNr,
EdbSegP seg 
)
983{
984 // Core function to check if a basetrack of the specific pattern matches the expectations
985 // for the desired quality cut (calculated on the estimations in CheckEdbPVRec(). ).
986 // Implementation for the Cuttype 0: Constant BT Density
987 // ---
988 // Note: since the eCutp1[i] values are calculated with
989 // this pattern->At() scheme labelling,
990 // its not necessarily guaranteed that seg->PID gives correct this
991 // number back. Thats why we have to give the PatternAtNr here again.
992 //
993 // And: it is not checked here if seg is contained in this specific
994 // pattern. It looks only for the quality cut!
995 if (gEDBDEBUGLEVEL>3) cout << "seg->W()* eCutp1[PatternAtNr] - eCutp0[PatternAtNr] = " << seg->W()* eCutp1[PatternAtNr] - eCutp0[PatternAtNr] << endl;
996
997 // Constant BT density cut:
998 if (seg->Chi2() >= seg->W()* eCutp1[PatternAtNr] - eCutp0[PatternAtNr]) return kFALSE;
999
1000 if (gEDBDEBUGLEVEL>3) cout <<"EdbShowPVRQuality::CheckSegmentQualityInPattern_ConstBTDens() Segment " << seg << " has passed ConstBTDens cut!" << endl;
1001 return kTRUE;
1002}
Float_t eCutp0[114]
Definition: EdbShowPVRQuality.h:74
Float_t eCutp1[114]
Definition: EdbShowPVRQuality.h:75

◆ CheckSegmentQualityInPattern_ConstQual()

Bool_t EdbShowPVRQuality::CheckSegmentQualityInPattern_ConstQual ( EdbPVRec ali,
Int_t  PatternAtNr,
EdbSegP seg 
)
1007{
1008 // Core function to check if a basetrack of the specific pattern matches the expectations
1009 // for the desired quality cut (calculated on the estimations in CheckEdbPVRec(). ).
1010 // Implementation for the Cuttype 1: Constant Quality.
1011 // ---
1012 // See comments in CheckSegmentQualityInPattern_ConstBTDens
1013 // Constant BT quality cut:
1015 if (agreementChi2>eAgreementChi2WDistCut[PatternAtNr]) return kFALSE;
1016
1017 if (gEDBDEBUGLEVEL>3) cout <<"EdbShowPVRQuality::CheckSegmentQualityInPattern_ConstQual() Segment " << seg << " has passed ConstQual cut!" << endl;
1018 return kTRUE;
1019}
Float_t eAgreementChi2CutMeanChi2
Definition: EdbShowPVRQuality.h:79
Float_t eAgreementChi2CutRMSW
Definition: EdbShowPVRQuality.h:82
Float_t eAgreementChi2CutRMSChi2
Definition: EdbShowPVRQuality.h:80
Float_t eAgreementChi2CutMeanW
Definition: EdbShowPVRQuality.h:81
Float_t eAgreementChi2WDistCut[114]
Definition: EdbShowPVRQuality.h:78

◆ ClassDef()

EdbShowPVRQuality::ClassDef ( EdbShowPVRQuality  ,
 
)

◆ CreateEdbPVRec()

void EdbShowPVRQuality::CreateEdbPVRec ( )
1025{
1026 // Creates the cleaned EdbPVRec object, containing only segments that
1027 // satisfied the cutcriteria.
1028 // Attention: the couples structure and the tracking structure will be lost,
1029 // this EdbPVRec object is only useful for the list of segments (shower reco).
1030 // DO NOT USE THIS ROUTINE FOR GENERAL I/O and/or GENERAL EdbPVRec operations!
1031
1032 cout << "----- void EdbShowPVRQuality::CreateEdbPVRec() -----" << endl;
1033 if (gEDBDEBUGLEVEL>2) {
1034 cout << "----- " << endl;
1035 cout << "----- This function makes out of the original eAli" << endl;
1036 cout << "----- a new EdbPVRec object having only those seg-" << endl;
1037 cout << "----- ments in it which satisfy the cutcriteria " << endl;
1038 cout << "----- determined in Execute_ConstantBTDensity, Execute_ConstantQuality" << endl;
1039 cout << "----- " << endl;
1040 cout << "----- WARNING: the couples structure and the tracking structure" << endl;
1041 cout << "----- will be lost, this PVR object is only useful for the" << endl;
1042 cout << "----- list of Segments (==ShowReco...) ... " << endl;
1043 cout << "----- DO NOT USE THIS ROUTINE FOR GENERAL I/O and/or EdbPVRec operations!" << endl;
1044 cout << "----- " << endl;
1045 cout << "CreateEdbPVRec() Mode 0:" << eCutMethodIsDone[0] << endl;
1046 cout << "CreateEdbPVRec() Mode 1:" << eCutMethodIsDone[1] << endl;
1047 cout << "----- " << endl;
1048 cout << "----- ----------------------------------------------" << endl;
1049 }
1050
1051 // Checks
1052 if (NULL==eAli_orig || eIsSource==kFALSE) {
1053 cout << "WARNING! EdbShowPVRQuality::CreateEdbPVRec NULL==eAli_orig || eIsSource==kFALSE return."<<endl;
1054 return;
1055 }
1056 // Checks
1057 if (eCutMethodIsDone[0]==kFALSE && eCutMethodIsDone[1]==kFALSE) {
1058 cout << "WARNING! EdbShowPVRQuality::CreateEdbPVRec eCutMethodIsDone[0]==kFALSE && eCutMethodIsDone[1]==kFALSE return."<<endl;
1059 return;
1060 }
1061
1062 // Make a new PVRec object anyway
1063 eAli_modified = new EdbPVRec();
1064
1065 // These two lines dont compile yet ... (???) ...
1066 // EdbScanCond* scancond = eAli_orig->GetScanCond();
1067 // eAli_modified->SetScanCond(*scancond);
1068
1069 Float_t agreementChi2;
1070
1071 // This makes pointer copies of patterns with segments list.
1072 // wARNING: the couples structure and the tracking structure
1073 // will be lost, this PVR object is only useful for the
1074 // list of Segments (==ShowReco...) ...
1075
1076 // Priority has the Execute_ConstantQuality cut.
1077 // If this was done we take this...:
1078
1079 // Loop over patterns
1080 for (int i = 0; i <eAli_orig->Npatterns(); i++ ) {
1081 EdbPattern* pat = eAli_orig->GetPattern(i);
1082 EdbPattern* pt= new EdbPattern();
1083 // SetPattern Values to the parent patterns:
1084 pt->SetID(pat->ID());
1085 pt->SetPID(pat->PID());
1086 pt->SetZ(pat->Z());
1087
1088 // Loop over segments
1089 for (int j = 0; j <pat->N(); j++ ) {
1090 EdbSegP* seg = pat->GetSegment(j);
1091
1092 // Put here the cut condition ...
1093 if (eCutMethodIsDone[0]==kTRUE && eCutMethodIsDone[1]==kFALSE) {
1094 // Constant BT density cut:
1095 if (seg->Chi2() >= seg->W()* eCutp1[i] - eCutp0[i]) continue;
1096 }
1097 else if (eCutMethodIsDone[1]==kTRUE) {
1098 // Constant Quality cut:
1100 if (agreementChi2>eAgreementChi2WDistCut[i]) continue;
1101 }
1102 else {
1103 // do nothing;
1104 }
1105
1106 // Add segment:
1107 pt->AddSegment(*seg);
1108 }
1110 }
1111
1112 eIsTarget=kTRUE;
1113
1114 //---------------------
1115 cout << "----- void EdbShowPVRQuality::CreateEdbPVRec()...Created new EdbPVR object at address " << eAli_modified << endl;
1116 //---------------------
1117 cout << "----- void EdbShowPVRQuality::CreateEdbPVRec()...done." << endl;
1118 return;
1119}
TPaveText * pt
Definition: Canv_SYSTEMATICS_ALLCOMBINED__RMSEnergy__vs__Energy__ELECTRON.C:160
Definition: EdbPVRec.h:148
int PID() const
Definition: EdbPattern.h:320
int ID() const
Definition: EdbPattern.h:319
void AddPattern(EdbPattern *pat)
Definition: EdbPattern.cxx:1707
Float_t Z() const
Definition: EdbPattern.h:84
EdbSegP * GetSegment(int i) const
Definition: EdbPattern.h:66
Bool_t eIsTarget
Definition: EdbShowPVRQuality.h:41
Bool_t eCutMethodIsDone[2]
Definition: EdbShowPVRQuality.h:48
EdbPVRec * eAli_modified
Definition: EdbShowPVRQuality.h:38

◆ Execute_ConstantBTDensity()

void EdbShowPVRQuality::Execute_ConstantBTDensity ( )
542{
543 // Execute the modified cut routines to achieve the basetrack density level,
544 // after application the specific cut on the segments of the specific plate (pattern).
545 // The Constant BT Density is defined by the number of BT/mm2 in the histogram.
546
547 for (int i=0; i<80; ++i) cout << "-";
548 cout << endl;
549 cout << "-";
550 cout << endl;
551 cout << "EdbShowPVRQuality::Execute_ConstantBTDensity " << endl;
552 cout << "-";
553 cout << endl;
554 for (int i=0; i<80; ++i) cout << "-";
555 cout << endl;
556
557 if (!eIsSource) return;
558 if (NULL==eAli_orig) return;
559
560 // Check the patterns of the EdbPVRec:
562 if (eAli_maxNpatterns>57) cout << " This tells us not yet if we do have one/two brick reconstruction done. A possibility could also be that the dataset was read with microtracks. Further investigation is needed! (On todo list)." << endl;
563 if (eAli_maxNpatterns>114) {
564 cout << "ERROR! EdbShowPVRQuality::Execute_ConstantBTDensity eAli_orig->Npatterns()= " << eAli_maxNpatterns << " is greater than possible basetrack data two bricks. This class does (not yet) work with this large number of patterns. Set maximum patterns to 114!!!." << endl;
566 }
567
568 TH1F* histPatternBTDensity = new TH1F("histPatternBTDensity","histPatternBTDensity",200,0,200);
569
570 // Loop over the patterns:
571 for (int i=0; i<eAli_maxNpatterns; i++) {
572 if (i>56) {
573 cout << "WARNING EdbShowPVRQuality::Execute_ConstantBTDensity() Your EdbPVRec object has more than 57 patterns! " << endl;
574 cout << "WARNING EdbShowPVRQuality::Execute_ConstantBTDensity() Check it! " << endl;
575 }
576 if (gEDBDEBUGLEVEL>2) cout << "Execute_ConstantBTDensity Doing Pattern " << i << endl;
577
578 // Now the condition loop:
579 // Loop over 20 steps a 0.15,0.145,0.14 ... down to 0.07
580 for (int l=0; l<20; l++) {
581 eHistYX->Reset(); // important to clean the histogram
582 eHistChi2W->Reset(); // important to clean the histogram
583 histPatternBTDensity->Reset(); // important to clean the histogram
584
586 Int_t npat=pat->N();
587 EdbSegP* seg=0;
588 // Loop over the segments.
589 for (int j=0; j<npat; j++) {
590 seg=(EdbSegP*)pat->At(j);
591 // Very important:
592 // For the data case, we assume the following:
593 // Data (MCEvt<0) will be taken for BTdensity calculation
594 // Sim (MCEvt>0) will NOT be taken for BTdensity calculation
595 // We take it ONLY and ONLY into account if it is especially wished
596 // by the user!
597 // Therefore (s)he needs to know how many Gauge Coupling Parameters
598 // in the Standard Model exist (at least)...
599 Bool_t result=kTRUE;
600 if (seg->MCEvt()>0) {
602 result = kTRUE;
603 // cout << "result = kTRUE !! " << endl;
604 }
605 else {
606 result = kFALSE;
607 }
608 }
609
610 if (gEDBDEBUGLEVEL>4) cout << "Doing segment " << j << " result for bool query is: " << result << endl;
611
612 // Main decision for segment to be kept or not (seg is of MC or data type).
613 if ( kFALSE == result ) continue;
614 // Constant BT density cut:
615 if (seg->Chi2() >= seg->W()* eCutp1[i] - eCutp0[i]) continue;
616
617 eHistYX->Fill(seg->Y(),seg->X());
618 eHistChi2W->Fill(seg->W(),seg->Chi2());
619 }
620
621 if (gEDBDEBUGLEVEL>2) cout << "I have filled the eHistYX Histogram. Entries = " << eHistYX->GetEntries() << endl;
622
623 int nbins=eHistYX->GetNbinsX()*eHistYX->GetNbinsY();
624 for (int k=1; k<nbins-1; k++) {
625 if (eHistYX->GetBinContent(k)==0) continue;
626 histPatternBTDensity->Fill(eHistYX->GetBinContent(k));
627 }
628
629 ePatternBTDensity_modified[i]=histPatternBTDensity->GetMean();
630 if (gEDBDEBUGLEVEL>2) cout <<"Execute_ConstantBTDensity Loop l= " << l << ": for the eCutp1[i] : " << eCutp1[i] << " we have a dens: " << ePatternBTDensity_modified[i] << endl;
631
632 // Now the condition check:
634 if (gEDBDEBUGLEVEL>2) cout << "Execute_ConstantBTDensity We reached the loop end due to good BT density level ... and break loop." << endl;
635 break;
636 }
637 else {
638 // Next step, tighten cut:
639 eCutp1[i] += -0.005;
640 }
641
642 } // of condition loop...
643
644 // Fill target profile histogram:
645 Float_t bincontentXY=histPatternBTDensity->GetMean();
646 eProfileBTdens_vs_PID_target->Fill(i,bincontentXY);
647
648 } // of Npattern loops..
649
650 eCutMethodIsDone[0]=kTRUE;
651
652 // Check if modified or original PVRec object should be returned:
654 eNeedModified=kTRUE;
655 cout << "----------void EdbShowPVRQuality::Execute_ConstantBTDensity() Check if modified or original PVRec object should be returned: " << endl;
656 cout << "----------void EdbShowPVRQuality::Execute_ConstantBTDensity() " << eNeedModified << endl;
657 }
658
659
660 // This will be commented when using in batch mode...
661 // For now its there for clarity reasons.
662 TCanvas* c1 = new TCanvas();
663 c1->Divide(2,2);
664 c1->cd(1);
665 eHistYX->DrawCopy("colz");
666 c1->cd(2);
667 eHistChi2W->DrawCopy("colz");
668 c1->cd(3);
669 histPatternBTDensity->DrawCopy("");
670 c1->cd(4);
671 eProfileBTdens_vs_PID_source->Draw("profileZ");
672 eProfileBTdens_vs_PID_source->GetXaxis()->SetRangeUser(0,eAli_maxNpatterns+2);
673 eProfileBTdens_vs_PID_target->SetLineStyle(2);
674 eProfileBTdens_vs_PID_target->Draw("profileZsame");
675 c1->cd();
676 histPatternBTDensity->Reset();
677 eHistYX->Reset();
678 eHistChi2W->Reset();
679
680
685
690
691 delete histPatternBTDensity;
692
693 cout << "----------void EdbShowPVRQuality::Execute_ConstantBTDensity() Cuts are done and saved to obtain desired BT density. " << endl;
694 cout << "----------void EdbShowPVRQuality::Execute_ConstantBTDensity() If you want to apply the cuts now, run the CreateEdbPVRec() function now. " << endl;
695 cout << "----------void EdbShowPVRQuality::Execute_ConstantBTDensity()....done." << endl;
696 return;
697}
Float_t ePatternBTDensity_modified[114]
Definition: EdbShowPVRQuality.h:53
Float_t eProfileBTdens_vs_PID_target_meanY
Definition: EdbShowPVRQuality.h:68
Bool_t eNeedModified
Definition: EdbShowPVRQuality.h:39
Float_t eProfileBTdens_vs_PID_target_meanX
Definition: EdbShowPVRQuality.h:68
Float_t eProfileBTdens_vs_PID_target_rmsX
Definition: EdbShowPVRQuality.h:69
TProfile * eProfileBTdens_vs_PID_target
Definition: EdbShowPVRQuality.h:67
Float_t eBTDensityLevel
Definition: EdbShowPVRQuality.h:47
Float_t eProfileBTdens_vs_PID_target_rmsY
Definition: EdbShowPVRQuality.h:69

◆ Execute_ConstantQuality()

void EdbShowPVRQuality::Execute_ConstantQuality ( )

______ now same code as in the function Execute_ConstantBTDensity ___________________

703{
704 // The cut method that compares the passing muon tracks with all scanned segments.
705 // This may help to improve, since it takes into account the actual segment quality,
706 // which varies from scan to scan anyway.
707 // Works for tracks passing the volume to extract their mean chi2/W. Then a distance
708 // measurement Sqrt{0.5*[ ((chi2-chi2mean)/chi2rms)^2 + ((W-Wmean)/Wrms)^2 )] } is
709 // build and the value of single basetracks is compared to this variable. Starting from
710 // dist_max = 3 going down until desired target level is achieved.
711 // If eAli->Tracks is there we take them from there.
712 // If not, we try if there is a file linked_tracks.root, and we take tracks from there.
713 // If that doesnt work either, nothing is done.
714
715 if (!eIsSource) return;
716
717 for (int i=0; i<80; ++i) cout << "-";
718 cout << endl;
719 cout << "-";
720 cout << endl;
721 cout << "EdbShowPVRQuality::Execute_ConstantBTDensity " << endl;
722 cout << "-";
723 cout << endl;
724 for (int i=0; i<80; ++i) cout << "-";
725 cout << endl;
726 cout << "---------- Works for tracks passing the volume to extract their mean chi2/W " << endl;
727 cout << "---------- If eAli->Tracks is there we take them from there." << endl;
728 cout << "---------- If not, we try if there is a file linked_tracks.root " << endl;
729 cout << "---------- we take tracks from there. " << endl;
730 cout << "---------- If that doesnt work either, nothing is done." << endl;
731 for (int i=0; i<80; ++i) cout << "-";
732 cout << endl;
733 cout << "-";
734 cout << endl;
735
736 cout << "EdbShowPVRQuality::Execute_ConstantQuality() eAli_orig.eTracks :" << eAli_orig->eTracks << endl;
737
738 Float_t meanChi2=0.5;
739 Float_t rmsChi2=0.2;
740 Float_t meanW=22;
741 Float_t rmsW=4;
742 Float_t agreementChi2=100;
743
744 // No eAli.Tracks ? Look for tracks in linked_track.root:
745 if (NULL == eAli_orig->eTracks) {
746 cout << "EdbShowPVRQuality::Execute_ConstantQuality() No eAli.Tracks. Look for tracks in linked_track.root" << endl;
747 TFile* trackfile = new TFile("linked_tracks.root");
748 trackfile->ls();
749 TTree* tracks = (TTree*)trackfile->Get("tracks");
750 if (NULL == tracks) {
751 cout << "EdbShowPVRQuality::Execute_ConstantQuality() No tracks in linked_track.root file. Return, leave eAli_orig unchanged and dont do any cleaning. You might try Execute_ConstantBTDensity instead. " << endl;
752 return;
753 }
754 // TH1F* h1; Attention the usage of _npl_ is disvavoured
755 // because this can cause problems.
756 tracks->Draw("nseg>>h(60,0,60)","","");
757 TH1F *h1 = (TH1F*)gPad->GetPrimitive("h");
758
759 // Short implementation of getting the last bin filled:
760 int lastfilledbin=0;
761 for (int k=1; k<h1->GetNbinsX()-1; k++) {
762 if (h1->GetBinContent(k)>0) lastfilledbin=k;
763 }
764
765 Float_t travellenghtpercentage=Float_t(lastfilledbin)/Float_t(eAli_maxNpatterns);
766
767 TString cutstring = TString(Form("nseg>=%d",int(h1->GetBinCenter(lastfilledbin-3)) ));
768 tracks->Draw("s.eChi2>>hChi2(100,0,2)",cutstring);
769 TH1F *hChi2 = (TH1F*)gPad->GetPrimitive("hChi2");
770
771 cout << "EdbShowPVRQuality::Execute_ConstantQuality() Found tracks in the file! Good. Info:" << endl;
772 cout << "EdbShowPVRQuality::Execute_ConstantQuality() The long tracks in the file have an average percentage length of:" << travellenghtpercentage << endl;
773 cout << "EdbShowPVRQuality::Execute_ConstantQuality() Mean(RMS) of Chi2 distribution of passing tracks: " << hChi2->GetMean() << "+-" << hChi2->GetRMS() << endl;
774
775 TCanvas* c1 = new TCanvas();
776 c1->cd();
777 tracks->Draw("s.eW>>hW(50,0,50)",cutstring);
778 TH1F *hW = (TH1F*)gPad->GetPrimitive("hW");
779 cout << "EdbShowPVRQuality::Execute_ConstantQuality() Mean(RMS) of W distribution of passing tracks: " << hW->GetMean() << "+-" << hW->GetRMS() << endl;
780
781 meanChi2=hChi2->GetMean();
782 rmsChi2=hChi2->GetRMS();
783 meanW=hW->GetMean();
784 rmsW=hW->GetRMS();
785
786 // Quick check if these values are reasonable:
787 if (TMath::Abs(meanChi2-0.5)>1 || TMath::Abs(meanW-22)>6) {
788 cout << "WARNING EdbShowPVRQuality::Execute_ConstantQuality() The tracks might have a strange distribution. You are urgently requested to check these manually!!!"<<endl;
789 }
790 else {
791 cout << "EdbShowPVRQuality::Execute_ConstantQuality() The tracks have reasonable values."<<endl;
792 }
793
794
795 // since the values for the passing tracks are
796 // calculated for the whole volume, we assume that the cutvalues
797 // are valid for all plates.
802 }
803
805
806 // Check the patterns of the EdbPVRec:
808 cout << "EdbShowPVRQuality::Execute_ConstantQuality eAli_orig->Npatterns()= " << eAli_maxNpatterns << endl;
809 if (eAli_maxNpatterns>57) cout << " This tells us not yet if we do have one/two brick reconstruction done. A possibility could also be that the dataset was read with microtracks. Further investigation is needed! (On todo list)." << endl;
810 if (eAli_maxNpatterns>114) {
811 cout << "ERROR! EdbShowPVRQuality::Execute_ConstantQuality eAli_orig->Npatterns()= " << eAli_maxNpatterns << " is greater than possible basetrack data two bricks. This class does (not yet) work with this large number of patterns. Set maximum patterns to 114!!!." << endl;
813 }
814
815 TH1F* histPatternBTDensity = new TH1F("histPatternBTDensity","histPatternBTDensity",100,0,200);
816 TH1F* histagreementChi2 = new TH1F("histagreementChi2","histagreementChi2",100,0,5);
817
818
819 cout << "Execute_ConstantQuality Loop over the patterns..." << endl;
820 for (int i=0; i<eAli_maxNpatterns; i++) {
821 if (i>56) {
822 cout << "ERROR EdbShowPVRQuality::Execute_ConstantQuality() Your EdbPVRec object has more than 57 patterns! " << endl;
823 cout << "ERROR EdbShowPVRQuality::Execute_ConstantQuality() Check it! " << endl;
824 }
825
826 if (gEDBDEBUGLEVEL>2) cout << "Execute_ConstantQuality Doing Pattern " << i << endl;
827
828 // Now the condition loop:
829 // Loop over 30 steps agreementChi2 step 0.05
830 for (int l=0; l<30; l++) {
831
832 if (gEDBDEBUGLEVEL>2) cout << "Execute_ConstantQuality Doing condition loop = " << l << endl;
833
834 eHistYX->Reset(); // important to clean the histogram
835 eHistChi2W->Reset(); // important to clean the histogram
836 histPatternBTDensity->Reset(); // important to clean the histogram
837
839 Int_t npat=pat->N();
840 EdbSegP* seg=0;
841
842 if (gEDBDEBUGLEVEL>2) cout << "Execute_ConstantQuality Loop over segments of pattern " << i << ",Number of segments= " << npat << endl;
843 for (int j=0; j<npat; j++) {
844 seg=(EdbSegP*)pat->At(j);
845
846 if (gEDBDEBUGLEVEL>4) cout << "Execute_ConstantQuality Doing segment= " << j << endl;
847 // Very important:
848 // For the data case, we assume the following:
849 // Data (MCEvt<0) will be taken for BTdensity calculation
850 // Sim (MCEvt>0) will NOT be taken for BTdensity calculation
851 // We take it ONLY and ONLY into account if it is especially wished
852 // by the user!
853 // Therefore (s)he needs to know how many Gauge Coupling Parameters
854 // in the Standard Model exist (at least)...
855 Bool_t result=kTRUE;
856 if (seg->MCEvt()>0) {
858 result = kTRUE;
859 // cout << "result = kTRUE !! " << endl;
860 }
861 else {
862 result = kFALSE;
863 }
864 }
865
866 if (gEDBDEBUGLEVEL>4) cout << "Doing segment " << j << " result for bool query is: " << result << endl;
867
868 // Main decision for segment to be kept or not (seg is of MC or data type).
869 if ( kFALSE == result ) continue;
870
871 // Change here to the quality with values obtained from the tracks.
872 // Constant BT quality cut:
873 agreementChi2=TMath::Sqrt(0.5 *( (seg->Chi2()-meanChi2)/rmsChi2)*((seg->Chi2()-meanChi2)/rmsChi2) + ((seg->W()-meanW)/rmsW)*((seg->W()-meanW)/rmsW) );
874 histagreementChi2->Fill(agreementChi2);
875
876 // Constant BT quality cut:
877 if (agreementChi2>eAgreementChi2WDistCut[i]) continue;
878
879 eHistYX->Fill(seg->Y(),seg->X());
880 eHistChi2W->Fill(seg->W(),seg->Chi2());
881 }
882
883 if (gEDBDEBUGLEVEL>2) cout << "I have filled the eHistYX Histogram. Entries = " << eHistYX->GetEntries() << endl;
884
885 Int_t nbins=eHistYX->GetNbinsX()*eHistYX->GetNbinsY();
886 Int_t nemptybinsXY=0;
887 for (int k=1; k<nbins-1; k++) {
888 if (eHistYX->GetBinContent(k)==0) {
889 ++nemptybinsXY;
890 continue;
891 }
892 histPatternBTDensity->Fill(eHistYX->GetBinContent(k));
893 }
894
895 // failsafe warning in case that there are many bins with zero content.
896 // for now we print a error message:
898
899 ePatternBTDensity_modified[i]=histPatternBTDensity->GetMean();
900
901 if (gEDBDEBUGLEVEL>2) cout <<"Execute_ConstantQuality Loop l= " << l << ": for the eAgreementChi2WDistCut : " << eAgreementChi2WDistCut[i] << " we have a dens: " << ePatternBTDensity_modified[i] << endl;
902
903 // Now the condition check:
905 if (gEDBDEBUGLEVEL>2) cout << "Execute_ConstantQuality We reached the loop end due to good BT density level ... and break loop." << endl;
906 // But dont forget to set values:
907 eCutDistChi2[i]=meanChi2;
908 eCutDistW[i]=meanW;
913 break;
914 }
915 else {
916 // Next step, tighten cut:
917 eAgreementChi2WDistCut[i]+= -0.05;
918 }
919
920 } // of condition loop...
921
922 // Fill target profile histogram:
923 Float_t bincontentXY=histPatternBTDensity->GetMean();
924 cout << "Execute_ConstantQuality histPatternBTDensity->GetMean()." << histPatternBTDensity->GetMean() << endl;
925 eProfileBTdens_vs_PID_target->Fill(i,bincontentXY);
926
927 } // of Npattern loops..
928 cout << "Execute_ConstantQuality Loop over the patterns...done." << endl;
929
930 // Check if modified or original PVRec object should be returned:
932 eNeedModified=kTRUE;
933 cout << "----------void EdbShowPVRQuality::Execute_ConstantBTDensity() Check if modified or original PVRec object should be returned: " << endl;
934 cout << "----------void EdbShowPVRQuality::Execute_ConstantBTDensity() " << eNeedModified << endl;
935 }
936
937 eCutMethodIsDone[1]=kTRUE;
938
939 // This will be commented when using in batch mode...
940 // For now its there for clarity reasons.
941 TCanvas* c1 = new TCanvas();
942 c1->Divide(2,2);
943 c1->cd(1);
944 eHistYX->DrawCopy("colz");
945 c1->cd(2);
946 eHistChi2W->DrawCopy("colz");
947 c1->cd(3);
948 histPatternBTDensity->DrawCopy("");
949 c1->cd(4);
950 eProfileBTdens_vs_PID_source->Draw("profileZ");
951 eProfileBTdens_vs_PID_source->GetXaxis()->SetRangeUser(0,eAli_maxNpatterns+2);
952 eProfileBTdens_vs_PID_target->SetLineStyle(2);
953 eProfileBTdens_vs_PID_target->Draw("profileZsame");
954 c1->cd();
955
960
965
966 histPatternBTDensity->Reset();
967 eHistYX->Reset();
968 eHistChi2W->Reset();
969 delete histPatternBTDensity;
970
971
972 cout << "----------void EdbShowPVRQuality::Execute_ConstantQuality() Cuts are done and saved to obtain desired BT density. " << endl;
973 cout << "----------void EdbShowPVRQuality::Execute_ConstantQuality() If you want to apply the cuts now, run the CreateEdbPVRec() function now. " << endl;
974
975 cout << "----------void EdbShowPVRQuality::Execute_ConstantQuality()....done." << endl;
976 return;
977}
TObjArray * eTracks
Definition: EdbPVRec.h:161
Float_t eCutDistChi2[114]
Definition: EdbShowPVRQuality.h:76
Float_t eCutDistW[114]
Definition: EdbShowPVRQuality.h:77
TTree * tracks
Definition: check_tr.C:19
TH1F * h1
Definition: energy.C:16

◆ FindFirstBinAbove()

Int_t EdbShowPVRQuality::FindFirstBinAbove ( TH1 *  hist,
Double_t  threshold,
Int_t  axis 
)
1579 {
1580 // The TH1 function FindFirstBinAbove is only implemented in
1581 // root version >= 5.24. But since many scanning labs use old root
1582 // versions persistently, I had to copy the TH1 function as a
1583 // memberfunction of EdbShowPVRQuality
1584 // Code taken from
1585 // http://root.cern.ch/root/html/src/TH1.cxx.html#biA7FC
1586 TAxis* ax=0;
1587 if (axis==1) ax = hist->GetXaxis();
1588 if (axis==2) ax = hist->GetYaxis();
1589 if (axis==3) ax = hist->GetZaxis();
1590 int nb = ax->GetNbins();
1591 for (Int_t i=0; i<=nb; i++) {
1592 if (hist->GetBinContent(i)>threshold) return i;
1593 }
1594 return 0;
1595}
void hist()
Definition: init.C:23

◆ FindLastBinAbove()

Int_t EdbShowPVRQuality::FindLastBinAbove ( TH1 *  hist,
Double_t  threshold,
Int_t  axis 
)
1599 {
1600 // The TH1 function FindFirstBinAbove is only implemented in
1601 // root version >= 5.24. But since many scanning labs use old root
1602 // versions persistently, I had to copy the TH1 function as a
1603 // memberfunction of EdbShowPVRQuality
1604 // Code taken from
1605 // http://root.cern.ch/root/html/src/TH1.cxx.html#biA7FC
1606 TAxis* ax=0;
1607 if (axis==1) ax = hist->GetXaxis();
1608 if (axis==2) ax = hist->GetYaxis();
1609 if (axis==3) ax = hist->GetZaxis();
1610 int nb = ax->GetNbins();
1611 for (Int_t i=nb; i>=1; i--) {
1612 if (hist->GetBinContent(i)>threshold) return i;
1613 }
1614 return 0;
1615}

◆ GetagreementChi2Cut() [1/2]

Float_t * EdbShowPVRQuality::GetagreementChi2Cut ( )
inline
196 {
198 }

◆ GetagreementChi2Cut() [2/2]

Float_t EdbShowPVRQuality::GetagreementChi2Cut ( Int_t  patNR)
inline
199 {
200 return eAgreementChi2WDistCut[patNR];
201 }

◆ GetagreementChi2CutMeanChi2()

Float_t EdbShowPVRQuality::GetagreementChi2CutMeanChi2 ( )
inline
183 {
185 }

◆ GetagreementChi2CutMeanW()

Float_t EdbShowPVRQuality::GetagreementChi2CutMeanW ( )
inline
189 {
191 }

◆ GetagreementChi2CutRMSChi2()

Float_t EdbShowPVRQuality::GetagreementChi2CutRMSChi2 ( )
inline
186 {
188 }

◆ GetagreementChi2CutRMSW()

Float_t EdbShowPVRQuality::GetagreementChi2CutRMSW ( )
inline
192 {
194 }

◆ GetBTDensityLevel()

Float_t EdbShowPVRQuality::GetBTDensityLevel ( )
inline
166 {
167 return eBTDensityLevel;
168 }

◆ GetBTDensityLevelCalcMethodMC()

Bool_t EdbShowPVRQuality::GetBTDensityLevelCalcMethodMC ( )
inline
106 {
108 }

◆ GetBTDensityLevelCalcMethodMCConfirmation()

Int_t EdbShowPVRQuality::GetBTDensityLevelCalcMethodMCConfirmation ( )
inline

◆ GetCutMethod()

Int_t EdbShowPVRQuality::GetCutMethod ( )
inline
159 {
160 return eCutMethod;
161 }
Int_t eCutMethod
Definition: EdbShowPVRQuality.h:46

◆ GetCutMethodIsDone()

Bool_t EdbShowPVRQuality::GetCutMethodIsDone ( Int_t  type)
inline
162 {
163 if (type>1) return 0;
164 return eCutMethodIsDone[type];
165 }
Int_t type
Definition: testBGReduction_By_ANN.C:15

◆ GetCutp0() [1/2]

Float_t * EdbShowPVRQuality::GetCutp0 ( )
inline
170 {
171 return eCutp0;
172 }

◆ GetCutp0() [2/2]

Float_t EdbShowPVRQuality::GetCutp0 ( Int_t  patNR)
inline
176 {
177 return eCutp0[patNR];
178 }

◆ GetCutp1() [1/2]

Float_t * EdbShowPVRQuality::GetCutp1 ( )
inline
173 {
174 return eCutp1;
175 }

◆ GetCutp1() [2/2]

Float_t EdbShowPVRQuality::GetCutp1 ( Int_t  patNR)
inline
179 {
180 return eCutp1[patNR];
181 }

◆ GetEdbPVRec() [1/3]

EdbPVRec * EdbShowPVRQuality::GetEdbPVRec ( )
inline
115 {
117 }
EdbPVRec * GetEdbPVRec()
Definition: EdbShowPVRQuality.h:115

◆ GetEdbPVRec() [2/3]

EdbPVRec * EdbShowPVRQuality::GetEdbPVRec ( Bool_t  NeedModified)
inline
127 {
128 cout << "Inline EdbPVRecType= " << NeedModified << endl;
129 if (NeedModified==1) {
130 return eAli_modified;
131 }
132 else {
133 return eAli_orig;
134 }
135 }

◆ GetEdbPVRec() [3/3]

EdbPVRec * EdbShowPVRQuality::GetEdbPVRec ( Int_t  EdbPVRecType)
inline
118 {
119 cout << "Inline EdbPVRecType= " << EdbPVRecType << endl;
120 if (EdbPVRecType==1) {
121 return eAli_modified;
122 }
123 else {
124 return eAli_orig;
125 }
126 }

◆ GetEdbPVRec_modified()

EdbPVRec * EdbShowPVRQuality::GetEdbPVRec_modified ( )
inline
140 {
141 return GetEdbPVRec(1);
142 }

◆ GetEdbPVRec_orig()

EdbPVRec * EdbShowPVRQuality::GetEdbPVRec_orig ( )
inline
137 {
138 return GetEdbPVRec(0);
139 }

◆ GetHistChi2W()

TH2F * EdbShowPVRQuality::GetHistChi2W ( )
inline
150 {
151 return eHistChi2W;
152 }

◆ GetHistYX()

TH2F * EdbShowPVRQuality::GetHistYX ( )
inline
153 {
154 return eHistYX;
155 }

◆ GetTracksFromLinkedTracksRootFile()

TObjArray * EdbShowPVRQuality::GetTracksFromLinkedTracksRootFile ( )
1620{
1621 cout << "EdbShowPVRQuality::GetTracksFromLinkedTracksRootFile()" << endl;
1622// cout << "___TODO___ STEP 1 : Check if linked_tracks.root exists " << endl;
1623// cout << "___TODO___ STEP 2 : Open linked_tracks " << endl;
1624// cout << "___TODO___ STEP 3 : Check if there are tracks inside " << endl;
1625// cout << "___TODO___ STEP 4 : Get Tracks as this TObjArray " << endl;
1626
1627 // Array containing EdbTrackP objects:
1628 TObjArray* trackarray = new TObjArray();
1629
1630 TFile* file = new TFile("linked_tracks.root");
1631 TTree* tr= (TTree*)file->Get("tracks");
1632 TClonesArray *seg= new TClonesArray("EdbSegP",60);
1633 EdbSegP* segment=0;
1634 EdbSegP* s2=0;
1635 EdbSegP* t=0;
1636
1637 int nentr = int(tr->GetEntries());
1638 // check if tracks has entries: if so then ok, otherwise return directly:
1639 if (nentr>0) {
1640 cout << "EdbShowPVRQuality::GetTracksFromLinkedTracksRootFile() " << nentr << " entries Total"<<endl;
1641 }
1642 else {
1643 cout << "EdbShowPVRQuality::GetTracksFromLinkedTracksRootFile() linked_tracks.root file has NO entries...Return NULL." << endl;
1644 return NULL;
1645 }
1646
1647 int nseg,n0,npl;
1648 tr->SetBranchAddress("t." , &t );
1649 tr->SetBranchAddress("s" , &seg );
1650 tr->SetBranchAddress("nseg" , &nseg );
1651 tr->SetBranchAddress("n0" , &n0 );
1652 tr->SetBranchAddress("npl" , &npl );
1653
1654 for (int i=0; i<nentr; i++ ) {
1655 tr->GetEntry(i);
1656 EdbTrackP* realTrack = new EdbTrackP(t);
1657 for (int j=0; j<nseg; j++ ) {
1658 s2=(EdbSegP*) seg->At(j);
1659 segment=(EdbSegP*)s2->Clone();
1660 realTrack->AddSegment(segment);
1661 }
1662 //realTrack->PrintNice();
1663 realTrack->SetNpl(); // Necessary otherwise these variables are not filled.
1664 realTrack->SetN0(); // Necessary otherwise these variables are not filled.
1665 trackarray->Add(realTrack);
1666 }
1667 delete tr;
1668 file->Close();
1669 delete file;
1670 cout << "EdbShowPVRQuality::GetTracksFromLinkedTracksRootFile()...done." << endl;
1671 return trackarray;
1672}
TTree * tr
Definition: Shower_E_FromShowerRoot.C:5
Definition: EdbPattern.h:113
void SetN0(int n0)
Definition: EdbPattern.h:161
void AddSegment(EdbSegP *s)
Definition: EdbPattern.h:214
void SetNpl(int npl)
Definition: EdbPattern.h:168
TTree * t
Definition: check_shower.C:4
EdbSegP * s2
Definition: tlg2couples.C:30
TFile * file
Definition: write_pvr.C:3

◆ Help()

void EdbShowPVRQuality::Help ( )
1186{
1187 // Basic Help Function.
1188
1189 cout << "----------------------------------------------" << endl;
1190 cout << "----- void EdbShowPVRQuality::Help() -----" << endl;
1191 cout << "-----" << endl;
1192 cout << "----- This Class helps you to determine the Quality Cut Plate by Plate" << endl;
1193 cout << "----- for suited BG level for shower reco." << endl;
1194 cout << "----- You use it like this:" << endl;
1195 cout << "----- EdbShowPVRQuality* QualityClass = new EdbShowPVRQuality(gAli)" << endl;
1196 cout << "----- where gAli is an EdbPVRec object pointer (usually the one" << endl;
1197 cout << "----- from the cp files). Then you can get the CutValues with" << endl;
1198 cout << "----- GetCutp0() (for each//all plate) back." << endl;
1199 cout << "-----" << endl;
1200 cout << "----- On TODO list: to be interfaced with the libShower reconstruction mode." << endl;
1201 cout << "-----" << endl;
1202 cout << "-----" << endl;
1203 cout << "-----" << endl;
1204 cout << "----- " << endl;
1205
1206 cout << "----- void EdbShowPVRQuality::Help() -----" << endl;
1207 cout << "----------------------------------------------" << endl;
1208}

◆ Init()

void EdbShowPVRQuality::Init ( void  )
protected

TEMPORARY, later the Standard Geometry should default be OPERA. For now we use a very large histogram that can contain both case x-y ranges. This takes slightly more memory but it shouldnt matter.

165{
166 // Basic Init function that creates objects in the memory which have
167 // to be created only ONE time per class instance.
168
169 Log(2,"EdbShowPVRQuality::Init","EdbShowPVRQuality::Init");
170
171 eProfileBTdens_vs_PID_source = new TProfile("eProfileBTdens_vs_PID_source","eProfileBTdens_vs_PID_source",114,-0.5,113.5,0,200);
172 eProfileBTdens_vs_PID_target = new TProfile("eProfileBTdens_vs_PID_target","eProfileBTdens_vs_PID_target",114,-0.5,113.5,0,200);
173
174 eHistChi2W = new TH2F("eHistChi2W","eHistChi2W",40,0,40,90,0,3);
175 eHistYX = new TH2F("eHistYX","eHistYX",100,0,1,100,0,1);
176
180// if (eHistGeometry==0) SetHistGeometry_OPERA();
181// cout << "EdbShowPVRQuality::Init() /// TEMPORARY SetHistGeometry_OPERAandMC " << endl;
183
184
185 Log(2,"EdbShowPVRQuality::Init","EdbShowPVRQuality::Init...done.");
186 return;
187}
void SetHistGeometry_OPERAandMC()
Definition: EdbShowPVRQuality.cxx:380

◆ Print()

void EdbShowPVRQuality::Print ( )
409{
410 Log(2,"EdbShowPVRQuality::Print","EdbShowPVRQuality::Print");
411 for (int i=0; i<80; ++i) cout << "-";
412 cout << endl;
413 cout << "-";
414 cout << endl;
415 cout << "EdbShowPVRQuality::Print()" << endl;
416 cout << "-";
417 cout << endl;
418 for (int i=0; i<80; ++i) cout << "-";
419 cout << endl;
420 cout << "-" << endl;
421 cout << "EdbShowPVRQuality::Print() --- GENERAL SETTINGS: ---" << endl;
422 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eCutMethod = " << eCutMethod << endl;
423 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eCutMethodIsDone[0] = " << eCutMethodIsDone[0] << endl;
424 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eCutMethodIsDone[1] = " << eCutMethodIsDone[1] << endl;
425 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eBTDensityLevel = " << eBTDensityLevel << endl;
426 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eBTDensityLevelCalcMethodMC = " << eBTDensityLevelCalcMethodMC << endl;
427 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eBTDensityLevelCalcMethodMCConfirmationNumber = " << eBTDensityLevelCalcMethodMCConfirmationNumber << endl;
428 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eHistGeometry = " << eHistGeometry << endl;
429 cout << "-" << endl;
430 cout << "EdbShowPVRQuality::Print() --- SETTINGS: Input data:" << endl;
431 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eAli_orig = " << eAli_orig << endl;
432 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eAli_maxNpatterns = " << eAli_maxNpatterns << endl;
433 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eIsSource = " << eIsSource << endl;
434 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eProfileBTdens_vs_PID_source_meanX = " << eProfileBTdens_vs_PID_source_meanX << endl;
435 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eProfileBTdens_vs_PID_source_meanY = " << eProfileBTdens_vs_PID_source_meanY << endl;
436 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eProfileBTdens_vs_PID_source_rmsX = " << eProfileBTdens_vs_PID_source_rmsX << endl;
437 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eProfileBTdens_vs_PID_source_rmsY = " << eProfileBTdens_vs_PID_source_rmsY << endl;
438 cout << "-" << endl;
439 cout << "EdbShowPVRQuality::Print() --- SETTINGS: Output data:" << endl;
440 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eAli_modified = " << eAli_modified << endl;
441 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eAli_maxNpatterns = " << eAli_maxNpatterns << endl;
442 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eIsTarget = " << eIsTarget << endl;
443 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eProfileBTdens_vs_PID_target_meanX = " << eProfileBTdens_vs_PID_target_meanX << endl;
444 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eProfileBTdens_vs_PID_target_meanY = " << eProfileBTdens_vs_PID_target_meanY << endl;
445 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eProfileBTdens_vs_PID_target_rmsX = " << eProfileBTdens_vs_PID_target_rmsX << endl;
446 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eProfileBTdens_vs_PID_target_rmsY = " << eProfileBTdens_vs_PID_target_rmsY << endl;
447 cout << "EdbShowPVRQuality::Print() --- " << setw(40) << "eNeedModified = " << eNeedModified << endl;
448 cout << "-" << endl;
449 cout << "EdbShowPVRQuality::Print()....done." << endl;
450 for (int i=0; i<80; ++i) cout << "-";
451 cout << endl;
452
453 Log(2,"EdbShowPVRQuality::Print","EdbShowPVRQuality::Print...done.");
454 return;
455}
Int_t eHistGeometry
Definition: EdbShowPVRQuality.h:42

◆ PrintCutType()

void EdbShowPVRQuality::PrintCutType ( )
396{
397 // Print PrintCutType data for this class.
398 cout << "EdbShowPVRQuality::Print() " << endl;
399 cout << "EdbShowPVRQuality::Print() eCutMethodIsDone[0] " << eCutMethodIsDone[0] << endl;
400 cout << "EdbShowPVRQuality::Print() eCutMethodIsDone[1] " << eCutMethodIsDone[1] << endl;
403 return;
404}
void PrintCutType1()
Definition: EdbShowPVRQuality.cxx:500
void PrintCutType0()
Definition: EdbShowPVRQuality.cxx:458

◆ PrintCutType0()

void EdbShowPVRQuality::PrintCutType0 ( )
459{
460 // Prints quality cut values for each plate of the original EdbPVRec object that were
461 // applied to achieve the basetrack density level, after application the specific cut
462 // on the segments of the specific plate (=pattern).
463
464 if (!eIsSource) return;
465
466 if (!eCutMethodIsDone[0]) {
467 cout << "--- EdbShowPVRQuality::PrintCutType0()----------" << endl;
468 cout << "--- This CutType has not been done. Run Execute_ConstantBTDensity() ---" << endl;
469 cout << "--- to see the results (return now). ---" << endl;
470 return;
471 }
472
473 cout << "----------void EdbShowPVRQuality::PrintCutType0()----------" << endl;
474 cout << "Pattern || Z() || Nseg || eCutp0[i] || eCutp1[i] || BTDensity_orig ||...|| Nseg_modified || BTDensity_modified ||"<< endl;
475
476 if (NULL==eAli_modified) {
478 cout << "WARNING eAli_modified==NULL ==>> Take eAli_orig instead. To build eAli_modified please run CreateEdbPVRec(eAli_orig)" << endl;
479 }
480
481 int Npat_orig = eAli_orig->Npatterns();
482 for (int i=0; i<Npat_orig; i++) {
483 EdbPattern* pat_orig = (EdbPattern*)eAli_orig->GetPattern(i);
484 Int_t npatO=pat_orig->N();
485 EdbPattern* pat_modified = (EdbPattern*)eAli_modified->GetPattern(i);
486 Int_t npatM=pat_modified->N();
487 //
488 cout << i;
489 cout << " ";
490 printf("%.1f %d %.3f %.3f %.2f",pat_orig->Z(),npatO, eCutp0[i], eCutp1[i] , ePatternBTDensity_orig[i]);
491 cout << " ... ";
492 printf("%.1f %d %.3f %.3f %.2f",pat_modified->Z(),npatM, eCutp0[i] ,eCutp1[i], ePatternBTDensity_modified[i]);
493 cout << endl;
494 }
495
496 return;
497}

◆ PrintCutType1()

void EdbShowPVRQuality::PrintCutType1 ( )
501{
502 // Prints quality cut values for each plate of the original EdbPVRec object that were
503 // applied to achieve the basetrack density level, after application the specific cut
504 // on the segments of the specific plate (=pattern).
505
506 if (!eIsSource) return;
507 if (!eCutMethodIsDone[1]) {
508 cout << "--- EdbShowPVRQuality::PrintCutType1()----------" << endl;
509 cout << "--- This CutType has not been done. Run Execute_ConstantQuality() ---" << endl;
510 cout << "--- to see the results (return now). ---" << endl;
511 return;
512 }
513
514 cout << "----------void EdbShowPVRQuality::PrintCutType1()----------" << endl;
515
516 cout << "Pattern || Z() || Nseg || BTDensity_orig || Chi2CutMeanChi2 || Chi2CutRMSChi2 || Chi2CutMeanW || Chi2CutRMSW || Chi2Cut[i] || BTDensity_modified ||"<< endl;
517
518 if (NULL==eAli_modified) {
520 cout << "WARNING eAli_modified==NULL ==>> Take eAli_orig instead. To build eAli_modified please run CreateEdbPVRec(eAli_orig)" << endl;
521 }
522
523 int Npat_orig = eAli_orig->Npatterns();
524 for (int i=0; i<Npat_orig; i++) {
525 EdbPattern* pat_orig = (EdbPattern*)eAli_orig->GetPattern(i);
526 Int_t npatO=pat_orig->N();
527 EdbPattern* pat_modified = (EdbPattern*)eAli_modified->GetPattern(i);
528 Int_t npatM=pat_modified->N();
529 cout << i;
530 cout << " ";
531 printf("%.1f %d %.2f %.2f %.2f %.2f %.2f %.2f",pat_orig->Z(),npatO, eAgreementChi2CutMeanChi2 , eAgreementChi2CutRMSChi2, eAgreementChi2CutMeanW , eAgreementChi2CutRMSW, eAgreementChi2WDistCut[i], ePatternBTDensity_orig[i]);
532 cout << " ... ";
533 printf("%.1f %d %.2f %.2f %.2f %.2f %.2f %.2f",pat_modified->Z(),npatM, eAgreementChi2CutMeanChi2 , eAgreementChi2CutRMSChi2, eAgreementChi2CutMeanW , eAgreementChi2CutRMSW, eAgreementChi2WDistCut[i], ePatternBTDensity_modified[i]);
534 cout << endl;
535 }
536 return;
537}

◆ Remove_DoubleBT()

EdbPVRec * EdbShowPVRQuality::Remove_DoubleBT ( EdbPVRec aliSource)
1215{
1216 // Remove double counted basetracks. This is a scanning induced error, at the
1217 // edges of the objective field, corrections are different, thus one
1218 // basetrack is seen differently in two different views.
1219 // In Bern data, we do see clearly this double bump, that is characterized by
1220 // _very_ close values of two BTs in X,Y,TX,TY (and also Chi2 and W).
1221 // Separation threshold values are 2microns in postion and 10 mrad in angle.
1222 // (again obtained from Bern data).
1223
1224 // Quick and Dirty implementation !
1225 cout << "EdbShowPVRQuality::Remove_DoubleBT()" << endl;
1226 cout << "EdbShowPVRQuality::Remove_DoubleBT() aliSource = " << aliSource << endl;
1227 cout << "EdbShowPVRQuality::Remove_DoubleBT() eAli_orig = " << eAli_orig << endl;
1228
1229 if (NULL==aliSource) {
1230 cout << "WARNING!----EdbShowPVRQuality::Remove_DoubleBT() Source EdbPVRec is NULL. Try to change to object eAli_orig: " << eAli_orig << endl;
1231 if (NULL==eAli_orig) {
1232 cout << "WARNING!----EdbShowPVRQuality::Remove_DoubleBT() Also eAli_orig EdbPVRec is NULL. Do nothing and return NULL pointer!" << endl;
1233 return NULL;
1234 }
1235 else {
1236 aliSource=eAli_orig;
1237 }
1238 }
1239
1240 // Make a new PVRec object anyway
1241 EdbPVRec* aliTarget = new EdbPVRec();
1242 aliTarget->Print();
1243
1244 Bool_t seg_seg_close=kFALSE;
1245 EdbSegP* seg=0;
1246 EdbSegP* seg1=0;
1247 Int_t NdoubleFoundSeg=0;
1248
1249 //gEDBDEBUGLEVEL=3;
1250 cout << "EdbShowPVRQuality::Remove_DoubleBT() Start looping now. Attention, this might take a while!" << endl;
1251
1252 for (int i = 0; i <aliSource->Npatterns(); i++ ) {
1253 if (gEDBDEBUGLEVEL==2) cout << "." << flush;
1254 if (gEDBDEBUGLEVEL>1) cout << "Looping over Ali_source->Pat()=" << i << ". Until now double Candidates found: " << NdoubleFoundSeg << endl;
1255
1256 // Create Target Pattern:
1257 EdbPattern* pat = aliSource->GetPattern(i);
1258 EdbPattern* pt= new EdbPattern();
1259 // SetPattern Values to the parent patterns:
1260 pt->SetID(pat->ID());
1261 pt->SetPID(pat->PID());
1262 pt->SetZ(pat->Z());
1263 // Helper Variable for cout purpose
1264 Int_t nPat=pat->N();
1265
1266 for (int j = 0; j < nPat-1; j++ ) {
1267 if (gEDBDEBUGLEVEL>2) if (j%(nPat/10)==0) cout << "10%more of loop done." << flush;
1268
1269 seg = pat->GetSegment(j);
1270 seg_seg_close=kFALSE;
1271
1272 for (int k = j+1; k <pat->N(); k++ ) {
1273 if (seg_seg_close) continue;
1274
1275 if (gEDBDEBUGLEVEL>3) cout << "Looping over pattern for segment pair nr=" << j << "," << k << endl;
1276 seg1 = pat->GetSegment(k);
1277
1278 // Here decide f.e. which segments to check...
1279 Bool_t toContinue=kFALSE;
1280
1281 // Skip already large distances
1282 if (TMath::Abs(seg->X()-seg1->X())>20.1) continue;
1283 if (TMath::Abs(seg->Y()-seg1->Y())>20.1) continue;
1284 if (TMath::Abs(seg->TX()-seg1->TX())>0.1) continue;
1285 if (TMath::Abs(seg->TY()-seg1->TY())>0.1) continue;
1286
1287 //
1288 // a) all BT parings within a square DeltaR and Delta Theta:
1289 //
1290 // b) all BT pairings within dT<->dR line:
1291 //
1292 Float_t dR_allowed=0;
1293 Float_t deltaR=0;
1294 Float_t deltaTheta=0;
1295 Float_t dR_allowed_m_deltaR=0;
1296
1297 deltaTheta=TMath::Sqrt( TMath::Power(seg->TX()-seg1->TX(),2)+TMath::Power(seg->TY()-seg1->TY(),2) );
1298 if (deltaTheta>0.1) continue;
1299
1300 deltaR=TMath::Sqrt( TMath::Power(seg->X()-seg1->X(),2)+TMath::Power(seg->Y()-seg1->Y(),2) );
1301
1302 // line defined by experimental values:
1303 // dR_allowed = 11(micron)/0.1(rad) * dTheta
1304 dR_allowed = 10.75/0.1*deltaTheta;
1305 dR_allowed_m_deltaR=TMath::Abs(deltaR-dR_allowed);
1306
1307 if (deltaR>3.0||deltaTheta>0.005) toContinue=kTRUE;
1308 if (dR_allowed_m_deltaR<1.0) toContinue=kFALSE;
1309 if (toContinue) continue;
1310
1311 if (gEDBDEBUGLEVEL>3) cout << "EdbPVRQuality::Remove_DoubleBT() Found incompatible segment for dR_dT line condition: " << endl;
1312
1313 // if (gEDBDEBUGLEVEL>3) cout << "EdbPVRQuality::Remove_DoubleBT() Do last check if both are MCEvt segments of different event number! " << endl;
1314 // if ((seg->MCEvt()!=seg1->MCEvt())&&seg->MCEvt()<0) continue;
1315 // I have doubts if I shall take MC events also out, but I guess Yes, because if
1316 // in data these small pairings appear, then they will fall also under the
1317 // category "fake doublet" and then be removed, even if they are real double BT
1318 // from a signal.
1319
1320
1321 ++NdoubleFoundSeg;
1322 seg_seg_close=kTRUE;
1323 // if (seg_seg_close) break;
1324 }
1325 if (seg_seg_close) continue;
1326
1327 // Add segment:
1328 if (gEDBDEBUGLEVEL>3) cout << "// Add segment:" << endl;
1329 pt->AddSegment(*seg);
1330 }
1331 if (gEDBDEBUGLEVEL>2) cout << "// Add Pattern:" << endl;
1332 aliTarget->AddPattern(pt);
1333 }
1334
1335 if (gEDBDEBUGLEVEL>1) aliSource->Print();
1336 if (gEDBDEBUGLEVEL>1) aliTarget->Print();
1337
1338 cout << "EdbShowPVRQuality::Remove_DoubleBT() Statistics: " << endl;
1339 cout << "EdbShowPVRQuality::Remove_DoubleBT() We found " << NdoubleFoundSeg << " double segments too close to each other to be different BTracks." << endl;
1340 cout << "EdbShowPVRQuality::Remove_DoubleBT()...done." << endl;
1341 return aliTarget;
1342}
void Print() const
Definition: EdbPattern.cxx:1693
Float_t TX() const
tangens = deltaX/deltaZ
Definition: EdbSegP.h:175
Float_t TY() const
tangens = deltaY/deltaZ
Definition: EdbSegP.h:176

◆ Remove_Passing()

EdbPVRec * EdbShowPVRQuality::Remove_Passing ( EdbPVRec aliSource)
1347{
1348 // Removes Passing Tracks from the EdbPVRec source object.
1349 // Unfortunately, there does Not Exist an easy function like
1350 // ->RemoveSegment from EdbPVRec object. That makes implementation complicated.
1351
1352 cout << "EdbShowPVRQuality::Remove_Passing()." << endl;
1353 cout << "EdbShowPVRQuality::Remove_Passing() aliSource = " << aliSource << endl;
1354 cout << "EdbShowPVRQuality::Remove_Passing() eAli_orig = " << eAli_orig << endl;
1355
1356 if (NULL==aliSource) {
1357 cout << "WARNING!----EdbShowPVRQuality::Remove_Passing() Source EdbPVRec is NULL. Try to change to object eAli_orig: " << eAli_orig << endl;
1358 if (NULL==eAli_orig) {
1359 cout << "WARNING!----EdbShowPVRQuality::Remove_Passing() Also eAli_orig EdbPVRec is NULL. Do nothing and return NULL pointer!" << endl;
1360 return NULL;
1361 }
1362 else {
1363 aliSource=eAli_orig;
1364 }
1365 }
1366
1367 // Get the tracks from the source object:
1368 TObjArray* Tracks=aliSource->eTracks;
1369
1370 // If eAli_source has no tracks, we return here and stop.
1371 if (NULL == Tracks) {
1372 cout << "WARNING!----EdbShowPVRQuality::Remove_Passing() NULL == eTracks." << endl;
1373 cout << "EdbShowPVRQuality::Remove_Passing() Read tracks from a linked_tracks.root file if there is any." << endl;
1375
1376 if (NULL == Tracks) {
1377 cout << "WARNING!----EdbShowPVRQuality::GetTracksFromLinkedTracksRootFile() NULL == eTracks." << endl;
1378 cout << "WARNING!----EdbShowPVRQuality::Remove_Passing() Reading tracks from a linked_tracks.root file failed." << endl;
1379 cout << "WARNING!----EdbShowPVRQuality::Remove_Passing() Return the original object back." << endl;
1380 return aliSource;
1381 }
1382 }
1383
1384 // Now the object array Tracks should be filled anyway, either from the source EdbPVR
1385 // or from the linked_tracks.root file:
1386 Int_t TracksN=Tracks->GetEntries();
1388 EdbSegP* trackseg;
1389
1390 if (gEDBDEBUGLEVEL>1) {
1391 cout << "EdbShowPVRQuality::Remove_Passing() aliSource/linked_tracks.root has Tracks N=" << TracksN << endl;
1392 }
1393
1394 // Make a new PVRec object anyway
1395 EdbPVRec* aliTarget = new EdbPVRec();
1396 Bool_t seg_in_eTracks=kFALSE;
1397 Int_t totallyRemovedSegments=0;
1398 Int_t totallyAddedSegments=0;
1399 Int_t aliSourceNPat=aliSource->Npatterns();
1400
1401 for (int i = 0; i < aliSourceNPat; i++ ) {
1402 if (gEDBDEBUGLEVEL>2) cout << "Looping over aliSource->Pat()=" << i << endl;
1403
1404 EdbPattern* pat = aliSource->GetPattern(i);
1405 EdbPattern* pt= new EdbPattern();
1406 // SetPattern Values to the parent patterns:
1407 pt->SetID(pat->ID());
1408 pt->SetPID(pat->PID());
1409 pt->SetZ(pat->Z());
1410
1411 for (int j = 0; j <pat->N(); j++ ) {
1412 EdbSegP* seg = pat->GetSegment(j);
1413 seg_in_eTracks=kFALSE;
1414
1415 if (gEDBDEBUGLEVEL>3) cout << "Checking segment j= " << j << " and loop over tracks to check correspondance." <<endl;
1416
1417 for (int ll = 0; ll<TracksN; ll++ ) {
1418 track=(EdbTrackP*)Tracks->At(ll);
1419 Int_t TrackSegN=track->N();
1420
1421 // Compare if this track is regarded as "passing"
1422 // when its less than5 plates than the original pvr object
1423 if (track->Npl()<aliSourceNPat-5) continue;
1424// cout << track->Npl() << " " << aliSourceNPat << endl;
1425
1426 // Here decide f.e. which tracks to check...
1427 // On cosmics it would be nice that f.e. Npl()>=Npatterns-4 ...
1428 // if ....()....
1429 // Since addresses of objects can vary, its better to compare them
1430 // by absolute positions.
1431
1432 for (int kk = 0; kk<TrackSegN; kk++ ) {
1433 if (seg_in_eTracks) continue;
1434 trackseg=track->GetSegment(kk);
1435 if (TMath::Abs(seg->Z()-trackseg->Z())>10.1) continue;
1436 if (TMath::Abs(seg->X()-trackseg->X())>5.1) continue;
1437 if (TMath::Abs(seg->Y()-trackseg->Y())>5.1) continue;
1438 if (TMath::Abs(seg->TX()-trackseg->TX())>0.05) continue;
1439 if (TMath::Abs(seg->TY()-trackseg->TY())>0.05) continue;
1440 //cout << "Found compatible segment!! " << endl;
1441 totallyRemovedSegments++;
1442 seg_in_eTracks=kTRUE;
1443 }
1444 if (seg_in_eTracks) continue;
1445 }
1446 if (seg_in_eTracks) continue;
1447 // Arrived here then the segmen has no equivalent in any segments
1448 // of the tracks array.
1449
1450 // Add segment:
1451 if (gEDBDEBUGLEVEL>3) cout << "EdbShowPVRQuality::Remove_Passing() Add segment:" << endl;
1452 pt->AddSegment(*seg);
1453 totallyAddedSegments++;
1454 }
1455 if (gEDBDEBUGLEVEL>2) cout << "EdbShowPVRQuality::Remove_Passing() TotallyAddedSegments (up to now) = " << totallyAddedSegments << " // Add Pattern:" << endl;
1456 aliTarget->AddPattern(pt);
1457 }
1458
1459 if (gEDBDEBUGLEVEL>1) aliSource->Print();
1460 if (gEDBDEBUGLEVEL>1) aliTarget->Print();
1461
1462
1463
1464 cout << "EdbShowPVRQuality::Remove_Passing() Totally removed segments= " << totallyRemovedSegments << endl;
1465 cout << "EdbShowPVRQuality::Remove_Passing()...done." << endl;
1466 return aliTarget;
1467}
Float_t Z() const
Definition: EdbSegP.h:153
TObjArray * GetTracksFromLinkedTracksRootFile()
Definition: EdbShowPVRQuality.cxx:1619
Definition: bitview.h:14

◆ Remove_Segment()

EdbPVRec * EdbShowPVRQuality::Remove_Segment ( EdbSegP seg,
EdbPVRec aliSource = NULL,
Int_t  Option = 0 
)
1518 {
1519 // Remove one segment (EdbSegP* seg type).
1520 // Return a new EdbPVRec* with previous element removed.
1521 cout << "EdbShowPVRQuality::Remove_Segment():" << endl;
1522 TObjArray* segArray= new TObjArray();
1523 segArray->Add(seg);
1524 Remove_SegmentArray(segArray,aliSource,Option);
1525 delete segArray;
1526 cout << "EdbShowPVRQuality::Remove_Segment()...done." << endl;
1527 return NULL;
1528}
EdbPVRec * Remove_SegmentArray(TObjArray *segarray, EdbPVRec *aliSource=NULL, Int_t Option=0)
Definition: EdbShowPVRQuality.cxx:1474

◆ Remove_SegmentArray()

EdbPVRec * EdbShowPVRQuality::Remove_SegmentArray ( TObjArray *  segarray,
EdbPVRec aliSource = NULL,
Int_t  Option = 0 
)
1475{
1476 // Every "Remove_XXY()" function will call Remove_SegmentArray() since the removal of
1477 // a segment array is easiest to implement. As said earlier, there is no easy function
1478 // to remove tracks from a EdbPVRec object directly. So instead we have to create a
1479 // new EdbPVRec object, fill with all basetracks from the old one, except where the
1480 // segments of the SegmentArray coincide with those which are in the source EdbPVRec.
1481
1482 cout << "EdbShowPVRQuality::Remove_SegmentArray():" << endl;
1483
1484 if (NULL==aliSource) {
1485 cout << "EdbShowPVRQuality::Remove_SegmentArray() ERROR aliSource=NULL pointer. I cannot" << endl;
1486 cout << "EdbShowPVRQuality::Remove_SegmentArray() ERROR do something here and I will return a NULL pointer now. Stop." << endl;
1487 return NULL;
1488 }
1489
1490 if (gEDBDEBUGLEVEL>2) {
1491 cout << "EdbShowPVRQuality::Remove_SegmentArray() INFO " << endl;
1492 cout << "EdbShowPVRQuality::Remove_SegmentArray() segarray = " << segarray << " with " << segarray->GetEntries() << " entries." << endl;
1493 cout << "EdbShowPVRQuality::Remove_SegmentArray() aliSource = " << aliSource << " with " << aliSource->Npatterns() << " patterns." << endl;
1494 cout << "EdbShowPVRQuality::Remove_SegmentArray() Option = " << Option << "." << endl;
1495 }
1496
1497
1498 cout << "EdbShowPVRQuality::Remove_SegmentArray()...done." << endl;
1499 return NULL;
1500}

◆ Remove_Track()

EdbPVRec * EdbShowPVRQuality::Remove_Track ( EdbTrackP track,
EdbPVRec aliSource = NULL,
Int_t  Option = 0 
)
1532 {
1533 // Remove one track (EdbTrackP* track type).
1534 // Return a new EdbPVRec* with previous element removed.
1535 cout << "EdbShowPVRQuality::Remove_Track():" << endl;
1536 TObjArray* segArray= new TObjArray();
1537 segArray= TrackToSegmentArray(track);
1538 Remove_SegmentArray(segArray,aliSource,Option);
1539 delete segArray;
1540 cout << "EdbShowPVRQuality::Remove_Track()...done." << endl;
1541 return NULL;
1542}
TObjArray * TrackToSegmentArray(EdbTrackP *track)
Definition: EdbShowPVRQuality.cxx:1546

◆ Remove_TrackArray()

EdbPVRec * EdbShowPVRQuality::Remove_TrackArray ( TObjArray *  trackArray,
EdbPVRec aliSource = NULL,
Int_t  Option = 0 
)
1504 {
1505 // Remove a whole track array (Array of EdbTrackP type).
1506 // Return a new EdbPVRec* with previous element removed.
1507 cout << "EdbShowPVRQuality::Remove_TrackArray():" << endl;
1508 TObjArray* segArray= new TObjArray();
1509 segArray= TrackArrayToSegmentArray(trackArray);
1510 Remove_SegmentArray(segArray,aliSource,Option);
1511 delete segArray;
1512 cout << "EdbShowPVRQuality::Remove_TrackArray()...done." << endl;
1513 return NULL;
1514}
TObjArray * TrackArrayToSegmentArray(TObjArray *trackArray)
Definition: EdbShowPVRQuality.cxx:1561

◆ Set0()

void EdbShowPVRQuality::Set0 ( )
protected
103{
104 // Reset Values
105
106 Log(2,"EdbShowPVRQuality::Set0","EdbShowPVRQuality::Set0");
107
110 eNeedModified=kFALSE;
111 eIsSource=kFALSE;
112 eIsTarget=kFALSE;
114 for (int i=0; i<2; i++) eCutMethodIsDone[i]=kFALSE;
115
116 // Default BT density level for which the standard cutroutine
117 // will be put:
118 eBTDensityLevel=20; // #BT/mm2
119
120 // Default BT density level will use only segments for
121 // data calculation .
124
125 // Reset Default Geometry: 0 OperaGeometry, 1: MC Geometry
127
128
129
130 eHistYX->Reset();
131 eHistChi2W->Reset();
132
133 for (int i=0; i<114; i++) {
136 eCutp1[i]=0.15;
137 eCutp0[i]=1.0; // Maximum Cut Value for const, BT dens
138 eAgreementChi2WDistCut[i]=3.0; // Maximum Cut Value for const, BT dens
139 }
140
141
150
151 // Default values for cosmics, taken from a brick data:
152 // (I cant remember which one, I hope it doesnt matter).
157
158 Log(2,"EdbShowPVRQuality::Set0","EdbShowPVRQuality::Set0...done.");
159 return;
160}

◆ SetBTDensityLevel()

void EdbShowPVRQuality::SetBTDensityLevel ( Float_t  BTDensityLevel)
inline
96 {
97 eBTDensityLevel=BTDensityLevel;
98 }

◆ SetBTDensityLevelCalcMethodMC()

void EdbShowPVRQuality::SetBTDensityLevelCalcMethodMC ( Bool_t  BTDensityLevelCalcMethodMC)
inline
100 {
101 eBTDensityLevelCalcMethodMC=BTDensityLevelCalcMethodMC;
102 }

◆ SetBTDensityLevelCalcMethodMCConfirmation()

void EdbShowPVRQuality::SetBTDensityLevelCalcMethodMCConfirmation ( Int_t  BTDensityLevelCalcMethodMCConfirmationNumber)
inline
103 {
104 eBTDensityLevelCalcMethodMCConfirmationNumber=BTDensityLevelCalcMethodMCConfirmationNumber;
105 }

◆ SetCutMethod()

void EdbShowPVRQuality::SetCutMethod ( Int_t  CutMethod)
192{
193 // Set Cut Method of the background reduction:
194 // 0: (default): Do cut based on the linear relation: seg.Chi2()<seg.eW*a-b
195 // 1: (testing): Do cut based on a chi2-variable that compares with passing tracks
196 // (if available), i.e. cosmics.
197
198 Log(2,"EdbShowPVRQuality::SetCutMethod","EdbShowPVRQuality::SetCutMethod");
199
200 eCutMethod=CutMethod;
201 cout << "EdbShowPVRQuality::SetCutMethod eCutMethod= " << eCutMethod << endl;
202 if (CutMethod>1) {
203 eCutMethod=0;
204 cout << "WARNING EdbShowPVRQuality::SetCutMethod eCutMethod invalid, Set back to default eCutMethod= " << eCutMethod << endl;
205 }
206
207 Log(2,"EdbShowPVRQuality::SetCutMethod","EdbShowPVRQuality::SetCutMethod...done.");
208 return;
209}

◆ SetEdbPVRec()

void EdbShowPVRQuality::SetEdbPVRec ( EdbPVRec Ali_orig)
inline
144 {
145 eAli_orig=Ali_orig;
146 eIsSource=kTRUE;
147 eAli_maxNpatterns=Ali_orig->Npatterns();
148 }

◆ SetHistGeometry_MC()

void EdbShowPVRQuality::SetHistGeometry_MC ( )
370{
371 // Set the geometry of the basetrack density evaluation using simulation case,
372 // MC ORFEO size conventions: x=-xmax..0..+xmax;y=-ymax..0..ymax).
373 // BinArea is 1mmx1mm
374 eHistYX->Reset();
375 eHistYX->SetBins(100,-50000,50000,100,-50000,50000);
376 //cout << "SetHistGeometry_MC :binwidth (micron)= " << eHistYX->GetBinWidth(1) << endl;
377 return;
378}

◆ SetHistGeometry_OPERA()

void EdbShowPVRQuality::SetHistGeometry_OPERA ( )
359{
360 // Set the geometry of the basetrack density evaluation using OPERA case,
361 // European Scanning System size conventions: x=0..125000;y=0..100000).
362 // BinArea is 1mmx1mm.
363 eHistYX->Reset();
364 eHistYX->SetBins(100,0,100000,120,0,120000);
365 //cout << "SetHistGeometry_OPERA :binwidth (micron)= " << eHistYX->GetBinWidth(1) << endl;
366 return;
367}

◆ SetHistGeometry_OPERAandMC()

void EdbShowPVRQuality::SetHistGeometry_OPERAandMC ( )
381{
382 // Set the geometry of the basetrack density evaluation covering MC and DATA case,
383 // size conventions: x=-125000..0..+125000;y=-125000..0..125000).
384 // BinArea is 1mmx1mm
385
386 Log(3,"EdbShowPVRQuality::SetHistGeometry_OPERAandMC","EdbShowPVRQuality::SetHistGeometry_OPERAandMC");
387 eHistYX->Reset();
388 eHistYX->SetBins(250,-125000,125000,250,-125000,125000);
389 if (gEDBDEBUGLEVEL>2) cout << "SetHistGeometry_OPERAandMC::binwidth (micron)= " << eHistYX->GetBinWidth(1) << endl;
390 Log(3,"EdbShowPVRQuality::SetHistGeometry_OPERAandMC","EdbShowPVRQuality::SetHistGeometry_OPERAandMC...done");
391 return;
392}

◆ TrackArrayToSegmentArray()

TObjArray * EdbShowPVRQuality::TrackArrayToSegmentArray ( TObjArray *  trackArray)
1561 {
1562 // Convert segments of tracks of a track array into an TObjArray containing segments.
1563 TObjArray* segArray= new TObjArray();
1564 Int_t nSegTotal=0;
1565 for (int j=0; j<trackArray->GetEntries(); j++) {
1566 EdbTrackP* track = (EdbTrackP*)trackArray->At(j);
1567 for (int i=0; i<track->N(); i++) {
1568 EdbSegP* seg = track->GetSegment(i);
1569 segArray->Add(seg);
1570 ++nSegTotal;
1571 }
1572 }
1573 if (gEDBDEBUGLEVEL>3) cout << "EdbShowPVRQuality::TrackArrayToSegmentArray() Totally added: " << nSegTotal << " segments from the track array."<<endl;
1574 return segArray;
1575}

◆ TrackToSegmentArray()

TObjArray * EdbShowPVRQuality::TrackToSegmentArray ( EdbTrackP track)
1546 {
1547 // Convert segments of a track into an TObjArray containing segments.
1548 TObjArray* segArray= new TObjArray();
1549 Int_t nSegTotal=0;
1550 for (int i=0; i<track->N(); i++) {
1551 EdbSegP* seg = track->GetSegment(i);
1552 segArray->Add(seg);
1553 ++nSegTotal;
1554 }
1555 if (gEDBDEBUGLEVEL>3) cout << "EdbShowPVRQuality::TrackToSegmentArray() Totally added: " << nSegTotal << " segments from the track."<<endl;
1556 return segArray;
1557}

Member Data Documentation

◆ eAgreementChi2CutMeanChi2

Float_t EdbShowPVRQuality::eAgreementChi2CutMeanChi2
private

◆ eAgreementChi2CutMeanW

Float_t EdbShowPVRQuality::eAgreementChi2CutMeanW
private

◆ eAgreementChi2CutRMSChi2

Float_t EdbShowPVRQuality::eAgreementChi2CutRMSChi2
private

◆ eAgreementChi2CutRMSW

Float_t EdbShowPVRQuality::eAgreementChi2CutRMSW
private

◆ eAgreementChi2WDistCut

Float_t EdbShowPVRQuality::eAgreementChi2WDistCut[114]
private

◆ eAli_maxNpatterns

Int_t EdbShowPVRQuality::eAli_maxNpatterns
private

◆ eAli_modified

EdbPVRec* EdbShowPVRQuality::eAli_modified
private

◆ eAli_orig

EdbPVRec* EdbShowPVRQuality::eAli_orig
private

◆ eBTDensityLevel

Float_t EdbShowPVRQuality::eBTDensityLevel
private

◆ eBTDensityLevelCalcMethodMC

Bool_t EdbShowPVRQuality::eBTDensityLevelCalcMethodMC
private

◆ eBTDensityLevelCalcMethodMCConfirmationNumber

Int_t EdbShowPVRQuality::eBTDensityLevelCalcMethodMCConfirmationNumber
private

◆ eCutDistChi2

Float_t EdbShowPVRQuality::eCutDistChi2[114]
private

◆ eCutDistW

Float_t EdbShowPVRQuality::eCutDistW[114]
private

◆ eCutMethod

Int_t EdbShowPVRQuality::eCutMethod
private

◆ eCutMethodIsDone

Bool_t EdbShowPVRQuality::eCutMethodIsDone[2]
private

◆ eCutp0

Float_t EdbShowPVRQuality::eCutp0[114]
private

◆ eCutp1

Float_t EdbShowPVRQuality::eCutp1[114]
private

◆ eHistChi2W

TH2F* EdbShowPVRQuality::eHistChi2W
private

◆ eHistGeometry

Int_t EdbShowPVRQuality::eHistGeometry
private

◆ eHistYX

TH2F* EdbShowPVRQuality::eHistYX
private

◆ eIsSource

Bool_t EdbShowPVRQuality::eIsSource
private

◆ eIsTarget

Bool_t EdbShowPVRQuality::eIsTarget
private

◆ eNeedModified

Bool_t EdbShowPVRQuality::eNeedModified
private

◆ ePatternBTDensity_modified

Float_t EdbShowPVRQuality::ePatternBTDensity_modified[114]
private

◆ ePatternBTDensity_orig

Float_t EdbShowPVRQuality::ePatternBTDensity_orig[114]
private

◆ eProfileBTdens_vs_PID_source

TProfile* EdbShowPVRQuality::eProfileBTdens_vs_PID_source
private

◆ eProfileBTdens_vs_PID_source_meanX

Float_t EdbShowPVRQuality::eProfileBTdens_vs_PID_source_meanX
private

◆ eProfileBTdens_vs_PID_source_meanY

Float_t EdbShowPVRQuality::eProfileBTdens_vs_PID_source_meanY
private

◆ eProfileBTdens_vs_PID_source_rmsX

Float_t EdbShowPVRQuality::eProfileBTdens_vs_PID_source_rmsX
private

◆ eProfileBTdens_vs_PID_source_rmsY

Float_t EdbShowPVRQuality::eProfileBTdens_vs_PID_source_rmsY
private

◆ eProfileBTdens_vs_PID_target

TProfile* EdbShowPVRQuality::eProfileBTdens_vs_PID_target
private

◆ eProfileBTdens_vs_PID_target_meanX

Float_t EdbShowPVRQuality::eProfileBTdens_vs_PID_target_meanX
private

◆ eProfileBTdens_vs_PID_target_meanY

Float_t EdbShowPVRQuality::eProfileBTdens_vs_PID_target_meanY
private

◆ eProfileBTdens_vs_PID_target_rmsX

Float_t EdbShowPVRQuality::eProfileBTdens_vs_PID_target_rmsX
private

◆ eProfileBTdens_vs_PID_target_rmsY

Float_t EdbShowPVRQuality::eProfileBTdens_vs_PID_target_rmsY
private

◆ maxX

Float_t EdbShowPVRQuality::maxX
private

◆ maxY

Float_t EdbShowPVRQuality::maxY
private

◆ minX

Float_t EdbShowPVRQuality::minX
private

◆ minY

Float_t EdbShowPVRQuality::minY
private

◆ NbinsX

Int_t EdbShowPVRQuality::NbinsX
private

◆ NbinsY

Int_t EdbShowPVRQuality::NbinsY
private

The documentation for this class was generated from the following files: