1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 /////////////////////////////////////////////////////////////////////////////////
20 // AliTRDCalibraVector
22 // This class is for the vector methode of the TRD calibration.
25 // R. Bailhache (R.Bailhache@gsi.de)
27 //////////////////////////////////////////////////////////////////////////////////////
29 #include <TGraphErrors.h>
32 #include <TObjArray.h>
35 #include <TDirectory.h>
41 #include "AliTRDCalibraVector.h"
42 #include "AliTRDCommonParam.h"
44 ClassImp(AliTRDCalibraVector)
46 //______________________________________________________________________________________
47 AliTRDCalibraVector::AliTRDCalibraVector()
49 ,fVectorPH(new TObjArray())
50 ,fPlaPH(new TObjArray())
51 ,fVectorCH(new TObjArray())
52 ,fPlaCH(new TObjArray())
53 ,fVectorPRF(new TObjArray())
54 ,fPlaPRF(new TObjArray())
60 // Default constructor
65 //______________________________________________________________________________________
66 AliTRDCalibraVector::AliTRDCalibraVector(const AliTRDCalibraVector &c)
68 ,fVectorPH(new TObjArray())
69 ,fPlaPH(new TObjArray())
70 ,fVectorCH(new TObjArray())
71 ,fPlaCH(new TObjArray())
72 ,fVectorPRF(new TObjArray())
73 ,fPlaPRF(new TObjArray())
84 //____________________________________________________________________________________
85 AliTRDCalibraVector::~AliTRDCalibraVector()
88 // AliTRDCalibraVector destructor
93 //_____________________________________________________________________________
94 Int_t AliTRDCalibraVector::SearchBin(Float_t value, Int_t i) const
102 Int_t fbinmax = (Int_t) value;
103 Int_t fNumberOfBin = -1;
109 fNumberOfBin = fNumberBinCharge;
116 fNumberOfBin = fNumberBinPRF;
120 if ((value >= fbinmax) ||
126 reponse = (Int_t) ((fNumberOfBin*(value-fbinmin)) / (fbinmax-fbinmin));
133 //_____________________________________________________________________________
134 Int_t AliTRDCalibraVector::SearchInVector(Int_t group, Int_t i) const
137 // Search if the calibration group "group" has already been
138 // initialised by a previous track in the vector
142 for (Int_t k = 0; k < (Int_t) fPlaCH->GetEntriesFast(); k++) {
143 if (((AliTRDPlace *) fPlaCH->At(k))->GetPlace() == group) {
151 for (Int_t k = 0; k < (Int_t) fPlaPH->GetEntriesFast(); k++) {
152 if (((AliTRDPlace *) fPlaPH->At(k))->GetPlace() == group) {
160 for (Int_t k = 0; k < (Int_t) fPlaPRF->GetEntriesFast(); k++) {
161 if (((AliTRDPlace *) fPlaPRF->At(k))->GetPlace() == group) {
172 //_____________________________________________________________________________
173 Int_t AliTRDCalibraVector::SearchInTreeVector(TObjArray *vectorplace, Int_t group) const
176 // Search if the calibration group "group" is present in the tree
179 for (Int_t k = 0; k < (Int_t) vectorplace->GetEntriesFast(); k++) {
180 if (((AliTRDPlace *) vectorplace->At(k))->GetPlace() == group) {
189 //_____________________________________________________________________________
190 Bool_t AliTRDCalibraVector::UpdateVectorCH(Int_t group, Float_t value)
193 // Fill the vector if a new calibration group "group" or update the
194 // values of the calibration group "group" if already here
198 Int_t bin = SearchBin(value,0);
200 if ((bin < 0) || (bin >= fNumberBinCharge)) {
205 Int_t place = SearchInVector(group,0);
209 AliTRDPlace *placegroup = new AliTRDPlace();
210 placegroup->SetPlace(group);
211 fPlaCH->Add((TObject *) placegroup);
213 AliTRDCTInfo *fCHInfo = new AliTRDCTInfo();
214 UShort_t *entries = new UShort_t[fNumberBinCharge];
216 for(Int_t k = 0; k < fNumberBinCharge; k++) {
222 fCHInfo->SetEntries(entries);
224 fVectorCH->Add((TObject *) fCHInfo);
226 // Group already exits
229 AliTRDCTInfo *fCHInfo = new AliTRDCTInfo();
231 fCHInfo = ((AliTRDCTInfo *) fVectorCH->At(place));
232 UShort_t *entries = fCHInfo->GetEntries();
236 fCHInfo->SetEntries(entries);
238 fVectorCH->AddAt((TObject *) fCHInfo,place);
245 //_____________________________________________________________________________
246 Bool_t AliTRDCalibraVector::UpdateVectorPRF(Int_t group, Float_t x, Float_t y)
249 // Fill the vector if a new calibration group "group" or update the
250 // values of the calibration group "group" if already here
254 Int_t bin = SearchBin(x,2);
256 if ((bin < 0) || (bin >= fNumberBinPRF)) {
261 Int_t place = SearchInVector(group,2);
266 AliTRDPlace *placegroup = new AliTRDPlace();
267 placegroup->SetPlace(group);
268 fPlaPRF->Add((TObject *) placegroup);
269 AliTRDPInfo *fPRFInfo = new AliTRDPInfo();
271 Float_t *sum = new Float_t[fNumberBinPRF];
272 Float_t *sumsquare = new Float_t[fNumberBinPRF];
273 UShort_t *entries = new UShort_t[fNumberBinPRF];
276 for (Int_t k = 0; k < fNumberBinPRF; k++) {
284 sumsquare[bin] += y*y;
288 fPRFInfo->SetSum(sum);
289 fPRFInfo->SetSumSquare(sumsquare);
290 fPRFInfo->SetEntries(entries);
293 fVectorPRF->Add((TObject *) fPRFInfo);
296 // Group already exits
299 AliTRDPInfo *fPRFInfo = new AliTRDPInfo();
301 fPRFInfo = (AliTRDPInfo *) fVectorPRF->At(place);
303 Float_t *sum = fPRFInfo->GetSum();
304 Float_t *sumsquare = fPRFInfo->GetSumSquare();
305 UShort_t *entries = fPRFInfo->GetEntries();
308 Double_t calcul = (((Double_t) fPRFInfo->GetEntries()[bin])
309 * ((Double_t) fPRFInfo->GetSum()[bin]) + (Double_t) y)
310 / (((Double_t) fPRFInfo->GetEntries()[bin]) + 1);
311 sum[bin] = (Float_t) calcul;
312 Double_t calculsquare = (((Double_t) fPRFInfo->GetSumSquare()[bin])
313 * ((Double_t) fPRFInfo->GetEntries()[bin]) + ((Double_t) y)*((Double_t) y))
314 / (((Double_t) fPRFInfo->GetEntries()[bin]) + 1);
315 sumsquare[bin] = (Float_t) calculsquare;
319 fPRFInfo->SetSum(sum);
320 fPRFInfo->SetSumSquare(sumsquare);
321 fPRFInfo->SetEntries(entries);
324 fVectorPRF->AddAt((TObject *) fPRFInfo,place);
332 //_____________________________________________________________________________
333 Bool_t AliTRDCalibraVector::UpdateVectorPH(Int_t group, Int_t time, Float_t value)
336 // Fill the vector if a new calibration group "group" or update
337 // the values of the calibration group "group" if already here
349 Int_t place = SearchInVector(group,1);
354 AliTRDPlace *placegroup = new AliTRDPlace();
355 placegroup->SetPlace(group);
356 fPlaPH->Add((TObject *) placegroup);
357 AliTRDPInfo *fPHInfo = new AliTRDPInfo();
359 Float_t *sum = new Float_t[fTimeMax];
360 Float_t *sumsquare = new Float_t[fTimeMax];
361 UShort_t *entries = new UShort_t[fTimeMax];
364 for (Int_t k = 0; k < fTimeMax; k++) {
372 sumsquare[bin] += value*value;
376 fPHInfo->SetSum(sum);
377 fPHInfo->SetSumSquare(sumsquare);
378 fPHInfo->SetEntries(entries);
381 fVectorPH->Add((TObject *) fPHInfo);
384 // Group already exits
387 AliTRDPInfo *fPHInfo = new AliTRDPInfo();
389 fPHInfo = (AliTRDPInfo *) fVectorPH->At(place);
391 Float_t *sum = fPHInfo->GetSum();
392 Float_t *sumsquare = fPHInfo->GetSumSquare();
393 UShort_t *entries = fPHInfo->GetEntries();
396 Double_t calcul = (((Double_t) fPHInfo->GetEntries()[bin])
397 * ((Double_t) fPHInfo->GetSum()[bin]) + (Double_t) value)
398 / (((Double_t) fPHInfo->GetEntries()[bin]) + 1);
399 sum[bin] = (Float_t) calcul;
400 Double_t calculsquare = ((((Double_t) fPHInfo->GetSumSquare()[bin])
401 * ((Double_t) fPHInfo->GetEntries()[bin]))
402 + (((Double_t) value) * ((Double_t)value)))
403 / (((Double_t) fPHInfo->GetEntries()[bin]) + 1);
404 sumsquare[bin] = (Float_t) calculsquare;
408 fPHInfo->SetSum(sum);
409 fPHInfo->SetSumSquare(sumsquare);
410 fPHInfo->SetEntries(entries);
413 fVectorPH->AddAt((TObject *) fPHInfo,place);
421 //_____________________________________________________________________________
422 TGraphErrors *AliTRDCalibraVector::ConvertVectorPHisto(Int_t place
423 , const Char_t * name) const
426 // Convert the PInfo in a TGraphErrors
429 AliTRDPInfo *pInfo = new AliTRDPInfo();
431 pInfo = ((AliTRDPInfo *) fVectorPH->At(place));
433 histo = ConvertVectorPHistoI((AliTRDPInfo *)pInfo, name);
439 //_____________________________________________________________________________
440 TGraphErrors *AliTRDCalibraVector::ConvertVectorPHistoI(AliTRDPInfo *pInfo
441 , const Char_t *name) const
444 // Convert the PInfo in a 1D grapherror, name must contains "PRF"
445 // if PRF calibration and not "PRF" for Vdrift calibration
449 const Char_t *pattern1 = "PRF";
452 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
454 AliInfo("Could not get CommonParam, take the default 10MHz");
457 sf = parCom->GetSamplingFrequency();
469 if (strstr(name,pattern1)) {
470 ntimes = fNumberBinPRF;
475 x = new Double_t[ntimes]; // Xaxis
476 y = new Double_t[ntimes]; // Mean
477 ex = new Double_t[ntimes]; // Nentries
478 ey = new Double_t[ntimes]; // Sum of square/nentries
481 if (!strstr(name,pattern1)) {
486 step = (1.0 - (-1.0)) / fNumberBinPRF;
487 min = -1.0 + step / 2.0;
491 for (Int_t k = 0; k < ntimes; k++) {
496 // Fill only if there is more than 0 something
497 if (pInfo->GetEntries()[k] > 0) {
498 ex[k] = pInfo->GetEntries()[k];
499 y[k] = pInfo->GetSum()[k];
500 ey[k] = pInfo->GetSumSquare()[k];
504 // Define the TGraphErrors
505 histo = new TGraphErrors(ntimes,x,y,ex,ey);
506 histo->SetTitle(name);
511 //_____________________________________________________________________________
512 TH1F *AliTRDCalibraVector::ConvertVectorCTHisto(Int_t place
513 , const Char_t * name) const
516 // Convert the CTInfo in a 1D histo
519 AliTRDCTInfo *cHInfo = new AliTRDCTInfo();
521 cHInfo = ((AliTRDCTInfo *) fVectorCH->At(place));
523 histo = ConvertVectorCTHistoI((AliTRDCTInfo *)cHInfo,(const Char_t *) name);
529 //_____________________________________________________________________________
530 TH1F *AliTRDCalibraVector::ConvertVectorCTHistoI(AliTRDCTInfo *cTInfo
531 , const Char_t * name) const
534 // Convert the CTInfo in a 1D histo
539 Int_t ntimes = fNumberBinCharge;
540 UShort_t *entries = cTInfo->GetEntries();
543 histo = new TH1F(name,name,fNumberBinCharge,0,300);
546 for (Int_t k = 0; k < ntimes; k++) {
547 histo->SetBinContent(k+1,entries[k]);
548 histo->SetBinError(k+1,TMath::Sqrt(TMath::Abs(entries[k])));
555 //_____________________________________________________________________________
556 TTree *AliTRDCalibraVector::ConvertVectorCTTreeHisto(TObjArray *vVectorCT
559 , const Char_t *nametitle) const
562 // Convert the vector in a tree with two branchs: the group number
563 // and the TH1F histo reconstructed from the vector
566 // Size of the things
567 Int_t ntotal = (Int_t) pPlaCT->GetEntriesFast();
569 AliInfo("nothing to write!");
570 TTree *treeCT = new TTree(name,nametitle);
574 // Variable of the tree
575 Int_t groupnumber = -1; // Group calibration
577 TObjArray vectorCT = *vVectorCT;
578 TObjArray plaCT = *pPlaCT;
581 TTree *treeCT = new TTree(name,nametitle);
582 treeCT->Branch("groupnumber",&groupnumber,"groupnumber/I");
583 treeCT->Branch("histo","TH1F",&histo,32000,0);
589 groupnumber = ((AliTRDPlace *) plaCT.At(0))->GetPlace();
591 histo = ConvertVectorCTHistoI(((AliTRDCTInfo *) vectorCT.At(0)),nome);
593 vectorCT.RemoveAt(0);
604 //_____________________________________________________________________________
605 TTree *AliTRDCalibraVector::ConvertVectorPTreeHisto(TObjArray *vVectorP
608 , const Char_t *nametitle) const
611 // Convert the vector in a tree with two branchs: the group number
612 // and the TGraphErrors histo reconstructed from the vector.
613 // The name must contain "PRF" for PRF calibration and not "PRF"
614 // for Vdrift calibration
617 // Size of the things
618 Int_t ntotal = (Int_t) pPlaP->GetEntriesFast();
620 AliInfo("nothing to write!");
621 TTree *treeP = new TTree(name,nametitle);
625 // Variable of the tree
626 Int_t groupnumber = -1; // Group calibration
627 TGraphErrors *histo = 0x0;
628 TObjArray vectorP = *vVectorP;
629 TObjArray plaP = *pPlaP;
632 TTree *treeP = new TTree(name,nametitle);
633 treeP->Branch("groupnumber",&groupnumber,"groupnumber/I");
634 treeP->Branch("histo","TGraphErrors",&histo,32000,0);
640 groupnumber = ((AliTRDPlace *) plaP.At(0))->GetPlace();
642 histo = ConvertVectorPHistoI((AliTRDPInfo *) vectorP.At(0),nome);
655 //_____________________________________________________________________________
656 TObjArray *AliTRDCalibraVector::ConvertTreeVector(TTree *tree) const
659 // Convert the branch groupnumber of the tree taken from
660 // TRD.calibration.root in case of vector method in a std::vector
665 TObjArray *vectorplace = new TObjArray();
667 // Variable of the tree
668 Int_t groupnumber = -1; // Group calibration
671 tree->SetBranchAddress("groupnumber",&groupnumber);
674 Int_t ntotal = tree->GetEntries();
675 for (Int_t k = 0; k < ntotal; k++) {
677 AliTRDPlace *placegroupnumber = new AliTRDPlace();
678 placegroupnumber->SetPlace(groupnumber);
679 vectorplace->Add((TObject *) placegroupnumber);
686 //_____________________________________________________________________________
687 Bool_t AliTRDCalibraVector::MergeVectorCT(TObjArray *vVectorCT2, TObjArray *pPlaCT2)
690 // Add the two vectors and place the result in the first
693 if (((Int_t) pPlaCT2->GetEntriesFast()) != ((Int_t) vVectorCT2->GetEntriesFast())){
694 AliInfo("VectorCT2 doesn't correspond to PlaCT2!");
699 for (Int_t k = 0; k < (Int_t) fPlaCH->GetEntriesFast(); k++) {
701 // Look if PlaCT1[k] it is also in the second vector
703 for (Int_t j = 0; j < (Int_t) pPlaCT2->GetEntriesFast(); j++) {
704 if (((AliTRDPlace *) pPlaCT2->At(j))->GetPlace() ==
705 ((AliTRDPlace *) fPlaCH->At(k))->GetPlace()) {
711 // If not in the second vector nothing to do
713 // If in the second vector
716 AliTRDCTInfo *fCTInfo = new AliTRDCTInfo();
717 UShort_t *entries = new UShort_t[fNumberBinCharge];
719 for (Int_t nu = 0; nu < fNumberBinCharge; nu++) {
720 entries[nu] = ((AliTRDCTInfo *) fVectorCH->At(((AliTRDPlace *) fPlaCH->At(k))->GetPlace()))->GetEntries()[nu]
721 + ((AliTRDCTInfo *) vVectorCT2->At(((AliTRDPlace *) fPlaCH->At(k))->GetPlace()))->GetEntries()[nu];
725 fCTInfo->SetEntries(entries);
727 // Nothing to do on PlaCT1
730 fVectorCH->AddAt((TObject *) fCTInfo,((AliTRDPlace *) fPlaCH->At(k))->GetPlace());
736 // And at the end the vector in CT2 but not in CH1
737 for (Int_t k = 0; k < (Int_t) pPlaCT2->GetEntriesFast(); k++) {
739 // Look if pPlaCT2[k] it is also in the second vector
741 for (Int_t j = 0; j < (Int_t) fPlaCH->GetEntriesFast(); j++) {
742 if (((AliTRDPlace *) fPlaCH->At(j))->GetPlace() == ((AliTRDPlace *) pPlaCT2->At(k))->GetPlace()) {
748 // If not in the first vector
751 AliTRDCTInfo *fCTInfo = new AliTRDCTInfo();
752 fCTInfo = ((AliTRDCTInfo *) vVectorCT2->At(((AliTRDPlace *) pPlaCT2->At(k))->GetPlace()));
755 fPlaCH->Add((TObject *) (pPlaCT2->At(k)));
756 fVectorCH->Add((TObject *) fCTInfo);
766 //_____________________________________________________________________________
767 Bool_t AliTRDCalibraVector::MergeVectorP(TObjArray *vVectorP2
772 // Add the two vectors and place the result in the first
775 if (((Int_t) pPlaP2->GetEntriesFast()) != ((Int_t) vVectorP2->GetEntriesFast())) {
776 AliInfo("VectorP2 doesn't correspond to PlaP2!");
783 for (Int_t k = 0; k < (Int_t) fPlaPH->GetEntriesFast(); k++) {
785 // Look if fPlaPH[k] it is also in the second vector
787 for (Int_t j = 0; j < (Int_t) pPlaP2->GetEntriesFast(); j++) {
788 if (((AliTRDPlace *) pPlaP2->At(j))->GetPlace() == ((AliTRDPlace *) fPlaPH->At(k))->GetPlace()) {
794 // If not in the second vector nothing to do
796 // If in the second vector
799 AliTRDPInfo *fPInfo = new AliTRDPInfo();
800 UShort_t *entries = new UShort_t[fTimeMax];
801 Float_t *sum = new Float_t[fTimeMax];
802 Float_t *sumsquare = new Float_t[fTimeMax];
804 for (Int_t nu = 0; nu < fTimeMax; nu++) {
806 entries[nu] = ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]
807 + ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu];
809 Double_t calcul = ((((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSum()[nu])
810 * ((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]))
811 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSum()[nu])
812 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu])))
813 / ((Double_t) fPInfo->GetEntries()[nu]);
815 sum[nu] = (Float_t) calcul;
817 Double_t calculsquare = ((((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSumSquare()[nu])
818 * ((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]))
819 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSumSquare()[nu])
820 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu])))
821 / ((Double_t) fPInfo->GetEntries()[nu]);
823 sumsquare[nu] = calculsquare;
829 fPInfo->SetSumSquare(sumsquare);
830 fPInfo->SetEntries(entries);
832 // Nothing to do on PlaCT1
834 // Update the vector VectorCT1
835 fVectorPH->AddAt((TObject *) fPInfo,((AliTRDPlace *) fPlaPH->At(k))->GetPlace());
841 // And at the end the vector in P2 but not in CH1
842 for (Int_t k = 0; k < (Int_t) pPlaP2->GetEntriesFast(); k++) {
844 // Look if PlaCT2[k] it is also in the second vector
846 for (Int_t j = 0; j < (Int_t) fPlaPH->GetEntriesFast(); j++) {
847 if (((AliTRDPlace *) fPlaPH->At(j))->GetPlace() == ((AliTRDPlace *) pPlaP2->At(k))->GetPlace()) {
853 // If not in the first vector
856 AliTRDPInfo *fPInfo = new AliTRDPInfo();
857 fPInfo = (AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) pPlaP2->At(k))->GetPlace());
859 // Add at the end of CH1
860 fPlaPH->Add(((TObject *) pPlaP2->At(k)));
861 fVectorPH->Add((TObject *) fPInfo);
872 for (Int_t k = 0; k < (Int_t) fPlaPRF->GetEntriesFast(); k++) {
874 // Look if fPlaPRF[k] it is also in the second vector
876 for (Int_t j = 0; j < (Int_t) pPlaP2->GetEntriesFast(); j++) {
877 if (((AliTRDPlace *) pPlaP2->At(j))->GetPlace() == ((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()) {
883 // If not in the second vector nothing to do
885 // If in the second vector
888 AliTRDPInfo *fPInfo = new AliTRDPInfo();
889 UShort_t *entries = new UShort_t[fNumberBinPRF];
890 Float_t *sum = new Float_t[fNumberBinPRF];
891 Float_t *sumsquare = new Float_t[fNumberBinPRF];
893 for (Int_t nu = 0; nu < fNumberBinPRF; nu++) {
895 entries[nu] = ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]
896 + ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu];
898 Double_t calcul = ((((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSum()[nu])
899 * ((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]))
900 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSum()[nu])
901 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu])))
902 / ((Double_t) fPInfo->GetEntries()[nu]);
904 sum[nu] = (Float_t) calcul;
906 Double_t calculsquare = ((((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSumSquare()[nu])
907 * ((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]))
908 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSumSquare()[nu])
909 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu])))
910 / ((Double_t) fPInfo->GetEntries()[nu]);
912 sumsquare[nu] = calculsquare;
918 fPInfo->SetSumSquare(sumsquare);
919 fPInfo->SetEntries(entries);
921 // Nothing to do on PlaCT1
923 // Update the vector VectorCT1
924 fVectorPRF->AddAt((TObject *) fPInfo,((AliTRDPlace *) fPlaPRF->At(k))->GetPlace());
930 // And at the end the vector in P2 but not in CH1
931 for (Int_t k = 0; k < (Int_t) pPlaP2->GetEntriesFast(); k++) {
933 // Look if PlaCT2[k] it is also in the second vector
935 for (Int_t j = 0; j < (Int_t) fPlaPRF->GetEntriesFast(); j++) {
936 if (((AliTRDPlace *) fPlaPRF->At(j))->GetPlace() == ((AliTRDPlace *) pPlaP2->At(k))->GetPlace()) {
942 // If not in the first vector
945 AliTRDPInfo *fPInfo = new AliTRDPInfo();
946 fPInfo = (AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) pPlaP2->At(k))->GetPlace());
948 // Add at the end of CH1
949 fPlaPRF->Add(((TObject *) pPlaP2->At(k)));
950 fVectorPRF->Add((TObject *) fPInfo);
962 //_____________________________________________________________________________
963 TTree *AliTRDCalibraVector::Sum2Trees(const Char_t *filename1
964 , const Char_t *filename2
965 , const Char_t *variablecali)
968 // It returns the sum of two trees with the name variablecali
969 // in the files filenam1 and filename2 equivalent of merging two 2D histos
970 // The name of the resulting tree is the same as the two input trees
971 // variablecali can be treeCH2d, treePH2d or treePRF2d
975 TChain *treeChain = new TChain(variablecali);
976 TObjArray *vectorplace = new TObjArray();
977 TObjArray *where = new TObjArray();
981 TFile *file1 = new TFile(filename1,"READ");
982 TTree *tree1 = (TTree *) file1->Get(variablecali);
987 vectorplace = ConvertTreeVector(tree1);
989 // Say where it is in tree 1
990 for (Int_t jui = 0; jui < (Int_t) vectorplace->GetEntriesFast(); jui++) {
991 AliTRDPlace *placejui = new AliTRDPlace();
992 placejui->SetPlace(jui);
993 TObjArray *chainplace = new TObjArray();
994 chainplace->Add((TObject *) placejui);
995 where->Add((TObject *) chainplace);
999 treeChain->Add(filename1);
1004 TFile *file2 = new TFile(filename2,"READ");
1005 TTree *tree2 = (TTree *) file2->Get(variablecali);
1010 TObjArray *vector2 = ConvertTreeVector(tree2);
1011 Int_t j = treeChain->GetEntries();
1013 for (Int_t jui = 0; jui < (Int_t) vector2->GetEntriesFast(); jui++) {
1014 // Search if already found
1015 Int_t place = SearchInTreeVector(vectorplace,((AliTRDPlace *) vector2->At(jui))->GetPlace());
1016 // Create a new element in the two std vectors
1018 AliTRDPlace *placejjui = new AliTRDPlace();
1019 placejjui->SetPlace((j+jui));
1020 TObjArray *chainplace = new TObjArray();
1021 chainplace->Add((TObject *) placejjui);
1022 vectorplace->Add((TObject *) (vector2->At(jui)));
1023 where->Add((TObject *) chainplace);
1025 // Update the element at the place "place" in the std vector whereinthechain
1027 AliTRDPlace *placejjui = new AliTRDPlace();
1028 placejjui->SetPlace((j+jui));
1029 TObjArray *chainplace = ((TObjArray *) where->At(place));
1030 chainplace->Add((TObject *) placejjui);
1031 where->AddAt((TObject *) chainplace,place);
1036 treeChain->Add(filename2);
1039 // Take care of the profile
1040 const Char_t *pattern = "P";
1043 if (!strstr(variablecali,pattern)) {
1045 // Ready to read the chain
1047 treeChain->SetBranchAddress("histo",&his);
1049 // Initialise the final tree
1051 TH1F *histsum = 0x0;
1053 tree = new TTree(variablecali,variablecali);
1054 tree->Branch("groupnumber",&group,"groupnumber/I");
1055 tree->Branch("histo","TH1F",&histsum,32000,0);
1058 if (treeChain->GetEntries() < 1) {
1062 for (Int_t h = 0; h < (Int_t) vectorplace->GetEntriesFast(); h++) {
1063 group = ((AliTRDPlace *) vectorplace->At(h))->GetPlace();
1064 TObjArray *chainplace = ((TObjArray *) where->At(h));
1065 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(0))->GetPlace());
1066 //Init for the first time
1068 histsum = new TH1F("","",his->GetXaxis()->GetNbins()
1069 ,his->GetXaxis()->GetBinLowEdge(1)
1070 ,his->GetXaxis()->GetBinUpEdge(his->GetXaxis()->GetNbins()));
1073 // Reset for each new group
1074 histsum->SetEntries(0.0);
1075 for (Int_t l = 0; l <= histsum->GetXaxis()->GetNbins(); l++) {
1076 histsum->SetBinContent(l,0.0);
1077 histsum->SetBinError(l,0.0);
1079 histsum->Add(his,1);
1080 if ((Int_t) chainplace->GetEntriesFast() > 1) {
1081 for (Int_t s = 1; s < (Int_t) chainplace->GetEntriesFast(); s++) {
1082 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(s))->GetPlace());
1083 histsum->Add(his,1);
1092 // Ready to read the chain
1093 TGraphErrors *his = 0x0;
1094 treeChain->SetBranchAddress("histo",&his);
1096 // Initialise the final tree
1098 TGraphErrors *histsum = 0x0;
1099 Double_t *xref = 0x0;
1101 tree = new TTree(variablecali,variablecali);
1102 tree->Branch("groupnumber",&group,"groupnumber/I");
1103 tree->Branch("histo","TGraphErrors",&histsum,32000,0);
1106 if (treeChain->GetEntries() < 1) {
1110 for (Int_t h = 0; h < (Int_t) vectorplace->GetEntriesFast(); h++) {
1112 group = ((AliTRDPlace *) vectorplace->At(h))->GetPlace();
1113 TObjArray *chainplace = ((TObjArray *) where->At(h));
1114 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(0))->GetPlace());
1115 //Init or reset for a new group
1116 Int_t nbins = his->GetN();
1118 x = new Double_t[nbins];
1121 ex = new Double_t[nbins];
1123 y = new Double_t[nbins];
1125 ey = new Double_t[nbins];
1127 for (Int_t lo = 0; lo < nbins; lo++) {
1134 histsum = new TGraphErrors(nbins,x,y,ex,ey);
1137 histsum = AddProfiles(his,histsum);
1138 if ((Int_t) chainplace->GetEntriesFast() > 1) {
1139 for (Int_t s = 1; s < (Int_t) chainplace->GetEntriesFast(); s++) {
1140 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(s))->GetPlace());
1141 histsum = AddProfiles(his,histsum);
1155 //_____________________________________________________________________________
1156 TGraphErrors *AliTRDCalibraVector::AddProfiles(TGraphErrors *hist1
1157 , TGraphErrors *hist2) const
1160 // In the case of the vectors method we use TGraphErrors for PH and PRF
1161 // to be able to add the them after
1162 // Here we add the TGraphErrors
1165 // First TGraphErrors
1166 Int_t nbins1 = hist1->GetN();
1167 Double_t *x1 = hist1->GetX();
1168 Double_t *ex1 = hist1->GetEX();
1169 Double_t *y1 = hist1->GetY();
1170 Double_t *ey1 = hist1->GetEY();
1172 TGraphErrors *rehist = new TGraphErrors(nbins1);
1174 // Second TGraphErrors
1175 Double_t *ex2 = hist2->GetEX();
1176 Double_t *y2 = hist2->GetY();
1177 Double_t *ey2 = hist2->GetEY();
1179 // Define the Variables for the new TGraphErrors
1185 for (Int_t k = 0; k < nbins1; k++) {
1186 Double_t nentries = 0.0;
1191 if ((ex2[k] == 0.0) &&
1195 if ((ex2[k] == 0.0) &&
1202 if ((ex2[k] > 0.0) &&
1209 if ((ex2[k] > 0.0) &&
1211 nentries = ex1[k] + ex2[k];
1212 y = ( y1[k]*ex1[k]+ y2[k]*ex2[k]) / nentries;
1213 ey = (ey1[k]*ex1[k]+ey2[k]*ex2[k]) / nentries;
1216 rehist->SetPoint(k,x,y);
1217 rehist->SetPointError(k,ex,ey);