2 This macro creates for each pattern four 2D histograms: distance between MChit and COG vs impact angle. This information are taken
3 from the output of compClusHits.C. Information about the pattern are taken from the database generated by BuildClTop.C.
4 The macro fits these histograms with gaussians along y axis, in order to study the dependence of dZ and dX wrt
5 the imact angles. Two 1D histograms with dZ and dX are generated as a projection on to the y axis and are used to calculate mean value-sigma
7 The output is a PDF with 4 2D histograms (mean-value of dZ(dX) vs alpha, sigma of dZ(dX) vs alpha,mean-value of dZ(dX) vs beta,
8 sigma of dZ(dX) vs beta) and 2 histograms, dZ and dX, for each pattern, a drawing of the pattern and the information from the fit and the DB and
9 information from DB and fits.
10 You can set the range of patterns you want to study by setting LowerLimit/UpperLimit. If the range coincide with the the whole ragne of patterns,
11 the DB is updated with the information from these fits.
12 It is possible to loop over few clusters instead of all the cluster setting ntotCl: usefull to check the macro. (-1 means all the clusters).
13 It is possilble toprint a PDF with the 2D histograms, pattern drawing and info: setting number2d histos you decide the number of pattern for which
14 you want 2D histograms.
18 #if !defined(__CINT__) || defined(__MAKECINT__)
19 #include "TObjArray.h"
32 #include "TStopwatch.h"
35 #include "../ITS/UPGRADE/AliITSUClusterPix.h"
36 #include "../ITS/UPGRADE/AliITSURecoLayer.h"
37 #include "../ITS/UPGRADE/AliITSURecoDet.h"
38 #include "../ITS/UPGRADE/AliITSUHit.h"
39 #include "../ITS/UPGRADE/AliITSUGeomTGeo.h"
40 #include "AliITSsegmentation.h"
41 #include "AliGeomManager.h"
47 TObjArray* pattDB=0; //it is an array with all the patterns in TBits format (forom the most to the least frequent)
48 TVectorF* pattFR=0; //frequency of the pattern in the database
49 TVectorF* xCentrePix=0; //coordinate of the centre of the pixel containing the COG for the down-left corner in fracion of the pitch
50 TVectorF* zCentrePix=0;
51 TVectorF* xCentreShift=0; //Difference between x coordinate fo COG and the centre of the pixel containing the COG, in fraction of the pitch
52 TVectorF* zCentreShift=0;
53 TVectorF* NPix=0; //Number of fired pixels
54 TVectorF* NRow=0; //Number of rows of the pattern
55 TVectorF* NCol=0; //Number of columns of the pattern
57 TArrayF DeltaZmean; //mean value of the difference between MChit and COG z coordinate in micron derived from a gaussian fit
58 TArrayF DeltaZmeanErr;
59 TArrayF DeltaXmean; //mean value of the difference between MChit and COG x coordinate in micron derived from a gaussian fit
60 TArrayF DeltaXmeanErr;
61 TArrayF DeltaZsigma; //sigma of the difference between MChit and COG z coordinate in micron derived from a gaussian fit
62 TArrayF DeltaZsigmaErr;
63 TArrayF DeltaXsigma; //sigma of the difference between MChit and COG x coordinate in micron derived from a gaussian fit
64 TArrayF DeltaXsigmaErr;
65 TArrayF Chi2x; //Chi2 of the gaussian fit over x
66 TArrayF Chi2z; //Chi2 of the aussian fut over z
75 Float_t pitchX; //pitch of the pixel in x direction
76 Float_t pitchZ; //pitch of the pixel in z direction
78 enum{kDXAlpha=0, kDZAlpha=1, kDXBeta=2, kDZBeta=3};//to select the kind of TH2 we want
80 void LoadDB(const char* fname);
81 void CreateHistos(TObjArray* harr);
82 TH2* GetHisto(Int_t pattID, Int_t kind, TObjArray* harr);
83 void FitHistos(Int_t pattID, Float_t Freq, Int_t nPix, Int_t nRow, Int_t nCol,
84 TObjArray* harr, TCanvas* canvas, TCanvas* canvasbis);
85 void UpDateDB(const char* outDBFile);
107 Float_t alpha; // alpha is the angle in y-radius plane in local frame
108 Float_t beta; // beta is the angle in xz plane, taken from z axis, growing counterclockwise
122 Int_t LowerLimit=0;//set the first cluster of the range of interest, INCLUDED (to update DB LowerLimit=0)
123 Int_t UpperLimit=4000; // EXCLUDED; -1 means number of Patterns. In this case the database is updated with the information about DeltaxMean, DeltaXerr, etc.
125 Int_t ntotCl=-1;//Set -1 to loop over all the clusters, otherwise set the number of clusters
126 Int_t number2dhisto=10;//Set the number pattern from 0 to number2dhisto-1 for which you want to print the TH2 histos (dZ or dX vs alpha or beta)
132 LoadDB("clusterTopology.root");
134 AliGeomManager::LoadGeometry("geometry.root");
135 AliITSUGeomTGeo* gm = new AliITSUGeomTGeo(kTRUE);
136 AliITSUClusterPix::SetGeom(gm);
137 const AliITSsegmentation* segm = gm->GetSegmentation(0);
138 pitchX = segm->Dpx(0)*10000; //pitch in X in micron
139 printf("pitchX: %f\n",pitchX);
140 pitchZ = segm->Dpz(0)*10000; //pitch in Z in micron
141 printf("pitchX: %f\n",pitchX);
144 nPatterns = (pattDB->GetEntriesFast());
147 UpperLimit=nPatterns;
150 nStudied=UpperLimit-LowerLimit+1;
152 printf("\n\n nPatterns %d\n\n",nPatterns);
154 DeltaZmean.Set(nPatterns+100);
155 DeltaZmeanErr.Set(nPatterns+100);
156 DeltaXmean.Set(nPatterns+100);
157 DeltaXmeanErr.Set(nPatterns+100);
158 DeltaZsigma.Set(nPatterns+100);
159 DeltaZsigmaErr.Set(nPatterns+100);
160 DeltaXsigma.Set(nPatterns+100);
161 DeltaXsigmaErr.Set(nPatterns+100);
162 Chi2x.Set(nPatterns+100);
163 Chi2z.Set(nPatterns+100);
164 NDFx.Set(nPatterns+100);
165 NDFz.Set(nPatterns+100);
167 printf("Loading input tree... ");
169 TFile *input = new TFile("clInfo.root","read");
170 TTree *clitsu = (TTree*) input->Get("clitsu");
173 ntotCl = (Int_t) clitsu->GetEntries();
175 TBranch *brdZ = clitsu->GetBranch("dZ");
176 brdZ->SetAddress(&Summ.dZ);
177 TBranch *brdX = clitsu->GetBranch("dX");
178 brdX->SetAddress(&Summ.dX);
179 TBranch *brpattID = clitsu->GetBranch("pattID");
180 brpattID->SetAddress(&Summ.pattID);
181 TBranch *brfreq = clitsu->GetBranch("freq");
182 brfreq->SetAddress(&Summ.freq);
183 TBranch *bralpha = clitsu->GetBranch("alpha");
184 bralpha->SetAddress(&Summ.alpha);
185 TBranch *brbeta = clitsu->GetBranch("beta");
186 brbeta->SetAddress(&Summ.beta);
187 printf("done!!\n\n");
189 CreateHistos(&histoArr);
190 printf("Starting to process clusters\n\n");
192 for(Int_t iCl=0; iCl<ntotCl; iCl++){
193 clitsu->GetEntry(iCl);
194 Int_t ID = Summ.pattID;
195 printf("Processing cluster %d... ",iCl);
196 if(ID<LowerLimit || ID>=UpperLimit) {
197 printf("skipped!\n\n");
200 GetHisto(ID,kDXAlpha,&histoArr)->Fill(Summ.alpha,Summ.dX);
201 GetHisto(ID,kDZAlpha,&histoArr)->Fill(Summ.alpha,Summ.dZ);
202 GetHisto(ID,kDXBeta,&histoArr)->Fill(Summ.beta,Summ.dX);
203 GetHisto(ID,kDZBeta,&histoArr)->Fill(Summ.beta,Summ.dZ);
206 printf("All clusters processed!\n");
210 TCanvas* cnv1 = new TCanvas("cnv1","cnv1",900,700);
211 TCanvas* cnv2 = new TCanvas("cnv2","cnv2",900,700);
213 printf("Printing PDF...\n\n");
216 cnv1->Print("angles.pdf[");
217 cnv2->Print("angles2D.pdf[");
219 for(Int_t i=LowerLimit; i<UpperLimit; i++){
221 FitHistos(i,(*pattFR)[i],(*NPix)[i],(*NRow)[i],(*NCol)[i],&histoArr, cnv1, cnv2);
224 cnv1->Print("angles.pdf]");
225 cnv2->Print("angles2D.pdf]");
227 printf("\n\ndone!!\n\n");
233 if(LowerLimit==0 && UpperLimit==nPatterns){
234 UpDateDB("clusterTopology2.root");
239 void LoadDB(const char* fname){ //load the data base built by BuildClTopDB.C
242 TFile* fl = TFile::Open(fname);
243 if(!fl){printf("Could not find %s",fname); exit(1);}
244 pattDB = (TObjArray*)fl->Get("TopDB");
245 pattFR = (TVectorF*)fl->Get("TopFreq");
246 xCentrePix =(TVectorF*)fl->Get("xCOG");
247 zCentrePix =(TVectorF*)fl->Get("zCOG");
248 xCentreShift =(TVectorF*)fl->Get("xShift");
249 zCentreShift =(TVectorF*)fl->Get("zShift");
250 NPix =(TVectorF*)fl->Get("NPix");
251 NCol =(TVectorF*)fl->Get("NCol");
252 NRow =(TVectorF*)fl->Get("NRow");
255 void CreateHistos(TObjArray* harr){ //creates four 2D histograms (dX or dZ vs alpha or beta) for each pattern included in the study range
257 printf("Creating histograms... ");
259 //for each pattern we crate 4 2D histograms corresponding the cathegories within the enumeration
260 if(!harr) harr=&histoArr;
262 for(Int_t i=0; i<nStudied; i++){
264 TH2* h0 = new TH2F(Form("hXalpha%d", i),
265 Form("#DeltaX vs #alpha for pattID %d", i),10,-0.1, 1.1* TMath::Pi()/2,100,-30,30);
267 harr->AddAtAndExpand(h0, i*4);//This histogram is the first of the 4-histos series
269 TH2* h1 = new TH2F(Form("hZalpha%d", i),
270 Form("#DeltaZ vs #alpha for pattID %d", i),10, -0.1, 1.1* TMath::Pi()/2,100,-30,30);
272 harr->AddAtAndExpand(h1, i*4+1);//This histogram is the second of the 4-histos series and so on...
274 TH2* h2 = new TH2F(Form("hXbeta%d", i),
275 Form("#DeltaX vs #beta for pattID %d", i),10,-0.1, 1.1* TMath::Pi()/2,100,-30,30);
277 harr->AddAtAndExpand(h2, i*4+2);
279 TH2* h3 = new TH2F(Form("hZbeta%d", i),
280 Form("#DeltaZ vs #beta for pattID %d", i),10,-0.1, 1.1* TMath::Pi()/2,100,-30,30);
282 harr->AddAtAndExpand(h3, i*4+3);
286 printf("done!!\n\n");
289 TH2* GetHisto(Int_t pattID, Int_t kind, TObjArray* harr){
291 if(!harr) harr=&histoArr ;
293 h=(TH2*)harr->At((pattID-LowerLimit)*4+kind);
295 printf("Unknown histo id=%d\n",pattID);
302 void FitHistos(Int_t pattID, Float_t Freq, Int_t nPix, Int_t nRow,
303 Int_t nCol, TObjArray* harr, TCanvas* canvas, TCanvas* canvasbis){ //Fit the 2D histograms with gaussian along the Y and creates two 1D histograms,
304 // with dX and dZ, bay a projection to the y axis and the fit with gaussian
306 static TF1* gs = new TF1("gs","gaus",-50,50);
308 if(!harr) harr=&histoArr;
310 //Creating the pads for the output
315 Float_t width = (1 - canvas->GetRightMargin() - canvas->GetLeftMargin() )/3;
317 TPad* titlepad = new TPad("titlepad","",canvas->GetLeftMargin(),0.825,1-canvas->GetRightMargin(),0.95);
318 titlepad->SetFillColor(kYellow);
320 TPad* pad1 = new TPad("pad1","",canvas->GetLeftMargin(),(0.8-canvas->GetBottomMargin())/2+canvas->GetBottomMargin(),canvas->GetLeftMargin()+width,0.8);
322 TPad* pad2 = new TPad("pad2","",canvas->GetLeftMargin()+width,(0.8-canvas->GetBottomMargin())/2+canvas->GetBottomMargin(),canvas->GetLeftMargin()+2*width,0.8);
324 TPad* pad3 = new TPad("pad3","",canvas->GetLeftMargin()+2*width,(0.8-canvas->GetBottomMargin())/2+canvas->GetBottomMargin(),1-canvas->GetRightMargin(),0.8);
326 TPad* pad4 = new TPad("pad3","",canvas->GetLeftMargin(),canvas->GetBottomMargin(),canvas->GetLeftMargin()+width,(0.8-canvas->GetBottomMargin())/2+canvas->GetBottomMargin());
328 TPad* pad5 = new TPad("pad5","",canvas->GetLeftMargin()+width,canvas->GetBottomMargin(),canvas->GetLeftMargin()+2*width,(0.8-canvas->GetBottomMargin())/2+canvas->GetBottomMargin());
330 TPad* pad6 = new TPad("pad6","",canvas->GetLeftMargin()+2*width,canvas->GetBottomMargin(),1-canvas->GetRightMargin(),(0.8-canvas->GetBottomMargin())/2+canvas->GetBottomMargin());
332 TPad* padPattern = new TPad("padPattern","",0.75,0.825,0.9,0.95);
335 //...................................................................................................
340 TH2* hXA = (TH2*)harr->At((pattID-LowerLimit)*4);
341 hXA->GetYaxis()->SetTitle("#DeltaX, #mum");
342 hXA->GetXaxis()->SetTitle("#alpha");
344 hXA->FitSlicesY(gs,0,-1,0,"qnr",&harrFXA);
345 TH1* hXA_1 = (TH1*)harrFXA.At(1);
347 hXA_1->SetDirectory(0);
348 hXA_1->GetYaxis()->SetTitle("#DeltaX_{#mu}(#DeltaZ_{#mu}), #mum");
349 TH1* hXA_2 = (TH1*)harrFXA.At(2);
351 hXA_2->SetDirectory(0);
352 hXA_2->GetYaxis()->SetTitle("#DeltaX_{#sigma}(#DeltaZ_{#sigma}), #mum");
354 //............................................................................................................
360 TH1* hdx = hXA->ProjectionY("hdx");
361 hdx->GetXaxis()->SetTitle("#DeltaX, #mum");
362 hdx->SetTitle(Form("#DeltaX for pattID %d",pattID));
365 gs->SetParameters(hdx->GetMaximum(),hdx->GetMean(),hdx->GetRMS());
366 if((hdx->GetEntries())<100) fitStatusX = hdx->Fit("gs","ql");
367 else fitStatusX = hdx->Fit("gs","q");
370 DeltaXmean[pattID]= gs->GetParameter(1);
371 DeltaXmeanErr[pattID]= gs->GetParError(1);
372 DeltaXsigma[pattID]= gs->GetParameter(2);
373 DeltaXsigmaErr[pattID]= gs->GetParError(2);
374 Chi2x[pattID] = gs->GetChisquare();
375 Int_t varNDFx = gs->GetNDF();
377 NDFx[pattID] = varNDFx;
383 DeltaXmean[pattID]=0;
384 DeltaXmeanErr[pattID]=0;
385 DeltaXsigma[pattID]=0;
386 DeltaXsigmaErr[pattID]=0;
394 //.............................................................................................................
400 TH2* hZA = (TH2*)harr->At((pattID-LowerLimit)*4+1);
401 hZA->GetYaxis()->SetTitle("#DeltaZ, #mum");
402 hZA->GetXaxis()->SetTitle("#alpha");
404 hZA->FitSlicesY(gs,0,-1,0,"qnr",&harrFZA);
405 TH1* hZA_1 = (TH1*)harrFZA.At(1);
406 hZA_1->SetMarkerColor(8);
407 hZA_1->SetLineColor(8);
409 TH1* hZA_2 = (TH1*)harrFZA.At(2);
410 hZA_2->SetMarkerColor(8);
411 hZA_2->SetLineColor(8);
414 Double_t maxA_1 = TMath::Max(hXA_1->GetMaximum(),hZA_1->GetMaximum());
415 Double_t minA_1 = TMath::Min(hXA_1->GetMinimum(),hZA_1->GetMinimum());
416 Double_t addtorangeA_1 = (maxA_1-minA_1)*0.1;
417 Double_t rangeA_1 = addtorangeA_1*10;
419 Double_t maxA_2 = TMath::Max(hXA_2->GetMaximum(),hZA_2->GetMaximum());
420 Double_t minA_2 = TMath::Min(hXA_2->GetMinimum(),hZA_2->GetMinimum());
421 Double_t addtorangeA_2 = (maxA_2-minA_2)*0.1;
422 Double_t rangeA_2 = addtorangeA_2*10;
424 hXA_1->SetMaximum(maxA_1+addtorangeA_1);
425 hXA_1->SetMinimum(minA_1-addtorangeA_1);
427 hXA_2->SetMaximum(maxA_2+addtorangeA_2);
428 hXA_2->SetMinimum(minA_2-addtorangeA_2);
433 TLegend* legA_1 = new TLegend(0.8*1.1* TMath::Pi()/2, maxA_1-0.2*rangeA_1, TMath::Pi()/2, maxA_1,"","");
434 legA_1->AddEntry(hXA_1, " #DeltaX_{#mu}", "l");
435 legA_1->AddEntry(hZA_1, " #DeltaZ_{#mu}", "l");
436 legA_1->SetFillColor(0);
437 legA_1->SetBorderSize(0);
442 TLegend* legA_2 = new TLegend(0.8*1.1* TMath::Pi()/2, maxA_2-0.2*rangeA_2, TMath::Pi()/2, maxA_2,"","");
443 legA_2->AddEntry(hXA_2, " #DeltaX_{#sigma}", "l");
444 legA_2->AddEntry((TObject*)0, "", "");
445 legA_2->AddEntry(hZA_2, " #DeltaZ_{#sigma}", "l");
446 legA_2->SetFillColor(0);
447 legA_2->SetBorderSize(0);
450 //...........................................................................................
457 TH1* hdz = hZA->ProjectionY("hdz");
458 hdz->GetXaxis()->SetTitle("#DeltaZ, #mum");
459 hdz->SetTitle(Form("#DeltaZ for pattID %d",pattID));
462 gs->SetParameters(hdz->GetMaximum(),hdz->GetMean(),hdz->GetRMS());
463 if((hdz->GetEntries())<100) fitStatusZ = hdz->Fit("gs","ql");
464 else fitStatusZ = hdz->Fit("gs","q");
467 DeltaZmean[pattID]= gs->GetParameter(1);
468 DeltaZmeanErr[pattID]= gs->GetParError(1);
469 DeltaZsigma[pattID]= gs->GetParameter(2);
470 DeltaZsigmaErr[pattID]= gs->GetParError(2);
471 Chi2z[pattID] = gs->GetChisquare();
472 Int_t varNDFz = gs->GetNDF();
474 NDFz[pattID] = varNDFz;
480 DeltaZmean[pattID]=0;
481 DeltaZmeanErr[pattID]=0;
482 DeltaZsigma[pattID]=0;
483 DeltaZsigmaErr[pattID]=0;
490 //.............................................................................................................
495 TH2* hXB = (TH2*)harr->At((pattID-LowerLimit)*4+2);
496 hXB->GetYaxis()->SetTitle("#DeltaX, #mum");
497 hXB->GetXaxis()->SetTitle("#beta");
499 hXB->FitSlicesY(gs,0,-1,0,"qnr",&harrFXB);
500 TH1* hXB_1 = (TH1*)harrFXB.At(1);
501 hXB_1->SetDirectory(0);
503 hXB_1->GetYaxis()->SetTitle("#DeltaX_{#mu}(#DeltaZ_{#mu}), #mum");
504 TH1* hXB_2 = (TH1*)harrFXB.At(2);
505 hXB_2->SetDirectory(0);
507 hXB_2->GetYaxis()->SetTitle("#DeltaX_{#sigma}(#DeltaZ_{#sigma}), #mum");
509 //............................................................................................
514 TH2* hZB = (TH2*)harr->At((pattID-LowerLimit)*4+3);
515 hZB->GetYaxis()->SetTitle("#DeltaZ, #mum");
516 hZB->GetXaxis()->SetTitle("#alpha");
518 hZB->FitSlicesY(gs,0,-1,0,"qnr",&harrFZB);
519 TH1* hZB_1 = (TH1*)harrFZB.At(1);
520 hZB_1->SetMarkerColor(8);
521 hZB_1->SetLineColor(8);
523 TH1* hZB_2 = (TH1*)harrFZB.At(2);
524 hZB_2->SetMarkerColor(8);
525 hZB_2->SetLineColor(8);
528 Double_t maxB_1 = TMath::Max(hXB_1->GetMaximum(),hZB_1->GetMaximum());
529 Double_t minB_1 = TMath::Min(hXB_1->GetMinimum(),hZB_1->GetMinimum());
530 Double_t addtorangeB_1 = (maxB_1-minB_1)*0.1;
531 Double_t rangeB_1 = 10*addtorangeB_1;
533 Double_t maxB_2 = TMath::Max(hXB_2->GetMaximum(),hZB_2->GetMaximum());
534 Double_t minB_2 = TMath::Min(hXB_2->GetMinimum(),hZB_2->GetMinimum());
535 Double_t addtorangeB_2 = (maxB_2-minB_2)*0.1;
536 Double_t rangeB_2 = 10*addtorangeB_2;
539 hXB_1->SetMaximum(maxB_1+addtorangeB_1);
540 hXB_1->SetMinimum(minB_1-addtorangeB_1);
542 hXB_2->SetMaximum(maxB_2+addtorangeB_2);
543 hXB_2->SetMinimum(minB_2-addtorangeB_2);
548 TLegend* legB_1 = new TLegend(0.8*1.1* TMath::Pi()/2, maxB_1-0.2*rangeB_1, TMath::Pi()/2, maxB_1,"","");
549 legB_1->AddEntry(hXA_1, " #DeltaX_{#mu}", "l");
550 legB_1->AddEntry((TObject*)0, "", "");
551 legB_1->AddEntry(hZA_1, " #DeltaZ_{#mu}", "l");
552 legB_1->SetFillColor(0);
553 legB_1->SetBorderSize(0);
558 TLegend* legB_2 = new TLegend(0.8*1.1* TMath::Pi()/2, maxB_2-0.2*rangeB_2, TMath::Pi()/2, maxB_2,"","");
559 legB_2->AddEntry(hXB_2, " #DeltaX_{#sigma}", "l");
560 legB_2->AddEntry((TObject*)0, "", "");
561 legB_2->AddEntry(hZB_2, " #DeltaZ_{#sigma}", "l");
562 legB_2->SetFillColor(0);
563 legB_2->SetBorderSize(0);
566 //cnv->Print("patt0.gif");
568 //..................................................................................................................
572 const char* text = Form("pattID: %d freq: %f nPix: %d nRow: %d nCol: %d xCOG (fraction): %0.2f + (%0.2f) zCOG (fraction): %0.2f + (%0.2f)",
573 pattID,Freq,nPix,nRow,nCol, (*xCentrePix)[pattID],(*xCentreShift)[pattID],(*zCentrePix)[pattID],
574 (*zCentreShift)[pattID]);
575 const char* text1 = Form("#DeltaX_{#mu}: (%.2f #pm %.2f) #mum #DeltaX_{#sigma}: (%.2f #pm %.2f) #mum #DeltaX_{MC-centre}: %.2f #mum #chi^{2}/NDF: %.2f / %d",
576 DeltaXmean[pattID],DeltaXmeanErr[pattID],DeltaXsigma[pattID],DeltaXsigmaErr[pattID], DeltaXmean[pattID]+(*xCentreShift)[pattID]*pitchX, Chi2x[pattID], NDFx[pattID]);
577 const char* text2 = Form("#DeltaZ_{#mu}: (%.2f #pm %.2f) #mum #DeltaZ_{#sigma}: (%.2f #pm %.2f) #mum #DeltaZ_{MC-centre}: %.2f #mum #chi^{2}/NDF: %.2f / %d",
578 DeltaZmean[pattID],DeltaZmeanErr[pattID],DeltaZsigma[pattID],DeltaZsigmaErr[pattID], DeltaZmean[pattID]+(*zCentreShift)[pattID]*pitchZ, Chi2z[pattID], NDFz[pattID]);
580 title.SetTextSize(0.1);
581 title.SetTextAlign(12);
583 title.DrawLatex(0.085,0.75,text);
584 title.DrawLatex(0.085,0.5,text1);
585 title.DrawLatex(0.085,0.25,text2);
588 //Drawing the pattern
590 TH2* patternIm = new TH2F("patternIm","", nCol,-0.5,nCol-0.5,nRow,-0.5,nRow-0.5);
591 for (Int_t ir=0; ir<nRow; ir++) {
592 for (Int_t ic=0; ic<nCol; ic++) {
593 if(((TBits*)(*pattDB)[pattID])->TestBitNumber(ir*nCol+ic)){
594 patternIm->Fill(ic,ir);
599 patternIm->SetStats(0);
601 patternIm->Draw("Acol");
603 //Drawing vertical lines on the pattern
605 TObjArray vertlinesarray;
606 for(int i=2; i<=nCol; i++){
607 TLine* vline = new TLine(patternIm->GetXaxis()->GetBinLowEdge(i),patternIm->GetYaxis()->GetBinLowEdge(1),
608 patternIm->GetXaxis()->GetBinLowEdge(i),patternIm->GetYaxis()->GetBinLowEdge(nRow+1));
609 vertlinesarray.AddAtAndExpand(vline,i);
612 for(int i=2; i<=nCol; i++){
613 TLine* v1 = (TLine*)vertlinesarray.At(i);
618 //Drawing horizontal lines on the pattern
620 TObjArray horlinesarray;
621 for(int i=2; i<=nRow; i++){
622 TLine* hline = new TLine(patternIm->GetXaxis()->GetBinLowEdge(1),patternIm->GetYaxis()->GetBinLowEdge(i),
623 patternIm->GetXaxis()->GetBinLowEdge(nCol+1),patternIm->GetYaxis()->GetBinLowEdge(i));
624 horlinesarray.AddAtAndExpand(hline,i);
627 for(int i=2; i<=nRow; i++){
628 TLine* hor1 = (TLine*)horlinesarray.At(i);
633 //Drawing MChit and COG
635 TMarker* COG = new TMarker((*zCentrePix)[pattID]+(*zCentreShift)[pattID]-0.5,
636 (*xCentrePix)[pattID]+(*xCentreShift)[pattID]-0.5,20);
639 TMarker* MChit = new TMarker((*zCentrePix)[pattID]+(*zCentreShift)[pattID]+(DeltaZmean[pattID])/pitchZ-0.5, //20 is the pixel pitch
640 (*xCentrePix)[pattID]+(*xCentreShift)[pattID]+(DeltaXmean[pattID])/pitchX-0.5,29);
641 MChit->SetMarkerColor(kWhite);
644 canvas->Print("angles.pdf");
647 if(pattID<number2dhisto){
650 TPad* titlepadbis = new TPad("titlepadbis","",canvas->GetLeftMargin(),0.825,1-canvas->GetRightMargin(),0.95);
651 titlepadbis->SetFillColor(kYellow);
653 TPad* pad1bis = new TPad("pad1bis","",canvasbis->GetLeftMargin(),(0.8-canvasbis->GetBottomMargin())/2+canvasbis->GetBottomMargin(),0.5,0.8);
655 TPad* pad2bis = new TPad("pad2bis","",0.5,(0.8-canvasbis->GetBottomMargin())/2+canvasbis->GetBottomMargin(),1-canvasbis->GetRightMargin(),0.8);
657 TPad* pad3bis = new TPad("pad3bis","",canvasbis->GetLeftMargin(),canvasbis->GetBottomMargin(),0.5,(0.8-canvasbis->GetBottomMargin())/2+canvasbis->GetBottomMargin());
659 TPad* pad4bis = new TPad("pad4bis","",0.5,canvasbis->GetBottomMargin(),1-canvasbis->GetRightMargin(),(0.8-canvasbis->GetBottomMargin())/2+canvasbis->GetBottomMargin());
661 TPad* padPatternbis = new TPad("padPatternbis","",0.75,0.825,0.9,0.95);
662 padPatternbis->Draw();
665 title.DrawLatex(0.085,0.75,text);
666 title.DrawLatex(0.085,0.5,text1);
667 title.DrawLatex(0.085,0.25,text2);
677 patternIm->Draw("Acol");
678 for(int i=2; i<=nCol; i++){
679 TLine* v1 = (TLine*)vertlinesarray.At(i);
683 for(int i=2; i<=nRow; i++){
684 TLine* hor1 = (TLine*)horlinesarray.At(i);
690 canvasbis->Print("angles2D.pdf");
698 horlinesarray.Delete();
699 vertlinesarray.Delete();
703 void UpDateDB(const char* outDBFile){
705 printf("\n\n\nStoring data in the DataBase\n\n\n");
707 TString outDBFileName = outDBFile;
708 if (outDBFileName.IsNull()) outDBFileName = "clusterTopology2.root";
710 TVectorF arrDeltaZmean(nPatterns);
711 TVectorF arrDeltaZmeanErr(nPatterns);
712 TVectorF arrDeltaXmean(nPatterns);
713 TVectorF arrDeltaXmeanErr(nPatterns);
714 TVectorF arrDeltaZsigma(nPatterns);
715 TVectorF arrDeltaZsigmaErr(nPatterns);
716 TVectorF arrDeltaXsigma(nPatterns);
717 TVectorF arrDeltaXsigmaErr(nPatterns);
718 TVectorF arrChi2x(nPatterns);
719 TVectorF arrNDFx(nPatterns);
720 TVectorF arrChi2z(nPatterns);
721 TVectorF arrNDFz(nPatterns);
723 for(Int_t ID=0; ID<nPatterns; ID++){
725 printf("Processing pattern %d... ", ID);
727 arrDeltaZmean[ID]=DeltaZmean[ID];
728 arrDeltaZmeanErr[ID]=DeltaZmeanErr[ID];
729 arrDeltaXmean[ID]=DeltaXmean[ID];
730 arrDeltaXmeanErr[ID]=DeltaXmeanErr[ID];
731 arrDeltaZsigma[ID]=DeltaZsigma[ID];
732 arrDeltaZsigmaErr[ID]=DeltaZsigmaErr[ID];
733 arrDeltaXsigma[ID]=DeltaXsigma[ID];
734 arrDeltaXsigmaErr[ID]=DeltaXsigmaErr[ID];
735 arrChi2x[ID]=Chi2x[ID];
736 arrNDFx[ID]=NDFx[ID];
737 arrChi2z[ID]=Chi2z[ID];
738 arrNDFz[ID]=NDFz[ID];
740 printf("done!!\n\n");
743 TFile* flDB = TFile::Open(outDBFileName.Data(), "recreate");
744 flDB->WriteObject(pattDB,"TopDB","kSingleKey");
745 flDB->WriteObject(NPix,"NPix","kSingleKey");
746 flDB->WriteObject(NRow,"NRow","kSingleKey");
747 flDB->WriteObject(NCol,"NCol","kSingleKey");
748 flDB->WriteObject(pattFR, "TopFreq","kSingleKey");
749 flDB->WriteObject(xCentrePix,"xCOG","kSingleKey");
750 flDB->WriteObject(xCentreShift,"xShift","kSingleKey");
751 flDB->WriteObject(zCentrePix,"zCOG","kSingleKey");
752 flDB->WriteObject(zCentreShift,"zShift","kSingleKey");
753 flDB->WriteObject(&arrDeltaZmean,"DeltaZmean","kSingleKey");
754 flDB->WriteObject(&arrDeltaZmeanErr,"DeltaZmeanErr","kSingleKey");
755 flDB->WriteObject(&arrDeltaXmean,"DeltaXmean","kSingleKey");
756 flDB->WriteObject(&arrDeltaXmeanErr,"DeltaXmeanErr","kSingleKey");
757 flDB->WriteObject(&arrDeltaZsigma,"DeltaZsigma","kSingleKey");
758 flDB->WriteObject(&arrDeltaZsigmaErr,"DeltaZsigmaErr","kSingleKey");
759 flDB->WriteObject(&arrDeltaXsigma,"DeltaXsigma","kSingleKey");
760 flDB->WriteObject(&arrDeltaXsigmaErr,"DeltaXsigmaErr","kSingleKey");
761 flDB->WriteObject(&arrChi2x,"Chi2x","kSingleKey");
762 flDB->WriteObject(&arrChi2z,"Chi2z","kSingleKey");
763 flDB->WriteObject(&arrNDFx,"NDFx","kSingleKey");
764 flDB->WriteObject(&arrNDFz,"NDFz","kSingleKey");
769 printf("\n\nDB Complete!!\n\n");