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 /////////////////////////////////////////////////////////////////////////////////
22 // This class is for the TRD calibration of the relative gain factor, the drift velocity,
23 // the time 0 and the pad response function.
24 // It can be used for the calibration per chamber but also per group of pads and eventually per pad.
25 // The user has to choose with the functions SetNz and SetNrphi the precision of the calibration.
31 <TR><TD><center>Nz</center></TD><TD><center> 0 </center></TD><TD><center> 1 </center></TD><TD><center> 2 </center></TD><TD><center> 3 </center></TD><TD><center> 4 </center></TD></TR>
32 <TR><TD><CENTER>group of row pads per detector</CENTER></TD><TD><CENTER>1</CENTER></TD><TD><CENTER>2</CENTER></TD><TD><CENTER>4</CENTER></TD><TD><CENTER>6(chamb2)<br> 8(others chambers)</CENTER></TD><TD><CENTER>12 (chamb2)<br> 16 (chamb0)</CENTER></TD></TR>
33 <TR><TD><CENTER>row pads per group</CENTER></TD><TD><CENTER>12 (chamb2)<br> 16 (chamb0)</CENTER></TD><TD><CENTER>6 (chamb2)<br> 8 (chamb0)</CENTER></TD><TD><CENTER>3 (chamb2)<br> 4 (chamb0)</CENTER></TD><TD><CENTER>2</CENTER></TD><TD><CENTER>1</CENTER></TD></TR>
34 <TR><TD><CENTER>~distance [cm]</CENTER></TD><TD><CENTER>106 (chamb2)<br> 130 (chamb0)</CENTER></TD><TD><CENTER>53 (chamb2)<br> 65 (chamb0)</CENTER></TD><TD><CENTER>26.5 (chamb2)<br> 32.5 (chamb0)</CENTER></TD><TD><CENTER>17 (chamb2)<br> 17 (chamb0)</CENTER></TD><TD><CENTER>9 (chamb2)<br> 9 (chamb0)</CENTER></TD></TR>
35 <CAPTION>In the z direction</CAPTION>
41 <TR><TD><center>Nrphi</center></TD><TD><center> 0 </center></TD><TD><center> 1 </center></TD><TD><center> 2 </center></TD><TD><center> 3 </center></TD><TD><center> 4 </center></TD><TD><center> 5 </center></TD><TD><center> 6 </center></TD></TR>
42 <TR><TD><CENTER>group of col pads per detector</CENTER></TD><TD><CENTER>1</CENTER></TD><TD><CENTER>2</CENTER></TD><TD><CENTER>4</CENTER></TD><TD><CENTER>8</CENTER></TD><TD><CENTER>16</CENTER></TD><TD><center>36</center></TD><TD><center>144</center></TD></TR>
43 <TR><TD><CENTER>col pads per group</CENTER></TD><TD><CENTER>144</CENTER></TD><TD><CENTER>72</CENTER></TD><TD><CENTER>36</CENTER></TD><TD><CENTER>18</CENTER></TD><TD><CENTER>9</CENTER></TD><TD><center>4</center></TD><TD><center>1</center></TD></TR>
44 <TR><TD><CENTER>~distance [cm]</CENTER></TD><TD><CENTER>113.4</CENTER></TD><TD><CENTER>56.7</CENTER></TD><TD><CENTER>25.3</CENTER></TD><TD><CENTER>14.3</CENTER></TD><TD><CENTER>7.25</CENTER></TD><TD><center>3.2</center></TD><TD><center>0.8</center></TD></TR>
45 <CAPTION>In the rphi direction</CAPTION>
52 // Fill histograms or vectors
53 //----------------------------
55 // 2D Histograms (Histo2d) or vectors (Vector2d), then converted in Trees, will be filled
56 // from RAW DATA in a run or from reconstructed TRD tracks during the offline tracking
57 // in the function "FollowBackProlongation" (AliTRDtracker)
58 // Per default the functions to fill are off.
62 Example of 2D histos for the relative gain (left) and the drift velocity (right) calibration of the sector 13 <br>
64 <img src="./gif/2dhisto.gif" width="600" height="350"><br>
69 // Fit the histograms to find the coefficients
70 //---------------------------------------------
72 // These 2D histograms or vectors (first converted in 1D histos) will be fitted
73 // if they have enough entries, otherwise the (default) value of the choosen database
74 // will be put. For the relative gain calibration the resulted factors will be globally
75 // normalized to the gain factors of the choosen database. It unables to precise
76 // previous calibration procedure.
77 // The function SetDebug enables the user to see:
78 // _fDebug = 0: nothing, only the values are written in the tree if wanted
79 // _fDebug = 1: a comparaison of the coefficients found and the default values
80 // in the choosen database.
81 // fCoef , histogram of the coefs as function of the calibration group number
82 // fDelta , histogram of the relative difference of the coef with the default
83 // value in the database as function of the calibration group number
84 // fError , dirstribution of this relative difference
85 // _fDebug = 2: only the fit of the choosen calibration group fFitVoir (SetFitVoir)
86 // _fDebug = 3: The coefficients in the choosen detector fDet (SetDet) as function of the
87 // pad row and col number
88 // _fDebug = 4; The coeffcicients in the choosen detector fDet (SetDet) like in the 3 but with
89 // also the comparaison histograms of the 1 for this detector
93 Example of fCoef for the relative gain calibration of the sector 13 <br>
95 <img src="./gif/coef.gif" width="400" height="460">
97 Example of fDelta (right) and fError (left) for the relative gain calibration of the sector 13 <br>
99 <img src="./gif/detlaerror.gif" width="550" height="380"><br>
105 // R. Bailhache (R.Bailhache@gsi.de)
107 //////////////////////////////////////////////////////////////////////////////////////
113 #include <TProfile2D.h>
116 #include <TGraphErrors.h>
117 #include <TObjArray.h>
125 #include <TStopwatch.h>
128 #include <TDirectory.h>
132 #include "AliCDBManager.h"
134 #include "AliRunLoader.h"
135 #include "AliLoader.h"
136 #include "AliRawReaderFile.h"
137 #include "AliRawReader.h"
139 #include "AliTRDCalibra.h"
140 #include "AliTRDcalibDB.h"
141 #include "AliTRDCommonParam.h"
142 #include "AliTRDmcmTracklet.h"
143 #include "AliTRDpadPlane.h"
144 #include "AliTRDcluster.h"
145 #include "AliTRDtrack.h"
146 #include "AliTRDdigit.h"
147 #include "AliTRDdigitsManager.h"
149 #include "AliTRDgeometry.h"
150 #include "./Cal/AliTRDCalROC.h"
151 #include "./Cal/AliTRDCalPad.h"
152 #include "./Cal/AliTRDCalDet.h"
153 #include "AliTRDrawData.h"
155 ClassImp(AliTRDCalibra)
157 AliTRDCalibra* AliTRDCalibra::fgInstance = 0;
158 Bool_t AliTRDCalibra::fgTerminated = kFALSE;
160 //_____________singleton implementation_________________________________________________
161 AliTRDCalibra *AliTRDCalibra::Instance()
164 // Singleton implementation
167 if (fgTerminated != kFALSE) {
171 if (fgInstance == 0) {
172 fgInstance = new AliTRDCalibra();
179 //______________________________________________________________________________________
180 void AliTRDCalibra::Terminate()
183 // Singleton implementation
184 // Deletes the instance of this class
187 fgTerminated = kTRUE;
189 if (fgInstance != 0) {
196 //______________________________________________________________________________________
197 AliTRDCalibra::AliTRDCalibra()
200 ,fMcmTracking(kFALSE)
201 ,fMcmCorrectAngle(kFALSE)
208 ,fCountRelativeScale(0)
209 ,fRelativeScaleAuto(kFALSE)
211 ,fThresholdClusterPRF1(0.0)
212 ,fThresholdClusterPRF2(0.0)
213 ,fCenterOfflineCluster(kFALSE)
219 ,fBeginFitCharge(0.0)
221 ,fMeanChargeOn(kFALSE)
222 ,fFitChargeBisOn(kFALSE)
241 ,fDetectorAliTRDtrack(kFALSE)
242 ,fChamberAliTRDtrack(-1)
243 ,fDetectorPreviousTrack(-1)
257 ,fScaleFitFactor(0.0)
278 // Default constructor
281 for (Int_t i = 0; i < 3; i++) {
286 for (Int_t k = 0; k < 3; k++) {
293 for (Int_t i = 0; i < 3; i++) {
294 fWriteCoef[i] = kFALSE;
299 for (Int_t k = 0; k < 3; k++) {
303 for (Int_t i = 0; i < 3; i++) {
312 //______________________________________________________________________________________
313 AliTRDCalibra::AliTRDCalibra(const AliTRDCalibra &c)
316 ,fMcmTracking(kFALSE)
317 ,fMcmCorrectAngle(kFALSE)
324 ,fCountRelativeScale(0)
325 ,fRelativeScaleAuto(kFALSE)
327 ,fThresholdClusterPRF1(0.0)
328 ,fThresholdClusterPRF2(0.0)
329 ,fCenterOfflineCluster(kFALSE)
335 ,fBeginFitCharge(0.0)
337 ,fMeanChargeOn(kFALSE)
338 ,fFitChargeBisOn(kFALSE)
357 ,fDetectorAliTRDtrack(kFALSE)
358 ,fChamberAliTRDtrack(-1)
359 ,fDetectorPreviousTrack(-1)
373 ,fScaleFitFactor(0.0)
399 //____________________________________________________________________________________
400 AliTRDCalibra::~AliTRDCalibra()
403 // AliTRDCalibra destructor
411 //_____________________________________________________________________________
412 void AliTRDCalibra::Destroy()
425 //_____________________________________________________________________________
426 void AliTRDCalibra::ClearHistos()
447 //_____________________________________________________________________________
448 void AliTRDCalibra::ClearTree()
473 //_____________________________________________________________________________
474 void AliTRDCalibra::Init()
477 // Init some default values
480 // How to fill the 2D
482 fThresholdClusterPRF1 = 2.0;
483 fThresholdClusterPRF2 = 20.0;
486 fNumberBinCharge = 100;
490 fWriteName = "TRD.calibration.root";
491 fWriteNameCoef = "TRD.coefficient.root";
495 fBeginFitCharge = 3.5;
500 // Internal variables
502 // Fill the 2D histos in the offline tracking
503 fDetectorPreviousTrack = -1;
504 fChamberAliTRDtrack = -1;
509 fNumberClusters = 18;
511 fNumberUsedCh[0] = 0;
512 fNumberUsedCh[1] = 0;
513 fNumberUsedPh[0] = 0;
514 fNumberUsedPh[1] = 0;
516 // Variables in the loop
517 for (Int_t k = 0; k < 4; k++) {
518 fChargeCoef[k] = 1.0;
519 fVdriftCoef[k] = 1.5;
522 for (Int_t i = 0; i < 2; i++) {
527 for (Int_t i = 0; i < 3; i++) {
539 // Local database to be changed
544 //____________Functions fit Online CH2d________________________________________
545 Bool_t AliTRDCalibra::FitCHOnline(TH2I *ch)
548 // Fit the 1D histos, projections of the 2D ch on the Xaxis, for each
549 // calibration group normalized the resulted coefficients (to 1 normally)
550 // and write the results in a tree
553 // Number of Xbins (detectors or groups of pads)
554 TAxis *xch = ch->GetXaxis();
555 Int_t nbins = xch->GetNbins();
556 TAxis *yph = ch->GetYaxis();
557 Int_t nybins = yph->GetNbins();
558 Double_t lowedge = xch->GetBinLowEdge(1);
559 Double_t upedge = xch->GetBinUpEdge(xch->GetNbins());
560 if (!InitFit(nbins,lowedge,upedge,0)) {
563 fStatisticMean = 0.0;
575 // Init fCountDet and fCount
576 InitfCountDetAndfCount(0);
578 // Beginning of the loop betwwen dect1 and dect2
579 for (Int_t idect = fDect1[0]; idect < fDect2[0]; idect++) {
581 TH1I *projch = (TH1I *) ch->ProjectionY("projch",idect+1,idect+1,(Option_t *)"e");
582 projch->SetDirectory(0);
584 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
585 UpdatefCountDetAndfCount(idect,0);
587 // Reconstruction of the row and pad group: rowmin, row max ...
588 ReconstructFitRowMinRowMax(idect, 0);
590 // Number of entries for this calibration group
591 Double_t nentries = 0.0;
592 for (Int_t k = 0; k < nybins; k++) {
593 nentries += ch->GetBinContent(ch->GetBin(idect+1,k+1));
599 // Rebin and statistic stuff
602 projch = ReBin((TH1I *) projch);
604 // This detector has not enough statistics or was off
605 if (nentries < fMinEntries) {
606 // Fill with the default infos
607 NotEnoughStatistic(idect,0);
615 // Statistics of the group fitted
616 AliInfo(Form("For the group number %d there are %f stats",idect,nentries));
617 fStatisticMean += nentries;
620 // Method Mean and fit
621 // idect is egal for fDebug = 0 and 2, only to fill the hist
622 FitCH((TH1 *) projch,(Int_t) (idect-fDect1[0]));
624 // idect is egal for fDebug = 0 and 2, only to fill the hist
625 if (fFitChargeBisOn) {
626 FitBisCH((TH1 *) projch,(Int_t) (idect-fDect1[0]));
629 // Visualise the detector for fDebug 3 or 4
630 // Here is the reconstruction of the pad and row group is used!
635 FillInfosFit(idect,0);
656 // 0 no plot, 1 and 4 error plot, 3 and 4 DB plot
667 if (fNumberFit > 0) {
668 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean/fNumberFit));
669 fStatisticMean = fStatisticMean / fNumberFit;
672 AliInfo("There is no fit!");
676 ConvertVectorFitCHTree();
685 //____________Functions fit Online CH2d________________________________________
686 Bool_t AliTRDCalibra::FitCHOnline()
689 // Reconstruct a 1D histo from the vectorCH for each calibration group,
690 // fit the histo, normalized the resulted coefficients (to 1 normally)
691 // and write the results in a tree
694 // Number of Xbins (detectors or groups of pads)
695 if (!InitFit(0,0,0,0)) {
698 fStatisticMean = 0.0;
702 // Init fCountDet and fCount
703 InitfCountDetAndfCount(0);
705 // Beginning of the loop between dect1 and dect2
706 for (Int_t idect = fDect1[0]; idect < fDect2[0]; idect++) {
708 // Search if the group is in the VectorCH
709 Int_t place = SearchInVector(idect,0);
717 AliTRDCTInfo *fCHInfo = new AliTRDCTInfo();
719 fCHInfo = ((AliTRDCTInfo *) fVectorCH->At(place));
720 projch = ConvertVectorCTHisto(fCHInfo,(const char *) name);
721 projch->SetDirectory(0);
725 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
726 UpdatefCountDetAndfCount(idect,0);
728 // Reconstruction of the row and pad group: rowmin, row max ...
729 ReconstructFitRowMinRowMax(idect,0);
732 Double_t nentries = 0.0;
734 for (Int_t k = 0; k < fNumberBinCharge; k++) {
735 nentries += projch->GetBinContent(k+1);
742 // Rebin and statistic stuff
746 projch = ReBin((TH1F *) projch);
749 // This detector has not enough statistics or was not found in VectorCH
752 (nentries < fMinEntries))) {
754 // Fill with the default infos
755 NotEnoughStatistic(idect,0);
766 // Statistic of the histos fitted
767 AliInfo(Form("For the group number %d there are %f stats",idect,nentries));
769 fStatisticMean += nentries;
771 // Method Mean and fit
772 // idect is egal for fDebug = 0 and 2, only to fill the hist
773 FitCH((TH1 *) projch,(Int_t) (idect-fDect1[0]));
776 // idect is egal for fDebug = 0 and 2, only to fill the hist
777 if (fFitChargeBisOn) {
778 FitBisCH((TH1 *) projch,(Int_t) (idect-fDect1[0]));
781 // Visualise the detector for fDebug 3 or 4
782 // Here is the reconstruction of the pad and row group is used!
788 FillInfosFit(idect,0);
809 // 0 no plot, 1 and 4 error plot, 3 and 4 DB plot
820 if (fNumberFit > 0) {
821 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean / fNumberFit));
822 fStatisticMean = fStatisticMean / fNumberFit;
825 AliInfo("There is no fit!");
829 ConvertVectorFitCHTree();
838 //____________Functions fit Online CH2d________________________________________
839 Bool_t AliTRDCalibra::FitCHOnline(TTree *tree)
842 // Look if the calibration group can be found in the tree, if yes take the
843 // histo, fit it, normalized the resulted coefficients (to 1 normally) and
844 // write the results in a tree
847 // Number of Xbins (detectors or groups of pads)
848 if (!InitFit(0,0,0,0)) {
851 fStatisticMean = 0.0;
863 // Init fCountDet and fCount
864 InitfCountDetAndfCount(0);
866 tree->SetBranchAddress("histo",&projch);
867 TObjArray *vectorplace = ConvertTreeVector(tree);
869 // Beginning of the loop between dect1 and dect2
870 for (Int_t idect = fDect1[0]; idect < fDect2[0]; idect++) {
872 //Search if the group is in the VectorCH
873 Int_t place = SearchInTreeVector(vectorplace,idect);
878 tree->GetEntry(place);
881 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
882 UpdatefCountDetAndfCount(idect,0);
884 // Reconstruction of the row and pad group: rowmin, row max ...
885 ReconstructFitRowMinRowMax(idect,0);
888 Double_t nentries = 0.0;
890 for (Int_t k = 0; k < projch->GetXaxis()->GetNbins(); k++) {
891 nentries += projch->GetBinContent(k+1);
898 // Rebin and statistic stuff
902 projch = ReBin((TH1F *) projch);
905 // This detector has not enough statistics or was not found in VectorCH
908 (nentries < fMinEntries))) {
910 // Fill with the default infos
911 NotEnoughStatistic(idect,0);
917 // Statistics of the group fitted
918 AliInfo(Form("For the group number %d there are %f stats",idect,nentries));
920 fStatisticMean += nentries;
922 // Method Mean and fit
923 // idect is egal for fDebug = 0 and 2, only to fill the hist
924 FitCH((TH1 *) projch,(Int_t) (idect-fDect1[0]));
927 // idect is egal for fDebug = 0 and 2, only to fill the hist
928 if (fFitChargeBisOn) {
929 FitBisCH((TH1 *) projch,(Int_t) (idect-fDect1[0]));
932 // Visualise the detector for fDebug 3 or 4
933 // Here is the reconstruction of the pad and row group is used!
939 FillInfosFit(idect,0);
955 // 0 no plot, 1 and 4 error plot, 3 and 4 DB plot
966 if (fNumberFit > 0) {
967 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean / fNumberFit));
968 fStatisticMean = fStatisticMean / fNumberFit;
971 AliInfo("There is no fit!");
975 ConvertVectorFitCHTree();
984 //________________functions fit Online PH2d____________________________________
985 Bool_t AliTRDCalibra::FitPHOnline(TProfile2D *ph)
988 // Take the 1D profiles (average pulse height), projections of the 2D PH
989 // on the Xaxis, for each calibration group
990 // Fit or use the slope of the average pulse height to reconstruct the
991 // drift velocity write the results in a tree
992 // A first calibration of T0 is also made using the same method (slope method)
995 // Number of Xbins (detectors or groups of pads)
996 TAxis *xph = ph->GetXaxis();
997 TAxis *yph = ph->GetYaxis();
998 Int_t nbins = xph->GetNbins();
999 Int_t nybins = yph->GetNbins();
1000 Double_t lowedge = xph->GetBinLowEdge(1);
1001 Double_t upedge = xph->GetBinUpEdge(xph->GetNbins());
1002 if (!InitFit(nbins,lowedge,upedge,1)) {
1005 fStatisticMean = 0.0;
1017 // Init fCountDet and fCount
1018 InitfCountDetAndfCount(1);
1020 // Beginning of the loop
1021 for (Int_t idect = fDect1[1]; idect < fDect2[1]; idect++) {
1023 TH1D *projph = (TH1D *) ph->ProjectionY("projph",idect+1,idect+1,(Option_t *) "e");
1024 projph->SetDirectory(0);
1026 // Number of entries for this calibration group
1027 Double_t nentries = 0;
1028 for (Int_t k = 0; k < nybins; k++) {
1029 nentries += ph->GetBinEntries(ph->GetBin(idect+1,k+1));
1035 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
1036 UpdatefCountDetAndfCount(idect,1);
1038 // Reconstruction of the row and pad group: rowmin, row max ...
1039 ReconstructFitRowMinRowMax(idect,1);
1041 // Rebin and statistic stuff
1042 // This detector has not enough statistics or was off
1043 if (nentries < fMinEntries) {
1045 // Fill with the default values
1046 NotEnoughStatistic(idect,1);
1057 // Statistics of the histos fitted
1058 AliInfo(Form("For the group number %d there are %f stats",idect,nentries));
1060 fStatisticMean += nentries;
1062 // Calcul of "real" coef
1063 CalculVdriftCoefMean(fCountDet[1],(Int_t) (idect - fDect1[1]));
1064 CalculT0CoefMean(fCountDet[1],(Int_t) (idect - fDect1[1]));
1066 // Method Mean and fit
1067 // idect is egal for fDebug = 0 and 2, only to fill the hist
1068 FitPente((TH1 *) projph,(Int_t) (idect - fDect1[1]));
1071 // idect is egal for fDebug = 0 and 2, only to fill the hist
1073 FitPH((TH1 *) projph,(Int_t) (idect - fDect1[1]));
1076 // Visualise the detector for fDebug 3 or 4
1077 // Here is the reconstruction of the pad and row group is used!
1083 // Fill the tree if end of a detector or only the pointer to the branch!!!
1084 FillInfosFit(idect,1);
1094 if ((fDebug == 1) ||
1101 // 0 no plot, 1 and 4 error plot, 3 and 4 DB plot
1102 if ((fDebug == 1) ||
1107 if ((fDebug == 4) ||
1114 if (fNumberFit > 0) {
1115 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean / fNumberFit));
1116 fStatisticMean = fStatisticMean / fNumberFit;
1119 AliInfo("There is no fit!");
1122 // Write the things!
1131 //____________Functions fit Online PH2d________________________________________
1132 Bool_t AliTRDCalibra::FitPHOnline()
1135 // Reconstruct the average pulse height from the vectorPH for each
1136 // calibration group
1137 // Fit or use the slope of the average pulse height to reconstruct the
1138 // drift velocity write the results in a tree
1139 // A first calibration of T0 is also made using the same method (slope method)
1142 // Number of Xbins (detectors or groups of pads)
1143 if (!InitFit(0,0,0,1)) {
1146 fStatisticMean = 0.0;
1150 // Init fCountDet and fCount
1151 InitfCountDetAndfCount(1);
1153 // Beginning of the loop
1154 for (Int_t idect = fDect1[1]; idect < fDect2[1]; idect++) {
1156 // Search if the group is in the VectorCH
1157 Int_t place = SearchInVector(idect,1);
1167 AliTRDPInfo *fPHInfo = new AliTRDPInfo();
1169 fPHInfo = (AliTRDPInfo *) fVectorPH->At(place);
1170 projph = CorrectTheError((TGraphErrors *) ConvertVectorPHisto(fPHInfo,(const char *) name));
1171 projph->SetDirectory(0);
1175 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
1176 UpdatefCountDetAndfCount(idect,1);
1178 // Reconstruction of the row and pad group: rowmin, row max ...
1179 ReconstructFitRowMinRowMax(idect,1);
1181 // Rebin and statistic stuff
1182 // This detector has not enough statistics or was off
1183 if ((place == -1) ||
1185 (fEntriesCurrent < fMinEntries))) {
1187 // Fill with the default values
1188 NotEnoughStatistic(idect,1);
1199 // Statistic of the histos fitted
1200 AliInfo(Form("For the group number %d there are %d stats",idect,fEntriesCurrent));
1202 fStatisticMean += fEntriesCurrent;
1204 // Calcul of "real" coef
1205 CalculVdriftCoefMean(fCountDet[1],(Int_t) (idect - fDect1[1]));
1206 CalculT0CoefMean(fCountDet[1],(Int_t) (idect - fDect1[1]));
1208 // Method Mean and fit
1209 // idect is egal for fDebug = 0 and 2, only to fill the hist
1210 FitPente((TH1 *) projph,(Int_t) (idect - fDect1[1]));
1213 // idect is egal for fDebug = 0 and 2, only to fill the hist
1215 FitPH((TH1 *) projph,(Int_t) (idect - fDect1[1]));
1218 // Visualise the detector for fDebug 3 or 4
1219 // Here is the reconstruction of the pad and row group is used!
1225 // Fill the tree if end of a detector or only the pointer to the branch!!!
1226 FillInfosFit(idect,1);
1236 if ((fDebug == 1) ||
1243 // 0 no plot, 1 and 4 error plot, 3 and 4 DB plot
1244 if ((fDebug == 1) ||
1249 if ((fDebug == 4) ||
1256 if (fNumberFit > 0) {
1257 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean / fNumberFit));
1258 fStatisticMean = fStatisticMean / fNumberFit;
1261 AliInfo("There is no fit!");
1264 // Write the things!
1265 if (fWriteCoef[1]) {
1273 //____________Functions fit Online PH2d________________________________________
1274 Bool_t AliTRDCalibra::FitPHOnline(TTree *tree)
1277 // Look if the calibration group can be found in the tree, if yes take the
1278 // histo, fit it, and write the results in a tree
1279 // A first calibration of T0 is also made using the same method (slope method)
1282 // Number of Xbins (detectors or groups of pads)
1283 if (!InitFit(0,0,0,1)) {
1286 fStatisticMean = 0.0;
1298 // Init fCountDet and fCount
1299 InitfCountDetAndfCount(1);
1300 TGraphErrors *projphtree = 0x0;
1301 tree->SetBranchAddress("histo",&projphtree);
1302 TObjArray *vectorplace = ConvertTreeVector(tree);
1304 // Beginning of the loop
1305 for (Int_t idect = fDect1[1]; idect < fDect2[1]; idect++) {
1307 // Search if the group is in the VectorCH
1308 Int_t place = SearchInTreeVector(vectorplace,idect);
1316 tree->GetEntry(place);
1317 projph = CorrectTheError(projphtree);
1320 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
1321 UpdatefCountDetAndfCount(idect,1);
1323 // Reconstruction of the row and pad group: rowmin, row max ...
1324 ReconstructFitRowMinRowMax(idect,1);
1326 // Rebin and statistic stuff
1327 // This detector has not enough statistics or was off
1330 (fEntriesCurrent < fMinEntries))) {
1332 // Fill with the default values
1333 NotEnoughStatistic(idect,1);
1344 // Statistics of the histos fitted
1345 AliInfo(Form("For the group number %d there are %d stats",idect,fEntriesCurrent));
1347 fStatisticMean += fEntriesCurrent;
1349 // Calcul of "real" coef
1350 CalculVdriftCoefMean(fCountDet[1],(Int_t) (idect - fDect1[1]));
1351 CalculT0CoefMean(fCountDet[1],(Int_t) (idect - fDect1[1]));
1353 // Method Mean and fit
1354 // idect is egal for fDebug = 0 and 2, only to fill the hist
1355 FitPente((TH1 *) projph,(Int_t) (idect - fDect1[1]));
1357 // idect is egal for fDebug = 0 and 2, only to fill the hist
1359 FitPH((TH1 *) projph,(Int_t) (idect - fDect1[1]));
1362 // Visualise the detector for fDebug 3 or 4
1363 // Here is the reconstruction of the pad and row group is used!
1369 // Fill the tree if end of a detector or only the pointer to the branch!!!
1370 FillInfosFit(idect,1);
1380 if ((fDebug == 1) ||
1387 // 0 no plot, 1 and 4 error plot, 3 and 4 DB plot
1388 if ((fDebug == 1) ||
1393 if ((fDebug == 4) ||
1400 if (fNumberFit > 0) {
1401 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean / fNumberFit));
1402 fStatisticMean = fStatisticMean / fNumberFit;
1405 AliInfo("There is no fit!");
1408 // Write the things!
1409 if (fWriteCoef[1]) {
1417 //____________Functions fit Online PRF2d_______________________________________
1418 Bool_t AliTRDCalibra::FitPRFOnline(TProfile2D *prf)
1421 // Take the 1D profiles (pad response function), projections of the 2D PRF
1422 // on the Xaxis, for each calibration group
1423 // Fit with a gaussian to reconstruct the sigma of the pad response function
1424 // write the results in a tree
1427 // Number of Xbins (detectors or groups of pads)
1428 TAxis *xprf = prf->GetXaxis();
1429 TAxis *yprf = prf->GetYaxis();
1430 Int_t nybins = yprf->GetNbins();
1431 Int_t nbins = xprf->GetNbins();
1432 Double_t lowedge = xprf->GetBinLowEdge(1);
1433 Double_t upedge = xprf->GetBinUpEdge(xprf->GetNbins());
1434 if (!InitFit(nbins,lowedge,upedge,2)) {
1437 fStatisticMean = 0.0;
1443 fVectorPRF->Clear();
1449 // Init fCountDet and fCount
1450 InitfCountDetAndfCount(2);
1452 // Beginning of the loop
1453 for (Int_t idect = fDect1[2]; idect < fDect2[2]; idect++) {
1455 TH1D *projprf = (TH1D *) prf->ProjectionY("projprf",idect+1,idect+1,(Option_t *) "e");
1456 projprf->SetDirectory(0);
1458 // Number of entries for this calibration group
1459 Double_t nentries = 0;
1460 for (Int_t k = 0; k < nybins; k++) {
1461 nentries += prf->GetBinEntries(prf->GetBin(idect+1,k+1));
1463 if(nentries > 0) fNumberEnt++;
1465 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
1466 UpdatefCountDetAndfCount(idect,2);
1468 // Reconstruction of the row and pad group: rowmin, row max ...
1469 ReconstructFitRowMinRowMax(idect,2);
1471 // Rebin and statistic stuff
1472 // This detector has not enough statistics or was off
1473 if (nentries < fMinEntries) {
1475 // Fill with the default values
1476 NotEnoughStatistic(idect,2);
1487 // Statistics of the histos fitted
1488 AliInfo(Form("For the group number %d there are %f stats",idect,nentries));
1490 fStatisticMean += nentries;
1492 // Calcul of "real" coef
1493 if ((fDebug == 1) ||
1495 CalculPRFCoefMean(fCountDet[2],(Int_t) (idect - fDect1[2]));
1498 // Method Mean and fit
1499 // idect is egal for fDebug = 0 and 2, only to fill the hist
1500 FitPRF((TH1 *) projprf,(Int_t) (idect - fDect1[2]));
1502 // Visualise the detector for fDebug 3 or 4
1503 // Here is the reconstruction of the pad and row group is used!
1508 // Fill the tree if end of a detector or only the pointer to the branch!!!
1509 FillInfosFit(idect,2);
1519 if ((fDebug == 1) ||
1525 // No plot, 1 and 4 error plot, 3 and 4 DB plot
1526 if ((fDebug == 1) ||
1530 if ((fDebug == 4) ||
1536 if (fNumberFit > 0) {
1537 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean / fNumberFit));
1538 fStatisticMean = fStatisticMean / fNumberFit;
1541 AliInfo("There is no fit!");
1544 // Write the things!
1545 if (fWriteCoef[2]) {
1553 //____________Functions fit Online PRF2d_______________________________________
1554 Bool_t AliTRDCalibra::FitPRFOnline(TTree *tree)
1557 // Look if the calibration group can be found in the tree, if yes take
1558 // the histo, fit it, and write the results in a tree
1561 // Number of Xbins (detectors or groups of pads)
1562 if (!InitFit(0,0,0,2)) {
1565 fStatisticMean = 0.0;
1571 fVectorPRF->Clear();
1577 // Init fCountDet and fCount
1578 InitfCountDetAndfCount(2);
1579 TGraphErrors *projprftree = 0x0;
1580 tree->SetBranchAddress("histo",&projprftree);
1581 TObjArray *vectorplace = ConvertTreeVector(tree);
1583 // Beginning of the loop
1584 for (Int_t idect = fDect1[2]; idect < fDect2[2]; idect++) {
1586 // Search if the group is in the VectorCH
1587 Int_t place = SearchInTreeVector(vectorplace,idect);
1590 TH1F *projprf = 0x0;
1595 tree->GetEntry(place);
1596 projprf = CorrectTheError(projprftree);
1599 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
1600 UpdatefCountDetAndfCount(idect,2);
1602 // Reconstruction of the row and pad group: rowmin, row max ...
1603 ReconstructFitRowMinRowMax(idect,2);
1605 // Rebin and statistic stuff
1606 // This detector has not enough statistics or was off
1607 if ((place == -1) ||
1609 (fEntriesCurrent < fMinEntries))) {
1611 // Fill with the default values
1612 NotEnoughStatistic(idect,2);
1623 // Statistics of the histos fitted
1624 AliInfo(Form("For the group number %d there are %d stats",idect,fEntriesCurrent));
1626 fStatisticMean += fEntriesCurrent;
1628 // Calcul of "real" coef
1629 if ((fDebug == 1) ||
1631 CalculPRFCoefMean(fCountDet[2],(Int_t) (idect - fDect1[2]));
1634 // Method Mean and fit
1635 // idect is egal for fDebug = 0 and 2, only to fill the hist
1636 FitPRF((TH1 *) projprf,(Int_t) (idect - fDect1[2]));
1637 // Visualise the detector for fDebug 3 or 4
1638 // Here is the reconstruction of the pad and row group is used!
1642 // Fill the tree if end of a detector or only the pointer to the branch!!!
1643 FillInfosFit(idect,2);
1653 if ((fDebug == 1) ||
1659 // No plot, 1 and 4 error plot, 3 and 4 DB plot
1660 if ((fDebug == 1) ||
1664 if ((fDebug == 4) ||
1670 if (fNumberFit > 0) {
1671 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean / fNumberFit));
1672 fStatisticMean = fStatisticMean / fNumberFit;
1675 AliInfo("There is no fit!");
1678 // Write the things!
1679 if (fWriteCoef[2]) {
1687 //____________Functions fit Online PRF2d_______________________________________
1688 Bool_t AliTRDCalibra::FitPRFOnline()
1691 // Reconstruct the 1D histo (pad response function) from the vectorPRD for
1692 // each calibration group
1693 // Fit with a gaussian to reconstruct the sigma of the pad response function
1694 // write the results in a tree
1697 // Number of Xbins (detectors or groups of pads)
1698 if (!InitFit(0,0,0,2)) {
1701 fStatisticMean = 0.0;
1705 // Init fCountDet and fCount
1706 InitfCountDetAndfCount(2);
1708 // Beginning of the loop
1709 for (Int_t idect = fDect1[2]; idect < fDect2[2]; idect++) {
1711 // Search if the group is in the VectorCH
1712 Int_t place = SearchInVector(idect,2);
1715 TH1F *projprf = 0x0;
1716 TString name("PRF");
1722 AliTRDPInfo *fPRFInfo = new AliTRDPInfo();
1724 fPRFInfo = (AliTRDPInfo *) fVectorPRF->At(place);
1725 projprf = CorrectTheError((TGraphErrors *) ConvertVectorPHisto(fPRFInfo,(const char *)name));
1726 projprf->SetDirectory(0);
1730 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
1731 UpdatefCountDetAndfCount(idect,2);
1733 // Reconstruction of the row and pad group: rowmin, row max ...
1734 ReconstructFitRowMinRowMax(idect,2);
1736 // Rebin and statistic stuff
1737 // This detector has not enough statistics or was off
1738 if ((place == -1) ||
1740 (fEntriesCurrent < fMinEntries))) {
1742 // Fill with the default values
1743 NotEnoughStatistic(idect,2);
1754 // Statistic of the histos fitted
1755 AliInfo(Form("For the group number %d there are %d stats", idect,fEntriesCurrent));
1757 fStatisticMean += fEntriesCurrent;
1759 // Calcul of "real" coef
1760 if ((fDebug == 1) ||
1762 CalculPRFCoefMean(fCountDet[2],(Int_t) (idect-fDect1[2]));
1765 // Method Mean and fit
1766 // idect is egal for fDebug = 0 and 2, only to fill the hist
1767 FitPRF((TH1 *) projprf,(Int_t) (idect-fDect1[2]));
1768 // Visualise the detector for fDebug 3 or 4
1769 // Here is the reconstruction of the pad and row group is used!
1773 // Fill the tree if end of a detector or only the pointer to the branch!!!
1774 FillInfosFit(idect,2);
1784 if ((fDebug == 1) ||
1790 // No plot, 1 and 4 error plot, 3 and 4 DB plot
1791 if ((fDebug == 1) ||
1795 if ((fDebug == 4) ||
1801 if (fNumberFit > 0) {
1802 AliInfo(Form("There is a mean statistic of: %d",(Int_t) fStatisticMean / fNumberFit));
1805 AliInfo("There is no fit!");
1808 // Write the things!
1809 if (fWriteCoef[2]) {
1817 //____________Functions for initialising the AliTRDCalibra in the code_________
1818 Bool_t AliTRDCalibra::Init2Dhistos()
1821 // For the offline tracking
1822 // This function will be called in the function AliReconstruction::Run()
1823 // Init the calibration mode (Nz, Nrphi), the 2D histograms if fHisto2d = kTRUE,
1828 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
1830 AliInfo("Could not get calibDB");
1835 fTimeMax = cal->GetNumberOfTimeBins();
1836 fSf = cal->GetSamplingFrequency();
1837 if (fRelativeScaleAuto) {
1841 fRelativeScale = 20;
1844 // Create the 2D histos corresponding to the pad groupCalibration mode
1847 AliInfo(Form("The pad calibration mode for the relative gain calibration: Nz %d, and Nrphi %d"
1851 // Calcul the number of Xbins
1853 ModePadCalibration(2,0);
1854 ModePadFragmentation(0,2,0,0);
1855 fDetChamb2[0] = fNfragZ[0] * fNfragRphi[0];
1857 AliInfo(Form("For the chamber 2: %d",fDetChamb2[0]));
1859 fNtotal[0] += 6 * 18 * fDetChamb2[0];
1860 ModePadCalibration(0,0);
1861 ModePadFragmentation(0,0,0,0);
1862 fDetChamb0[0] = fNfragZ[0] * fNfragRphi[0];
1864 AliInfo(Form("For the other chamber 0: %d",fDetChamb0[0]));
1866 fNtotal[0] += 6 * 4 * 18 * fDetChamb0[0];
1867 AliInfo(Form("Total number of Xbins: %d",fNtotal[0]));
1869 // Create the 2D histo
1871 CreateCH2d(fNtotal[0]);
1874 fVectorCH = new TObjArray();
1875 fPlaCH = new TObjArray();
1879 fAmpTotal = new Float_t[TMath::Max(fDetChamb2[0],fDetChamb0[0])];
1880 for (Int_t k = 0; k < TMath::Max(fDetChamb2[0],fDetChamb0[0]); k++) {
1888 AliInfo(Form("The pad calibration mode for the drift velocity calibration: Nz %d, and Nrphi %d"
1892 // Calcul the number of Xbins
1894 ModePadCalibration(2,1);
1895 ModePadFragmentation(0,2,0,1);
1896 fDetChamb2[1] = fNfragZ[1]*fNfragRphi[1];
1898 AliInfo(Form("For the chamber 2: %d",fDetChamb2[1]));
1900 fNtotal[1] += 6 * 18 * fDetChamb2[1];
1901 ModePadCalibration(0,1);
1902 ModePadFragmentation(0,0,0,1);
1903 fDetChamb0[1] = fNfragZ[1] * fNfragRphi[1];
1905 AliInfo(Form("For the chamber 0: %d",fDetChamb0[1]));
1907 fNtotal[1] += 6 * 4 * 18 * fDetChamb0[1];
1908 AliInfo(Form("Total number of Xbins: %d",fNtotal[1]));
1910 // Create the 2D histo
1912 CreatePH2d(fNtotal[1]);
1915 fVectorPH = new TObjArray();
1916 fPlaPH = new TObjArray();
1920 fPHPlace = new Short_t[fTimeMax];
1921 for (Int_t k = 0; k < fTimeMax; k++) {
1924 fPHValue = new Float_t[fTimeMax];
1925 for (Int_t k = 0; k < fTimeMax; k++) {
1933 AliInfo(Form("The pad calibration mode for the PRF calibration: Nz %d, and Nrphi %d"
1937 // Calcul the number of Xbins
1939 ModePadCalibration(2,2);
1940 ModePadFragmentation(0,2,0,2);
1941 fDetChamb2[2] = fNfragZ[2] * fNfragRphi[2];
1943 AliInfo(Form("For the chamber 2: %d",fDetChamb2[2]));
1945 fNtotal[2] += 6 * 18 * fDetChamb2[2];
1946 ModePadCalibration(0,2);
1947 ModePadFragmentation(0,0,0,2);
1948 fDetChamb0[2] = fNfragZ[2] * fNfragRphi[2];
1950 AliInfo(Form("For the chamber 0: %d",fDetChamb0[2]));
1952 fNtotal[2] += 6 * 4 * 18 * fDetChamb0[2];
1953 AliInfo(Form("Total number of Xbins: %d",fNtotal[2]));
1955 // Create the 2D histo
1957 CreatePRF2d(fNtotal[2]);
1960 fVectorPRF = new TObjArray();
1961 fPlaPRF = new TObjArray();
1970 //____________Functions for filling the histos in the code_____________________
1972 //____________Offine tracking in the AliTRDtracker_____________________________
1973 Bool_t AliTRDCalibra::ResetTrack()
1976 // For the offline tracking
1977 // This function will be called in the function
1978 // AliTRDtracker::FollowBackPropagation() at the beginning
1979 // Reset the parameter to know we have a new TRD track
1982 fDetectorAliTRDtrack = kFALSE;
1987 //____________Offline tracking in the AliTRDtracker____________________________
1988 Bool_t AliTRDCalibra::UpdateHistograms(AliTRDcluster *cl, AliTRDtrack *t)
1991 // For the offline tracking
1992 // This function will be called in the function
1993 // AliTRDtracker::FollowBackPropagation() in the loop over the clusters
1995 // Fill the 2D histos or the vectors with the info of the clusters at
1996 // the end of a detectors if the track is "good"
1999 // Get the parameter object
2000 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
2002 AliInfo("Could not get CommonParam");
2006 // Get the parameter object
2007 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
2009 AliInfo("Could not get calibDB");
2013 // Localisation of the detector
2014 Int_t detector = cl->GetDetector();
2015 Int_t chamber = GetChamber(detector);
2016 Int_t plane = GetPlane(detector);
2018 // Fill the infos for the previous clusters if not the same
2019 // detector anymore or if not the same track
2020 if (((detector != fDetectorPreviousTrack) || (!fDetectorAliTRDtrack)) &&
2021 (fDetectorPreviousTrack != -1)) {
2025 // If the same track, then look if the previous detector is in
2026 // the same plane, if yes: not a good track
2027 if (fDetectorAliTRDtrack &&
2028 (GetPlane(detector) <= GetPlane(fDetectorPreviousTrack))) {
2029 fGoodTrack = kFALSE;
2032 // Fill only if the track doesn't touch a masked pad or doesn't
2033 // appear in the middle (fGoodTrack)
2038 FillTheInfoOfTheTrackCH();
2043 FillTheInfoOfTheTrackPH();
2046 } // if a good track
2050 } // Fill at the end the charge
2052 // Calcul the position of the detector
2053 if (detector != fDetectorPreviousTrack) {
2054 LocalisationDetectorXbins(detector);
2057 // Reset the good track for the PRF
2058 Bool_t good = kTRUE;
2060 // Localisation of the cluster
2061 Double_t pos[3] = { 0.0, 0.0, 0.0 };
2062 pos[0] = cl->GetX();
2063 pos[1] = cl->GetY();
2064 pos[2] = cl->GetZ();
2065 Int_t time = cl->GetLocalTimeBin();
2067 // Reset the detector
2068 fDetectorPreviousTrack = detector;
2069 fDetectorAliTRDtrack = kTRUE;
2071 // Position of the cluster
2072 AliTRDpadPlane *padplane = parCom->GetPadPlane(plane,chamber);
2073 Int_t row = padplane->GetPadRowNumber(pos[2]);
2074 Double_t offsetz = padplane->GetPadRowOffset(row,pos[2]);
2075 Double_t offsettilt = padplane->GetTiltOffset(offsetz);
2076 Int_t col = padplane->GetPadColNumber(pos[1] + offsettilt,offsetz);
2078 // See if we are not near a masked pad
2079 if (!IsPadOn(detector,col,row)) {
2081 fGoodTrack = kFALSE;
2085 if (!IsPadOn(detector,col-1,row)) {
2086 fGoodTrack = kFALSE;
2092 if (!IsPadOn(detector,col+1,row)) {
2093 fGoodTrack = kFALSE;
2098 // Row of the cluster and position in the pad groups
2099 Int_t posr[3] = { 0, 0, 0 };
2100 if ((fCH2dOn) && (fNnZ[0] != 0)) {
2101 posr[0] = (Int_t) row / fNnZ[0];
2103 if ((fPH2dOn) && (fNnZ[1] != 0)) {
2104 posr[1] = (Int_t) row / fNnZ[1];
2106 if ((fPRF2dOn) && (fNnZ[2] != 0)) {
2107 posr[2] = (Int_t) row / fNnZ[2];
2110 // Col of the cluster and position in the pad groups
2111 Int_t posc[3] = { 0, 0, 0 };
2112 if ((fCH2dOn) && (fNnRphi[0] != 0)) {
2113 posc[0] = (Int_t) col / fNnRphi[0];
2115 if ((fPH2dOn) && (fNnRphi[1] != 0)) {
2116 posc[1] = (Int_t) col / fNnRphi[1];
2118 if ((fPRF2dOn) && (fNnRphi[2] != 0)) {
2119 posc[2] = (Int_t) col / fNnRphi[2];
2122 // Charge in the cluster
2123 // For the moment take the abs
2124 Float_t q = TMath::Abs(cl->GetQ());
2125 Short_t *signals = cl->GetSignals();
2127 // Correction due to the track angle
2128 Float_t correction = 1.0;
2129 Float_t normalisation = 6.67;
2130 if ((q >0) && (t->GetNdedx() > 0)) {
2131 correction = t->GetClusterdQdl((t->GetNdedx() - 1)) / (q * normalisation);
2134 // Fill the fAmpTotal with the charge
2137 fAmpTotal[(Int_t) (posc[0]*fNfragZ[0]+posr[0])] += q * correction;
2140 fAmpTotal[(Int_t) (posc[0]*fNfragZ[0]+posr[0])] += ((Float_t) signals[3]) * correction;
2144 // Fill the fPHPlace and value
2146 fPHPlace[time] = posc[1]*fNfragZ[1]+posr[1];
2148 fPHValue[time] = q * correction;
2151 fPHValue[time] = ((Float_t) signals[3]) * correction;
2155 // Fill direct the PRF
2156 if ((fPRF2dOn) && (good)) {
2158 Float_t yminus = 0.0;
2159 Float_t xcenter = 0.0;
2160 Float_t ycenter = 0.0;
2162 Bool_t echec = kFALSE;
2164 if ((cl->From3pad()) && (!cl->IsUsed())) {
2166 // Center 3 balanced
2167 if ((((Float_t) signals[3]) > fThresholdClusterPRF2) &&
2168 (((Float_t) signals[2]) > fThresholdClusterPRF2) &&
2169 (((Float_t) signals[4]) > fThresholdClusterPRF2) &&
2170 (((Float_t) signals[1]) < fThresholdClusterPRF1) &&
2171 (((Float_t) signals[5]) < fThresholdClusterPRF1) &&
2172 (q * correction > 130.0)) {
2173 // Col correspond to signals[3]
2174 if (fCenterOfflineCluster) {
2175 xcenter = cl->GetCenter();
2178 // Security of the denomiateur is 0
2179 if ((((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) /
2180 ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))) != 1.0) {
2181 xcenter = 0.5 * (TMath::Log((Float_t) (((Float_t) signals[4]) / ((Float_t) signals[2]))))
2182 / (TMath::Log(((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3])))
2183 / ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))));
2189 if ((xcenter > -0.5) && (xcenter < 0.5)) {
2190 ycenter = (Float_t) (((Float_t) signals[3])
2191 / (((Float_t) signals[2]) + ((Float_t) signals[3]) + (((Float_t) signals[4]))));
2192 yminus = (Float_t) (((Float_t) signals[2])
2193 / (((Float_t) signals[2]) + ((Float_t) signals[3]) + (((Float_t) signals[4]))));
2194 ymax = (Float_t) (((Float_t) signals[4])
2195 / (((Float_t) signals[2]) + ((Float_t) signals[3]) + (((Float_t) signals[4]))));
2196 if ((ycenter > 0.485) &&
2197 (TMath::Abs(((Float_t) signals[2]) + ((Float_t) signals[3]) + (((Float_t) signals[4])) - q) < 8.0)) {
2203 // Fill only if it is in the drift region!
2204 if ((((Float_t) (((Float_t) time) / fSf)) > 0.3) && (echec)) {
2206 fPRF2d->Fill((fXbins[2]+posc[2]*fNfragZ[2]+posr[2]+0.5),xcenter,ycenter);
2207 if (xcenter < 0.0) {
2208 fPRF2d->Fill((fXbins[2]+posc[2]*fNfragZ[2]+posr[2]+0.5),-(xcenter+1.0),yminus);
2210 if (xcenter > 0.0) {
2211 fPRF2d->Fill((fXbins[2]+posc[2]*fNfragZ[2]+posr[2]+0.5),1.0-xcenter,ymax);
2215 UpdateVectorPRF(fXbins[2]+posc[2]*fNfragZ[2]+posr[2],xcenter,ycenter);
2216 if (xcenter < 0.0) {
2217 UpdateVectorPRF(fXbins[2]+posc[2]*fNfragZ[2]+posr[2],-(xcenter+1.0),yminus);
2219 if (xcenter > 0.0) {
2220 UpdateVectorPRF(fXbins[2]+posc[2]*fNfragZ[2]+posr[2],1.0-xcenter,ymax);
2223 } // If in the drift region
2233 //____________Online trackling in AliTRDtrigger________________________________
2234 Bool_t AliTRDCalibra::UpdateHistogramcm(AliTRDmcmTracklet *trk)
2238 // This function will be called in the function AliTRDtrigger::TestTracklet
2239 // before applying the pt cut on the tracklets
2240 // Fill the infos for the tracklets fTrkTest if the tracklets is "good"
2243 // Localisation of the Xbins involved
2244 Int_t idect = trk->GetDetector();
2245 LocalisationDetectorXbins(idect);
2247 // Get the parameter object
2248 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
2250 AliInfo("Could not get calibDB");
2257 // Row of the tracklet and position in the pad groups
2258 Int_t row = trk->GetRow();
2259 Int_t posr[3] = { 0, 0, 0 };
2260 if ((fCH2dOn) && (fNnZ[0] != 0)) {
2261 posr[0] = (Int_t) row / fNnZ[0];
2263 if ((fPH2dOn) && (fNnZ[1] != 0)) {
2264 posr[1] = (Int_t) row / fNnZ[1];
2266 if ((fPRF2dOn) && (fNnZ[2] != 0)) {
2267 posr[2] = (Int_t) row / fNnZ[2];
2270 // Eventuelle correction due to track angle in z direction
2271 Float_t correction = 1.0;
2272 if (fMcmCorrectAngle) {
2273 Float_t z = trk->GetRowz();
2274 Float_t r = trk->GetTime0();
2275 correction = r / TMath::Sqrt((r*r+z*z));
2278 // Boucle sur les clusters
2279 // Condition on number of cluster: don't come from the middle of the detector
2280 if (trk->GetNclusters() >= fNumberClusters) {
2282 for (Int_t icl = 0; icl < trk->GetNclusters(); icl++) {
2284 Float_t amp[3] = { 0.0, 0.0, 0.0 };
2285 Int_t time = trk->GetClusterTime(icl);
2286 Int_t col = trk->GetClusterCol(icl);
2288 amp[0] = trk->GetClusterADC(icl)[0] * correction;
2289 amp[1] = trk->GetClusterADC(icl)[1] * correction;
2290 amp[2] = trk->GetClusterADC(icl)[2] * correction;
2292 if ((amp[0] < 0.0) ||
2298 // Col of cluster and position in the pad groups
2299 Int_t posc[3] = { 0, 0, 0 };
2300 if ((fCH2dOn) && (fNnRphi[0] != 0)) {
2301 posc[0] = (Int_t) col / fNnRphi[0];
2303 if ((fPH2dOn) && (fNnRphi[1] != 0)) {
2304 posc[1] = (Int_t) col / fNnRphi[1];
2306 if ((fPRF2dOn) && (fNnRphi[2] != 0)) {
2307 posc[2] = (Int_t) col / fNnRphi[2];
2310 // See if we are not near a masked pad
2311 Bool_t good = kTRUE;
2312 if (!IsPadOn(idect,col,row)) {
2313 fGoodTrack = kFALSE;
2318 if (!IsPadOn(idect,col-1,row)) {
2319 fGoodTrack = kFALSE;
2325 if (!IsPadOn(idect,col+1,row)) {
2326 fGoodTrack = kFALSE;
2333 fPHPlace[time] = posc[1] * fNfragZ[1] + posr[1];
2338 fAmpTotal[(Int_t) (posc[0]*fNfragZ[0]+posr[0])] += (Float_t) (amp[0]+amp[1]+amp[2]);
2341 fPHValue[time] = (Float_t) (amp[0]+amp[1]+amp[2]);
2346 fAmpTotal[(Int_t) (posc[0]*fNfragZ[0]+posr[0])] += (Float_t) amp[1];
2349 fPHValue[time] = amp[1];
2354 if (fPRF2dOn && good) {
2355 if ((amp[0] > fThresholdClusterPRF2) &&
2356 (amp[1] > fThresholdClusterPRF2) &&
2357 (amp[2] > fThresholdClusterPRF2) &&
2358 ((amp[0]+amp[1]+amp[2]) > 130.0)) {
2359 // Security of the denomiateur is 0
2360 if ((((Float_t) (((Float_t) amp[1]) * ((Float_t) amp[1])))
2361 / ((Float_t) (((Float_t) amp[0]) * ((Float_t) amp[2])))) != 1.0) {
2362 Float_t xcenter = 0.5 * (TMath::Log(amp[2] / amp[0]))
2363 / (TMath::Log((amp[1]*amp[1]) / (amp[0]*amp[2])));
2364 Float_t ycenter = amp[1] / (amp[0] + amp[1] + amp[2]);
2365 if ((xcenter > -0.5) &&
2367 (ycenter > 0.485)) {
2368 Float_t yminus = amp[0] / (amp[0]+amp[1]+amp[2]);
2369 Float_t ymax = amp[2] / (amp[0]+amp[1]+amp[2]);
2370 // Fill only if it is in the drift region!
2371 if (((Float_t) time / fSf) > 0.3) {
2373 fPRF2d->Fill((fXbins[2]+posc[2]*fNfragZ[2]+posr[2]+0.5),xcenter,ycenter);
2374 if (xcenter < 0.0) {
2375 fPRF2d->Fill((fXbins[2]+posc[2]*fNfragZ[2]+posr[2]+0.5),-(xcenter+1.0),yminus);
2377 if (xcenter > 0.0) {
2378 fPRF2d->Fill((fXbins[2]+posc[2]*fNfragZ[2]+posr[2]+0.5),(1.0-xcenter),ymax);
2382 UpdateVectorPRF((fXbins[2]+posc[2]*fNfragZ[2]+posr[2]),xcenter,ycenter);
2383 if (xcenter < 0.0) {
2384 UpdateVectorPRF(fXbins[2]+posc[2]*fNfragZ[2]+posr[2],-(xcenter+1.0),yminus);
2386 if (xcenter > 0.0) {
2387 UpdateVectorPRF(fXbins[2]+posc[2]*fNfragZ[2]+posr[2],(1.0-xcenter),ymax);
2396 } // Boucle clusters
2399 if (fCH2dOn && fGoodTrack) {
2400 FillTheInfoOfTheTrackCH();
2404 if (fPH2dOn && fGoodTrack) {
2405 FillTheInfoOfTheTrackPH();
2408 } // Condition on number of clusters
2414 //____________Functions for seeing if the pad is really okey___________________
2416 //_____________________________________________________________________________
2417 Bool_t AliTRDCalibra::SetModeCalibrationFromTObject(TObject *object, Int_t i)
2420 // Set fNz[i] and fNrphi[i] of the AliTRDCalibra::Instance()
2421 // corresponding to the given TObject
2424 const char *nametitle = object->GetTitle();
2427 const Char_t *patternz0 = "Nz0";
2428 const Char_t *patternz1 = "Nz1";
2429 const Char_t *patternz2 = "Nz2";
2430 const Char_t *patternz3 = "Nz3";
2431 const Char_t *patternz4 = "Nz4";
2432 const Char_t *patternrphi0 = "Nrphi0";
2433 const Char_t *patternrphi1 = "Nrphi1";
2434 const Char_t *patternrphi2 = "Nrphi2";
2435 const Char_t *patternrphi3 = "Nrphi3";
2436 const Char_t *patternrphi4 = "Nrphi4";
2437 const Char_t *patternrphi5 = "Nrphi5";
2438 const Char_t *patternrphi6 = "Nrphi6";
2441 UShort_t testrphi = 0;
2444 if (strstr(nametitle,patternz0)) {
2448 if (strstr(nametitle,patternz1)) {
2452 if (strstr(nametitle,patternz2)) {
2456 if (strstr(nametitle,patternz3)) {
2460 if (strstr(nametitle,patternz4)) {
2466 if (strstr(nametitle,patternrphi0)) {
2470 if (strstr(nametitle,patternrphi1)) {
2474 if (strstr(nametitle,patternrphi2)) {
2478 if (strstr(nametitle,patternrphi3)) {
2482 if (strstr(nametitle,patternrphi4)) {
2486 if (strstr(nametitle,patternrphi5)) {
2490 if (strstr(nametitle,patternrphi6)) {
2495 // Look if all is okey
2508 //_____________________________________________________________________________
2509 Bool_t AliTRDCalibra::IsPadOn(Int_t detector, Int_t col, Int_t row) const
2512 // Look in the choosen database if the pad is On.
2513 // If no the track will be "not good"
2516 // Get the parameter object
2517 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
2519 AliInfo("Could not get calibDB");
2524 Int_t colmcm = (Int_t) col / npads;
2526 if (!cal->IsChamberInstalled(detector) ||
2527 cal->IsChamberMasked(detector) ||
2528 cal->IsMCMMasked(detector,colmcm,row) ||
2529 cal->IsPadMasked(detector,col,row)) {
2538 //____________Functions for plotting the 2D____________________________________
2540 //_____________________________________________________________________________
2541 void AliTRDCalibra::Plot2d()
2544 // Plot the 2D histos
2559 //____________Writing the 2D___________________________________________________
2561 //_____________________________________________________________________________
2562 Bool_t AliTRDCalibra::Write2d()
2565 // Write the 2D histograms or the vectors converted in trees in the file
2566 // "TRD.calibration.root"
2569 TFile *fout = TFile::Open(fWriteName,"RECREATE");
2570 // Check if the file could be opened
2571 if (!fout || !fout->IsOpen()) {
2572 AliInfo("No File found!");
2575 AliInfo(Form("Numbertrack: %d Numberusedch[0]: %d, Numberusedch[1]: %d Numberusedph[0]: %d, Numberusedph[1]: %d"
2580 ,fNumberUsedPh[1]));
2582 TStopwatch stopwatch;
2586 if ((fCH2dOn ) && (fWrite[0])) {
2588 fout->WriteTObject(fCH2d);
2595 TTree *treeCH2d = ConvertVectorCTTreeHisto(fVectorCH,fPlaCH,"treeCH2d",(const char *) name);
2596 fout->WriteTObject(treeCH2d);
2599 if ((fPH2dOn ) && (fWrite[1])) {
2601 fout->WriteTObject(fPH2d);
2608 TTree *treePH2d = ConvertVectorPTreeHisto(fVectorPH,fPlaPH,"treePH2d",(const char *) name);
2609 fout->WriteTObject(treePH2d);
2612 if ((fPRF2dOn ) && (fWrite[2])) {
2614 fout->WriteTObject(fPRF2d);
2621 TTree *treePRF2d = ConvertVectorPTreeHisto(fVectorPRF,fPlaPRF,"treePRF2d",(const char *) name);
2622 fout->WriteTObject(treePRF2d);
2628 AliInfo(Form("Execution time Write2d: R:%.2fs C:%.2fs"
2629 ,stopwatch.RealTime(),stopwatch.CpuTime()));
2635 //_____________________________________________________________________________
2636 AliTRDCalDet *AliTRDCalibra::CreateDetObjectTree(TTree *tree, Int_t i)
2639 // It creates the AliTRDCalDet object from the tree of the coefficient
2640 // for the calibration i (i != 2)
2641 // It takes the mean value of the coefficients per detector
2642 // This object has to be written in the database
2645 // Create the DetObject
2646 AliTRDCalDet *object = 0x0;
2648 object = new AliTRDCalDet("ChamberGainFactor","GainFactor (detector value)");
2651 object = new AliTRDCalDet("ChamberVdrift","TRD drift velocities (detector value)");
2654 object = new AliTRDCalDet("ChamberT0","T0 (detector value)");
2658 Int_t detector = -1;
2659 Float_t values[2304];
2660 tree->SetBranchAddress("detector",&detector);
2662 tree->SetBranchAddress("gainPad",values);
2665 tree->SetBranchAddress("vdrift" ,values);
2668 tree->SetBranchAddress("t0" ,values);
2671 // For calculating the mean
2674 Int_t numberofentries = tree->GetEntries();
2676 if (numberofentries != 540) {
2677 AliInfo("The tree is not complete");
2680 for (Int_t det = 0; det < numberofentries; ++det) {
2681 tree->GetEntry(det);
2682 if (GetChamber(detector) == 2) {
2689 for (Int_t k = 0; k < nto; k++) {
2690 mean += TMath::Abs(values[k]) / nto;
2692 object->SetValue(detector,mean);
2699 //_____________________________________________________________________________
2700 TObject *AliTRDCalibra::CreatePadObjectTree(TTree *tree, Int_t i
2701 , AliTRDCalDet *detobject)
2704 // It Creates the AliTRDCalPad object from the tree of the
2705 // coefficient for the calibration i (i != 2)
2706 // You need first to create the object for the detectors,
2707 // where the mean value is put.
2708 // This object has to be written in the database
2711 // Create the DetObject
2712 AliTRDCalPad *object = 0x0;
2714 object = new AliTRDCalPad("GainFactor","GainFactor (local variations)");
2717 object = new AliTRDCalPad("LocalVdrift","TRD drift velocities (local variations)");
2720 object = new AliTRDCalPad("LocalT0","T0 (local variations)");
2724 Int_t detector = -1;
2725 Float_t values[2304];
2726 tree->SetBranchAddress("detector",&detector);
2728 tree->SetBranchAddress("gainPad",values);
2731 tree->SetBranchAddress("vdrift" ,values);
2734 tree->SetBranchAddress("t0" ,values);
2739 Int_t numberofentries = tree->GetEntries();
2741 if (numberofentries != 540) {
2742 AliInfo("The tree is not complete");
2745 for (Int_t det = 0; det < numberofentries; ++det) {
2746 tree->GetEntry(det);
2747 AliTRDCalROC *calROC = object->GetCalROC(detector);
2748 mean = detobject->GetValue(detector);
2752 Int_t rowMax = calROC->GetNrows();
2753 Int_t colMax = calROC->GetNcols();
2754 for (Int_t row = 0; row < rowMax; ++row) {
2755 for (Int_t col = 0; col < colMax; ++col) {
2757 calROC->SetValue(col,row,TMath::Abs(values[(Int_t) (col*rowMax+row)])/mean);
2767 //_____________________________________________________________________________
2768 TObject *AliTRDCalibra::CreatePadObjectTree(TTree *tree)
2771 // It Creates the AliTRDCalPad object from the tree of the
2772 // coefficient for the calibration PRF (i = 2)
2773 // This object has to be written in the database
2776 // Create the DetObject
2777 AliTRDCalPad *object = new AliTRDCalPad("PRFWidth","PRFWidth");
2780 Int_t detector = -1;
2781 Float_t values[2304];
2782 tree->SetBranchAddress("detector",&detector);
2783 tree->SetBranchAddress("width" ,values);
2786 Int_t numberofentries = tree->GetEntries();
2788 if (numberofentries != 540) {
2789 AliInfo("The tree is not complete");
2792 for (Int_t det = 0; det < numberofentries; ++det) {
2793 tree->GetEntry(det);
2794 AliTRDCalROC *calROC = object->GetCalROC(detector);
2795 Int_t rowMax = calROC->GetNrows();
2796 Int_t colMax = calROC->GetNcols();
2797 for (Int_t row = 0; row < rowMax; ++row) {
2798 for (Int_t col = 0; col < colMax; ++col) {
2799 calROC->SetValue(col,row,TMath::Abs(values[(Int_t) (col*rowMax+row)]));
2808 //_____________________________________________________________________________
2809 void AliTRDCalibra::SetRelativeScale(Float_t RelativeScale)
2812 // Set the factor that will divide the deposited charge
2813 // to fit in the histo range [0,300]
2816 if (RelativeScale > 0.0) {
2817 fRelativeScale = RelativeScale;
2820 AliInfo("RelativeScale must be strict positif!");
2825 //_____________________________________________________________________________
2826 void AliTRDCalibra::SetNz(Int_t i, Short_t Nz)
2829 // Set the mode of calibration group in the z direction for the parameter i
2837 AliInfo("You have to choose between 0 and 4");
2842 //_____________________________________________________________________________
2843 void AliTRDCalibra::SetNrphi(Int_t i, Short_t Nrphi)
2846 // Set the mode of calibration group in the rphi direction for the parameter i
2854 AliInfo("You have to choose between 0 and 6");
2859 //_____________________________________________________________________________
2860 void AliTRDCalibra::SetPeriodeFitPH(Int_t periodeFitPH)
2863 // Set FitPH if 1 then each detector will be fitted
2866 if (periodeFitPH > 0) {
2867 fFitPHPeriode = periodeFitPH;
2870 AliInfo("periodeFitPH must be higher than 0!");
2875 //_____________________________________________________________________________
2876 void AliTRDCalibra::SetBeginFitCharge(Float_t beginFitCharge)
2879 // The fit of the deposited charge distribution begins at
2880 // histo->Mean()/beginFitCharge
2881 // You can here set beginFitCharge
2884 if (beginFitCharge > 0) {
2885 fBeginFitCharge = beginFitCharge;
2888 AliInfo("beginFitCharge must be strict positif!");
2893 //_____________________________________________________________________________
2894 void AliTRDCalibra::SetT0Shift(Float_t t0Shift)
2897 // The t0 calculated with the maximum positif slope is shift from t0Shift
2898 // You can here set t0Shift
2905 AliInfo("t0Shift must be strict positif!");
2910 //_____________________________________________________________________________
2911 void AliTRDCalibra::SetRangeFitPRF(Float_t rangeFitPRF)
2914 // The fit of the PRF is from -rangeFitPRF to rangeFitPRF
2915 // You can here set rangeFitPRF
2918 if ((rangeFitPRF > 0) &&
2919 (rangeFitPRF <= 1.0)) {
2920 fRangeFitPRF = rangeFitPRF;
2923 AliInfo("rangeFitPRF must be between 0 and 1.0");
2928 //_____________________________________________________________________________
2929 void AliTRDCalibra::SetRebin(Short_t rebin)
2932 // Rebin with rebin time less bins the Ch histo
2933 // You can set here rebin that should divide the number of bins of CH histo
2938 AliInfo("You have to be sure that fRebin divides fNumberBinCharge used!");
2941 AliInfo("You have to choose a positiv value!");
2946 //_____________________________________________________________________________
2947 TTree *AliTRDCalibra::Sum2Trees(const Char_t *filename1
2948 , const Char_t *filename2
2949 , const Char_t *variablecali)
2952 // It returns the sum of two trees with the name variablecali
2953 // in the files filenam1 and filename2 equivalent of merging two 2D histos
2954 // The name of the resulting tree is the same as the two input trees
2955 // variablecali can be treeCH2d, treePH2d or treePRF2d
2959 TChain *treeChain = new TChain(variablecali);
2960 TObjArray *vectorplace = new TObjArray();
2961 TObjArray *where = new TObjArray();
2965 TFile *file1 = new TFile(filename1,"READ");
2966 TTree *tree1 = (TTree *) file1->Get(variablecali);
2971 vectorplace = ConvertTreeVector(tree1);
2973 // Say where it is in tree 1
2974 for (Int_t jui = 0; jui < (Int_t) vectorplace->GetEntriesFast(); jui++) {
2975 AliTRDPlace *placejui = new AliTRDPlace();
2976 placejui->SetPlace(jui);
2977 TObjArray *chainplace = new TObjArray();
2978 chainplace->Add((TObject *) placejui);
2979 where->Add((TObject *) chainplace);
2983 treeChain->Add(filename1);
2988 TFile *file2 = new TFile(filename2,"READ");
2989 TTree *tree2 = (TTree *) file2->Get(variablecali);
2994 TObjArray *vector2 = ConvertTreeVector(tree2);
2995 Int_t j = treeChain->GetEntries();
2997 for (Int_t jui = 0; jui < (Int_t) vector2->GetEntriesFast(); jui++) {
2998 // Search if already found
2999 Int_t place = SearchInTreeVector(vectorplace,((AliTRDPlace *) vector2->At(jui))->GetPlace());
3000 // Create a new element in the two std vectors
3002 AliTRDPlace *placejjui = new AliTRDPlace();
3003 placejjui->SetPlace((j+jui));
3004 TObjArray *chainplace = new TObjArray();
3005 chainplace->Add((TObject *) placejjui);
3006 vectorplace->Add((TObject *) (vector2->At(jui)));
3007 where->Add((TObject *) chainplace);
3009 // Update the element at the place "place" in the std vector whereinthechain
3011 AliTRDPlace *placejjui = new AliTRDPlace();
3012 placejjui->SetPlace((j+jui));
3013 TObjArray *chainplace = ((TObjArray *) where->At(place));
3014 chainplace->Add((TObject *) placejjui);
3015 where->AddAt((TObject *) chainplace,place);
3020 treeChain->Add(filename2);
3023 // Take care of the profile
3024 const Char_t *pattern = "P";
3027 if (!strstr(variablecali,pattern)) {
3029 // Ready to read the chain
3031 treeChain->SetBranchAddress("histo",&his);
3033 // Initialise the final tree
3035 TH1F *histsum = 0x0;
3037 tree = new TTree(variablecali,variablecali);
3038 tree->Branch("groupnumber",&group,"groupnumber/I");
3039 tree->Branch("histo","TH1F",&histsum,32000,0);
3042 if (treeChain->GetEntries() < 1) {
3046 for (Int_t h = 0; h < (Int_t) vectorplace->GetEntriesFast(); h++) {
3047 group = ((AliTRDPlace *) vectorplace->At(h))->GetPlace();
3048 TObjArray *chainplace = ((TObjArray *) where->At(h));
3049 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(0))->GetPlace());
3050 //Init for the first time
3052 histsum = new TH1F("","",his->GetXaxis()->GetNbins()
3053 ,his->GetXaxis()->GetBinLowEdge(1)
3054 ,his->GetXaxis()->GetBinUpEdge(his->GetXaxis()->GetNbins()));
3057 // Reset for each new group
3058 histsum->SetEntries(0.0);
3059 for (Int_t l = 0; l <= histsum->GetXaxis()->GetNbins(); l++) {
3060 histsum->SetBinContent(l,0.0);
3061 histsum->SetBinError(l,0.0);
3063 histsum->Add(his,1);
3064 if ((Int_t) chainplace->GetEntriesFast() > 1) {
3065 for (Int_t s = 1; s < (Int_t) chainplace->GetEntriesFast(); s++) {
3066 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(s))->GetPlace());
3067 histsum->Add(his,1);
3076 // Ready to read the chain
3077 TGraphErrors *his = 0x0;
3078 treeChain->SetBranchAddress("histo",&his);
3080 // Initialise the final tree
3082 TGraphErrors *histsum = 0x0;
3083 Double_t *xref = 0x0;
3085 tree = new TTree(variablecali,variablecali);
3086 tree->Branch("groupnumber",&group,"groupnumber/I");
3087 tree->Branch("histo","TGraphErrors",&histsum,32000,0);
3090 if (treeChain->GetEntries() < 1) {
3094 for (Int_t h = 0; h < (Int_t) vectorplace->GetEntriesFast(); h++) {
3096 group = ((AliTRDPlace *) vectorplace->At(h))->GetPlace();
3097 TObjArray *chainplace = ((TObjArray *) where->At(h));
3098 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(0))->GetPlace());
3099 //Init or reset for a new group
3100 Int_t nbins = his->GetN();
3102 x = new Double_t[nbins];
3105 ex = new Double_t[nbins];
3107 y = new Double_t[nbins];
3109 ey = new Double_t[nbins];
3111 for (Int_t lo = 0; lo < nbins; lo++) {
3118 histsum = new TGraphErrors(nbins,x,y,ex,ey);
3121 histsum = AddProfiles(his,histsum);
3122 if ((Int_t) chainplace->GetEntriesFast() > 1) {
3123 for (Int_t s = 1; s < (Int_t) chainplace->GetEntriesFast(); s++) {
3124 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(s))->GetPlace());
3125 histsum = AddProfiles(his,histsum);
3139 //____________Function fill 2D for the moment out of the code__________________
3141 //____________Function fill 2D all objects from digits_________________________
3142 Bool_t AliTRDCalibra::Create2DDiSimOnline(Int_t iev1, Int_t iev2)
3145 // Only for simulations, after the simulation, create the 2D histos
3146 // from the digits stored in the file "TRD.Digits.root"
3147 // Only for CH and PH
3150 const Int_t kNplan = 6;
3151 const Int_t kNcham = 5;
3153 // RunLoader and so on
3155 delete gAlice->GetRunLoader();
3160 AliRunLoader *rl = AliRunLoader::Open("galice.root");
3166 gAlice = rl->GetAliRun();
3171 // Import the Trees for the event nEvent in the file
3172 rl->LoadKinematics();
3176 AliLoader *loader = rl->GetLoader("TRDLoader");
3178 AliInfo("No TRDLLoader found!");
3182 // Get the pointer to the TRD detector
3183 AliTRD *trd = (AliTRD *) gAlice->GetDetector("TRD");
3185 AliInfo("No TRD detector found");
3189 // Get the pointer to the geometry object
3190 AliTRDgeometry *geo;
3192 geo = trd->GetGeometry();
3195 AliInfo("No TRD geometry found");
3200 AliCDBManager *man = AliCDBManager::Instance();
3202 AliInfo("Could not get CDB Manager");
3206 // Get the parameter object
3207 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
3209 AliInfo("Could not get CommonParam");
3212 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
3214 AliInfo("Could not get calibDB");
3219 fTimeMax = cal->GetNumberOfTimeBins();
3220 fSf = (Float_t) cal->GetSamplingFrequency();
3221 if (fRelativeScaleAuto) {
3225 if (fRelativeScale <= 0.0) {
3226 AliInfo("You have to set the relativescale factor per hand!");
3231 // Create the 2D histos corresponding to the pad group calibration mode
3234 AliInfo(Form("We will fill the CH2d histo with the pad calibration mode: Nz %d, and Nrphi %d"
3238 // Calcul the number of Xbins
3240 ModePadCalibration(2,0);
3241 ModePadFragmentation(0,2,0,0);
3242 fDetChamb2[0] = fNfragZ[0] * fNfragRphi[0];
3243 fNtotal[0] += 6 * 18 * fDetChamb2[0];
3244 ModePadCalibration(0,0);
3245 ModePadFragmentation(0,0,0,0);
3246 fDetChamb0[0] = fNfragZ[0] * fNfragRphi[0];
3247 fNtotal[0] += 6 * 4 * 18 * fDetChamb0[0];
3248 AliInfo(Form("Total number of Xbins: %d",fNtotal[0]));
3250 // Create the 2D histo
3252 CreateCH2d(fNtotal[0]);
3255 fVectorCH = new TObjArray();
3256 fPlaCH = new TObjArray();
3263 AliInfo(Form("We will fill the PH2d histo with the pad calibration mode: Nz %d, and Nrphi %d"
3267 // Calcul the number of Xbins
3269 ModePadCalibration(2,1);
3270 ModePadFragmentation(0,2,0,1);
3271 fDetChamb2[1] = fNfragZ[1] * fNfragRphi[1];
3272 fNtotal[1] += 6 * 18 * fDetChamb2[1];
3273 ModePadCalibration(0,1);
3274 ModePadFragmentation(0,0,0,1);
3275 fDetChamb0[1] = fNfragZ[1] * fNfragRphi[1];
3276 fNtotal[1] += 6 * 4 * 18 * fDetChamb0[1];
3277 AliInfo(Form("Total number of Xbins: %d",fNtotal[1]));
3279 // Create the 2D histo
3281 CreatePH2d(fNtotal[1]);
3284 fVectorPH = new TObjArray();
3285 fPlaPH = new TObjArray();
3290 loader->LoadDigits();
3291 AliInfo("LoadDigits ");
3292 AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager();
3294 //iev2 egal to the max if 0
3296 iev2 = rl->GetNumberOfEvents();
3297 AliInfo(Form("Total number of events: %d",iev2));
3301 for (Int_t ievent = iev1; ievent < iev2; ievent++) {
3302 AliInfo(Form("Process event %d",ievent));
3303 rl->GetEvent(ievent);
3304 if (!loader->TreeD()) {
3305 AliInfo("loader Loading Digits ... ");
3306 loader->LoadDigits();
3308 digitsManager->ReadDigits(loader->TreeD());
3309 AliInfo("digitsManager Read Digits Done");
3310 // Read the digits from the file
3311 if (!(digitsManager->ReadDigits(loader->TreeD()))) {
3316 for (Int_t iSect = 0; iSect < 18; iSect++) {
3317 for (Int_t iChamb = 0; iChamb < kNcham; iChamb++) {
3318 for (Int_t iPlane = 0; iPlane < kNplan; iPlane++) {
3320 // A little geometry:
3321 Int_t iDet = geo->GetDetector(iPlane,iChamb,iSect);
3322 Int_t rowMax = parCom->GetRowMax(iPlane,iChamb,iSect);
3323 Int_t colMax = parCom->GetColMax(iPlane);
3325 // Variables for the group
3326 LocalisationDetectorXbins(iDet);
3328 // In the cas of charge
3330 amptotal = new Float_t[fNfragRphi[0]*fNfragZ[0]];
3332 for (Int_t k = 0; k < fNfragRphi[0]*fNfragZ[0]; k++) {
3337 // Loop through the detector pixel
3338 for (Int_t time = 0; time < fTimeMax; time++) {
3339 for (Int_t col = 0; col < colMax; col++) {
3340 for (Int_t row = 0; row < rowMax; row++) {
3342 // Amplitude and position in pad group
3343 AliTRDdigit *digit = digitsManager->GetDigit(row,col,time,iDet);
3344 Int_t amp = digit->GetAmp();
3345 Int_t posr[2] = {0,0};
3346 Int_t posc[2] = {0,0};
3349 posr[0] = (Int_t) row / fNnZ[0];
3352 (fNnRphi[0] != 0)) {
3353 posc[0] = (Int_t) col / fNnRphi[0];
3357 posr[1] = (Int_t) row / fNnZ[1];
3360 (fNnRphi[1] != 0)) {
3361 posc[1] = (Int_t) col / fNnRphi[1];
3366 if (amp < fThresholdDigit) {
3369 amptotal[(Int_t) (posc[0]*fNfragZ[0]+posr[0])] += amp;
3373 fPH2d->Fill((fXbins[1]+posc[1]*fNfragZ[1]+posr[1])+0.5,(Float_t) time/fSf,(Double_t) amp);
3376 UpdateVectorPH((fXbins[1]+posc[1]*fNfragZ[1]+posr[1]),time,(Double_t) amp);
3389 // If automatic scale
3390 if ((fCountRelativeScale < 100) && (fRelativeScaleAuto)) {
3391 // Take only the one zone track
3392 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3393 if ((fCountRelativeScale < 100) && (amptotal[k] > 2.0)) {
3394 fRelativeScale += amptotal[k]*0.014*0.01;
3395 fCountRelativeScale++;
3400 // We fill the CH2d after having scale with the first 100
3401 if ((fCountRelativeScale >= 100) && (fRelativeScaleAuto)) {
3403 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3405 (amptotal[k] > 0.0)) {
3406 fCH2d->Fill(fXbins[0]+k+0.5,amptotal[k]/fRelativeScale);
3409 (amptotal[k] > 0.0)) {
3410 UpdateVectorCH(fXbins[0]+k ,amptotal[k]/fRelativeScale);
3415 // No relative salce
3416 if (!fRelativeScaleAuto) {
3417 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3419 (amptotal[k] > 0.0)) {
3420 fCH2d->Fill(fXbins[0]+k+0.5, amptotal[k]/fRelativeScale);
3423 (amptotal[k] > 0.0)) {
3424 UpdateVectorCH(fXbins[0]+k, amptotal[k]/fRelativeScale);
3437 loader->UnloadDigits();
3442 if (fPH2dOn && fHisto2d) {
3445 if (fCH2dOn && fHisto2d) {
3450 if (fWrite[0] || fWrite[1]) {
3452 TFile *fout = TFile::Open(fWriteName,"RECREATE");
3453 // Check if the file could be opened
3454 if (!fout || !fout->IsOpen()) {
3455 AliInfo("<No File found!");
3459 if (fCH2dOn && fHisto2d && fWrite[0]) {
3460 fout->WriteTObject(fCH2d);
3462 if (fPH2dOn && fHisto2d && fWrite[1]) {
3463 fout->WriteTObject(fPH2d);
3466 if (fVector2d && fCH2dOn && fWrite[0]) {
3471 TTree *treeCH2d = ConvertVectorCTTreeHisto(fVectorCH,fPlaCH,"treeCH2d",(const char *) name);
3472 fout->WriteTObject(treeCH2d);
3475 if (fVector2d && fPH2dOn && fWrite[1]) {
3480 TTree *treePH2d = ConvertVectorPTreeHisto(fVectorPH,fPlaPH,"treePH2d",(const char *) name);
3481 fout->WriteTObject(treePH2d);
3492 //____________Function fill 2D all objects from Raw Data_______________________
3493 Bool_t AliTRDCalibra::Create2DRaDaOnline(Int_t iev1, Int_t iev2)
3496 // After having written the RAW DATA in the current directory, create the
3497 // 2D histos from these RAW DATA
3498 // Only for CH and PH
3501 const Int_t kNplan = 6;
3502 const Int_t kNcham = 5;
3503 TString dirname(".");
3506 AliCDBManager *man = AliCDBManager::Instance();
3508 AliInfo("Could not get CDB Manager");
3512 // Get the parameter object
3513 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
3515 AliInfo("Could not get CommonParam");
3519 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
3521 AliInfo("Could not get calibDB");
3526 fTimeMax = cal->GetNumberOfTimeBins();
3527 fSf = (Float_t) cal->GetSamplingFrequency();
3528 if (fRelativeScaleAuto) {
3532 if (fRelativeScale <= 0.0) {
3533 AliInfo("You have to set the relativescale factor per hand!");
3538 // Create the 2D histo corresponding to the pad group calibration mode
3541 AliInfo(Form("We will fill the CH2d histo with the pad calibration mode: Nz %d, and Nrphi %d"
3545 // Calcul the number of Xbins
3547 ModePadCalibration(2,0);
3548 ModePadFragmentation(0,2,0,0);
3549 fDetChamb2[0] = fNfragZ[0] * fNfragRphi[0];
3550 fNtotal[0] += 6 * 18 * fDetChamb2[0];
3551 ModePadCalibration(0,0);
3552 ModePadFragmentation(0,0,0,0);
3553 fDetChamb0[0] = fNfragZ[0] * fNfragRphi[0];
3554 fNtotal[0] += 6 * 4 * 18 * fDetChamb0[0];
3555 AliInfo(Form("Total number of Xbins: %d",fNtotal[0]));
3557 // Create the 2D histo
3559 CreateCH2d(fNtotal[0]);
3562 fVectorCH = new TObjArray();
3563 fPlaCH = new TObjArray();
3570 AliInfo(Form("We will fill the PH2d histo with the pad calibration mode: Nz %d, and Nrphi %d"
3574 // Calcul the number of Xbins
3576 ModePadCalibration(2,1);
3577 ModePadFragmentation(0,2,0,1);
3578 fDetChamb2[1] = fNfragZ[1] * fNfragRphi[1];
3579 fNtotal[1] += 6 * 18 * fDetChamb2[1];
3580 ModePadCalibration(0,1);
3581 ModePadFragmentation(0,0,0,1);
3582 fDetChamb0[1] = fNfragZ[1] * fNfragRphi[1];
3583 fNtotal[1] += 6 * 4 * 18 * fDetChamb0[1];
3584 AliInfo(Form("Total number of Xbins: %d",fNtotal[1]));
3586 // Create the 2D histo
3588 CreatePH2d(fNtotal[1]);
3591 fVectorPH = new TObjArray();
3592 fPlaPH = new TObjArray();
3597 AliTRDrawData *rawdata = new AliTRDrawData();
3598 AliInfo("AliTRDrawData object created ");
3601 for (Int_t ievent = iev1; ievent < iev2; ievent++) {
3604 AliRawReaderFile *readerfile = new AliRawReaderFile(dirname,ievent);
3606 AliInfo("No readerfile found!");
3610 AliTRDdigitsManager *digitsManager = rawdata->Raw2Digits((AliRawReader *) readerfile);
3611 if (!digitsManager) {
3612 AliInfo("No DigitsManager done!");
3616 // Loop on detectors
3617 for (Int_t iSect = 0; iSect < 18; iSect++) {
3618 for (Int_t iPlane = 0; iPlane < kNplan; iPlane++) {
3619 for (Int_t iChamb = 0; iChamb < kNcham; iChamb++) {
3621 // A little geometry:
3622 Int_t iDet = AliTRDgeometry::GetDetector(iPlane,iChamb,iSect);
3623 Int_t rowMax = parCom->GetRowMax(iPlane,iChamb,iSect);
3624 Int_t colMax = parCom->GetColMax(iPlane);
3626 // Variables for the group
3627 LocalisationDetectorXbins(iDet);
3629 // In the cas of charge
3631 amptotal = new Float_t[fNfragRphi[0]*fNfragZ[0]];
3633 for (Int_t k = 0; k < fNfragRphi[0]*fNfragZ[0]; k++) {
3638 // Loop through the detector pixel
3639 for (Int_t time = 0; time < fTimeMax; time++) {
3640 for (Int_t col = 0; col < colMax; col++) {
3641 for (Int_t row = 0; row < rowMax; row++) {
3643 // Amplitude and position of the digit
3644 AliTRDdigit *digit = digitsManager->GetDigit(row,col,time,iDet);
3645 Int_t amp = digit->GetAmp();
3646 Int_t posr[2] = { 0, 0 };
3647 Int_t posc[2] = { 0, 0 };
3650 posr[0] = (Int_t) row / fNnZ[0];
3653 (fNnRphi[0] != 0)) {
3654 posc[0] = (Int_t) col / fNnRphi[0];
3658 posr[1] = (Int_t) row / fNnZ[1];
3661 (fNnRphi[1] != 0)) {
3662 posc[1] = (Int_t) col / fNnRphi[1];
3667 if (amp < fThresholdDigit) {
3670 amptotal[(Int_t) (posc[0]*fNfragZ[0]+posr[0])] += amp;
3675 fPH2d->Fill((fXbins[1]+posc[1]*fNfragZ[1]+posr[1])+0.5,(Float_t)time/fSf,amp);
3678 UpdateVectorPH(fXbins[1]+posc[1]*fNfragZ[1]+posr[1],time,amp);
3690 // If automatic scale
3691 if ((fCountRelativeScale < 100) && (fRelativeScaleAuto)) {
3692 // Take only the one zone track
3693 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3694 if ((fCountRelativeScale < 100) && (amptotal[k] > 2.0)) {
3695 fRelativeScale += amptotal[k] * 0.014 * 0.01;
3696 fCountRelativeScale++;
3701 // We fill the CH2d after having scale with the first 100
3702 if ((fCountRelativeScale >= 100) && (fRelativeScaleAuto)) {
3704 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3705 if (fHisto2d && (amptotal[k] > 0.0)) {
3706 fCH2d->Fill(fXbins[0]+k+0.5,amptotal[k]/fRelativeScale);
3708 if (fVector2d && (amptotal[k] > 0.0)) {
3709 UpdateVectorCH(fXbins[0]+k, amptotal[k]/fRelativeScale);
3714 // No relative salce
3715 if (!fRelativeScaleAuto) {
3716 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3718 (amptotal[k] > 0.0)) {
3719 fCH2d->Fill(fXbins[0]+k+0.5,amptotal[k]/fRelativeScale);
3722 (amptotal[k] > 0.0)) {
3723 UpdateVectorCH(fXbins[0]+k, amptotal[k]/fRelativeScale);
3736 delete digitsManager;
3742 if (fPH2dOn && fHisto2d) {
3745 if (fCH2dOn && fHisto2d) {
3750 if (fWrite[0] || fWrite[1]) {
3752 TFile *fout = TFile::Open(fWriteName,"UPDATE");
3753 // Check if the file could be opened
3754 if (!fout || !fout->IsOpen()) {
3755 AliInfo("<No File found!");
3759 if (fCH2dOn && fHisto2d && fWrite[0]) {
3760 fout->WriteTObject(fCH2d);
3762 if (fPH2dOn && fHisto2d && fWrite[1]) {
3763 fout->WriteTObject(fPH2d);
3766 if (fVector2d && fCH2dOn && fWrite[0]) {
3771 TTree *treeCH2d = ConvertVectorCTTreeHisto(fVectorCH,fPlaCH,"treeCH2d",(const Char_t *) name);
3772 fout->WriteTObject(treeCH2d);
3775 if (fVector2d && fPH2dOn && fWrite[1]) {
3780 TTree *treePH2d = ConvertVectorPTreeHisto(fVectorPH,fPlaPH,"treePH2d",(const Char_t *) name);
3781 fout->WriteTObject(treePH2d);
3790 //____________Pad Calibration Public___________________________________________
3792 //____________Define the number of pads per group for one detector and one calibration
3793 void AliTRDCalibra::ModePadCalibration(Int_t iChamb, Int_t i)
3796 // Definition of the calibration mode
3797 // from Nz and Nrphi, the number of row and col pads per calibration groups are setted
3804 if ((fNz[i] == 0) && (iChamb == 2)) {
3807 if ((fNz[i] == 0) && (iChamb != 2)) {
3810 if ((fNz[i] == 1) && (iChamb == 2)) {
3813 if ((fNz[i] == 1) && (iChamb != 2)) {
3816 if ((fNz[i] == 2) && (iChamb == 2)) {
3819 if ((fNz[i] == 2) && (iChamb != 2)) {
3829 if (fNrphi[i] == 0) {
3832 if (fNrphi[i] == 1) {
3835 if (fNrphi[i] == 2) {
3838 if (fNrphi[i] == 3) {
3841 if (fNrphi[i] == 4) {
3844 if (fNrphi[i] == 5) {
3847 if (fNrphi[i] == 6) {
3853 //____________Define the number of pad groups in one detector for one calibration
3854 Bool_t AliTRDCalibra::ModePadFragmentation(Int_t iPlane,Int_t iChamb, Int_t iSect, Int_t i)
3857 // Definition of the calibration mode
3858 // From the number of row and col pads per calibration groups the
3859 // number of calibration groups are setted
3865 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
3867 AliInfo("Could not get CommonParam Manager");
3871 // A little geometry:
3872 Int_t rowMax = parCom->GetRowMax(iPlane,iChamb,iSect);
3873 Int_t colMax = parCom->GetColMax(iPlane);
3875 // The fragmentation
3877 fNfragZ[i] = (Int_t) rowMax / fNnZ[i];
3880 if (fNnRphi[i] != 0) {
3881 fNfragRphi[i] = (Int_t) colMax / fNnRphi[i];
3888 //____________Protected Functions______________________________________________
3889 //____________Create the 2D histo to be filled online__________________________
3892 //_____________________________________________________________________________
3893 void AliTRDCalibra::CreatePRF2d(Int_t nn)
3896 // Create the 2D histos
3904 fPRF2d = new TProfile2D("PRF2d",(const Char_t *) name
3905 ,nn,0,nn,fNumberBinPRF,-1.0,1.0);
3906 fPRF2d->SetXTitle("Det/pad groups");
3907 fPRF2d->SetYTitle("Position x/W [pad width units]");
3908 fPRF2d->SetZTitle("Q_{i}/Q_{total}");
3909 fPRF2d->SetStats(0);
3913 //_____________________________________________________________________________
3914 void AliTRDCalibra::CreatePH2d(Int_t nn)
3917 // Create the 2D histos
3925 fPH2d = new TProfile2D("PH2d",(const Char_t *) name
3927 ,-0.5/fSf,(Float_t) (fTimeMax-0.5)/fSf);
3928 fPH2d->SetXTitle("Det/pad groups");
3929 fPH2d->SetYTitle("time [#mus]");
3930 fPH2d->SetZTitle("<PH> [a.u.]");
3935 //_____________________________________________________________________________
3936 void AliTRDCalibra::CreateCH2d(Int_t nn)
3939 // Create the 2D histos
3947 fCH2d = new TH2I("CH2d",(const Char_t *) name
3948 ,nn,0,nn,fNumberBinCharge,0,300);
3949 fCH2d->SetXTitle("Det/pad groups");
3950 fCH2d->SetYTitle("charge deposit [a.u]");
3951 fCH2d->SetZTitle("counts");
3957 //____________Offine tracking in the AliTRDtracker_____________________________
3958 void AliTRDCalibra::FillTheInfoOfTheTrackCH()
3961 // For the offline tracking or mcm tracklets
3962 // This function will be called in the functions UpdateHistogram...
3963 // to fill the info of a track for the relativ gain calibration
3966 Int_t nb = 0; // Nombre de zones traversees
3967 Int_t fd = -1; // Premiere zone non nulle
3970 // See if the track goes through different zones
3971 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3972 if (fAmpTotal[k] > 0.0) {
3980 // If automatic scale
3981 if ((fCountRelativeScale < 100) && (fRelativeScaleAuto)) {
3982 // Take only the one zone track
3984 fRelativeScale += fAmpTotal[fd] * 0.014 * 0.01;
3985 fCountRelativeScale++;
3989 // We fill the CH2d after having scale with the first 100
3990 if ((fCountRelativeScale >= 100) && (fRelativeScaleAuto)) {
3991 // Case of track with only one zone
3994 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
3997 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4000 // Case of track with two zones
4002 // Two zones voisines sinon rien!
4003 if ((fAmpTotal[fd] > 0.0) &&
4004 (fAmpTotal[fd+1] > 0.0)) {
4005 // One of the two very big
4006 if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+1]) {
4008 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
4011 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4014 if (fAmpTotal[fd+1] > fProcent*fAmpTotal[fd]) {
4016 fCH2d->Fill(fXbins[0]+fd+1.5,fAmpTotal[fd+1]/fRelativeScale);
4019 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd+1]/fRelativeScale);
4026 // Fill with no automatic scale
4027 if (!fRelativeScaleAuto) {
4028 // Case of track with only one zone
4032 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
4035 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4038 // Case of track with two zones
4040 // Two zones voisines sinon rien!
4042 if ((fAmpTotal[fd] > 0.0) &&
4043 (fAmpTotal[fd+1] > 0.0)) {
4044 // One of the two very big
4045 if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+1]) {
4047 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
4050 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4054 if (fAmpTotal[fd+1] > fProcent*fAmpTotal[fd]) {
4056 fCH2d->Fill(fXbins[0]+fd+1.5,fAmpTotal[fd+1]/fRelativeScale);
4059 UpdateVectorCH(fXbins[0]+fd+1,fAmpTotal[fd+1]/fRelativeScale);
4065 if (fNfragZ[0] > 1) {
4066 if (fAmpTotal[fd] > 0.0) {
4067 if ((fd+fNfragZ[0]) < (fNfragZ[0]*fNfragRphi[0])) {
4068 if (fAmpTotal[fd+fNfragZ[0]] > 0.0) {
4069 // One of the two very big
4070 if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+fNfragZ[0]]) {
4072 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
4075 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4079 if (fAmpTotal[fd+fNfragZ[0]] > fProcent*fAmpTotal[fd]) {
4081 fCH2d->Fill(fXbins[0]+fd+fNfragZ[0]+0.5,fAmpTotal[fd+fNfragZ[0]]/fRelativeScale);
4085 UpdateVectorCH(fXbins[0]+fd+fNfragZ[0],fAmpTotal[fd+fNfragZ[0]]/fRelativeScale);
4098 //____________Offine tracking in the AliTRDtracker_____________________________
4099 void AliTRDCalibra::ResetfVariables()
4102 // Reset values of fAmpTotal, fPHValue and fPHPlace for
4103 // the updateHistogram... functions
4106 // Reset the good track
4109 // Reset the fAmpTotal where we put value
4111 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
4116 // Reset the fPHValue
4118 for (Int_t k = 0; k < fTimeMax; k++) {
4126 //____________Offine tracking in the AliTRDtracker_____________________________
4127 void AliTRDCalibra::FillTheInfoOfTheTrackPH()
4130 // For the offline tracking or mcm tracklets
4131 // This function will be called in the functions UpdateHistogram...
4132 // to fill the info of a track for the drift velocity calibration
4135 Int_t nb = 1; // Nombre de zones traversees 1, 2 ou plus de 3
4136 Int_t fd1 = -1; // Premiere zone non nulle
4137 Int_t fd2 = -1; // Deuxieme zone non nulle
4138 Int_t k1 = -1; // Debut de la premiere zone
4139 Int_t k2 = -1; // Debut de la seconde zone
4141 // See if the track goes through different zones
4142 for (Int_t k = 0; k < fTimeMax; k++) {
4143 if (fPHValue[k] > 0.0) {
4148 if (fPHPlace[k] != fd1) {
4154 if (fPHPlace[k] != fd2) {
4162 // Case of track with only one zone
4165 for (Int_t i = 0; i < fTimeMax; i++) {
4166 if (fPHValue[i] > 0.0) {
4168 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4171 AliInfo(Form("WRITE nb %d ,place final: %d, fPHPlace[i]: %d, i: %d, fPHValue[i]: %f"
4172 ,nb,fXbins[1]+fPHPlace[i],fPHPlace[i],i,fPHValue[i]));
4175 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4180 // Case of track with two zones
4182 // Two zones voisines sinon rien!
4184 if ((fd1 == fd2+1) ||
4186 // One of the two fast all the think
4187 if (k2 > (k1+fDifference)) {
4189 for (Int_t i = k1; i < k2; i++) {
4190 if (fPHValue[i] > 0.0) {
4192 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4195 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4200 if ((k2+fDifference) < fTimeMax) {
4202 for (Int_t i = k2; i < fTimeMax; i++) {
4203 if (fPHValue[i] > 0.0) {
4205 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4208 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4214 // Two zones voisines sinon rien!
4215 if (fNfragZ[1] > 1) {
4217 if ((fd1+fNfragZ[1]) < (fNfragZ[1]*fNfragRphi[1])) {
4218 if (fd2 == (fd1+fNfragZ[1])) {
4219 // One of the two fast all the think
4220 if (k2 > (k1+fDifference)) {
4222 for (Int_t i = k1; i < k2; i++) {
4223 if (fPHValue[i] > 0.0) {
4225 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4228 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4233 if ((k2+fDifference) < fTimeMax) {
4235 for (Int_t i = k2; i < fTimeMax; i++) {
4236 if (fPHValue[i] > 0.0) {
4238 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4241 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4248 // Two zones voisines sinon rien!
4250 if ((fd1 - fNfragZ[1]) >= 0) {
4251 if (fd2 == (fd1 - fNfragZ[1])) {
4252 // One of the two fast all the think
4253 if (k2 > (k1 + fDifference)) {
4255 for (Int_t i = k1; i < k2; i++) {
4256 if (fPHValue[i] > 0.0) {
4258 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4261 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4266 if ((k2+fDifference) < fTimeMax) {
4268 for (Int_t i = k2; i < fTimeMax; i++) {
4269 if (fPHValue[i] > 0.0) {
4271 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4274 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4287 //____________Set the pad calibration variables for the detector_______________
4288 Bool_t AliTRDCalibra::LocalisationDetectorXbins(Int_t detector)
4291 // For the detector calcul the first Xbins and set the number of row
4292 // and col pads per calibration groups, the number of calibration
4293 // groups in the detector.
4296 // first Xbins of the detector
4298 CalculXBins(detector,0);
4301 CalculXBins(detector,1);
4304 CalculXBins(detector,2);
4307 // fragmentation of idect
4308 for (Int_t i = 0; i < 3; i++) {
4309 ModePadCalibration((Int_t) GetChamber(detector),i);
4310 ModePadFragmentation((Int_t) GetPlane(detector)
4311 , (Int_t) GetChamber(detector)
4312 , (Int_t) GetSector(detector),i);
4319 //____________Plot the 2D histos filled Online_________________________________
4321 //_____________________________________________________________________________
4322 void AliTRDCalibra::PlotPH2d()
4325 // Plot the 2D histo
4328 TCanvas *cph2d = new TCanvas("cph2d","",50,50,600,800);
4330 fPH2d->Draw("LEGO");
4334 //_____________________________________________________________________________
4335 void AliTRDCalibra::PlotCH2d()
4338 // Plot the 2D histos
4341 TCanvas *cch2d = new TCanvas("cch2d","",50,50,600,800);
4343 fCH2d->Draw("LEGO");
4347 //_____________________________________________________________________________
4348 void AliTRDCalibra::PlotPRF2d()
4351 // Plot the 2D histos
4354 TCanvas *cPRF2d = new TCanvas("cPRF2d","",50,50,600,800);
4356 fPRF2d->Draw("LEGO");
4360 //____________Fit______________________________________________________________
4362 //____________Create histos if fDebug == 1 or fDebug >= 3______________________
4364 //_____________________________________________________________________________
4365 void AliTRDCalibra::CreateFitHistoPH(Int_t nbins, Double_t low, Double_t high)
4368 // Create the histos for fDebug = 1 and fDebug = 4 (Fit functions)
4371 // Histograms to store the coef
4372 fCoefVdrift[0] = new TH1F("coefvdrift0" ,"",nbins,low ,high);
4373 fCoefVdrift[1] = new TH1F("coefvdrift1" ,"",nbins,low ,high);
4374 fCoefVdrift[2] = new TH1F("coefvdrift2" ,"",nbins,low ,high);
4376 // Histograms for Debug
4377 fDeltaVdrift[0] = new TH1F("deltavdrift0","",nbins,low ,high);
4378 fDeltaVdrift[1] = new TH1F("deltavdrift1","",nbins,low ,high);
4379 fErrorVdrift[0] = new TH1I("errorvdrift0","",300 ,-0.5,0.5);
4380 fErrorVdrift[1] = new TH1I("errorvdrift1","",300 ,-0.5,0.5);
4382 fCoefVdrift[0]->SetXTitle("Det/pad groups");
4383 fCoefVdrift[0]->SetYTitle("Vdrift [cm/#mus]");
4384 fCoefVdrift[1]->SetXTitle("Det/pad groups");
4385 fCoefVdrift[1]->SetYTitle("Vdrift [cm/#mus]");
4386 fCoefVdrift[2]->SetXTitle("Det/pad groups");
4387 fCoefVdrift[2]->SetYTitle("Vdrift [cm/#mus]");
4389 fDeltaVdrift[0]->SetXTitle("Det/pad groups");
4390 fDeltaVdrift[0]->SetYTitle("#Deltav/v_{sim}");
4391 fDeltaVdrift[1]->SetXTitle("Det/pad groups");
4392 fDeltaVdrift[1]->SetYTitle("#Deltav/v_{sim}");
4394 fErrorVdrift[0]->SetXTitle("#Deltav/v_{sim}");
4395 fErrorVdrift[0]->SetYTitle("counts");
4396 fErrorVdrift[1]->SetXTitle("#Deltav/v_{sim}");
4397 fErrorVdrift[1]->SetYTitle("counts");
4399 fCoefVdrift[0]->SetStats(0);
4400 fCoefVdrift[1]->SetStats(0);
4401 fCoefVdrift[2]->SetStats(0);
4402 fDeltaVdrift[0]->SetStats(0);
4403 fDeltaVdrift[1]->SetStats(0);
4404 fErrorVdrift[0]->SetStats(0);
4405 fErrorVdrift[1]->SetStats(0);
4407 fCoefVdrift[0]->SetMarkerColor(6);
4408 fCoefVdrift[0]->SetMarkerStyle(26);
4409 fCoefVdrift[0]->SetLineColor(6);
4410 fCoefVdrift[1]->SetMarkerColor(2);
4411 fCoefVdrift[1]->SetMarkerStyle(24);
4412 fCoefVdrift[1]->SetLineColor(2);
4413 fCoefVdrift[2]->SetLineColor(4);
4415 fDeltaVdrift[1]->SetMarkerColor(2);
4416 fDeltaVdrift[1]->SetMarkerStyle(24);
4417 fDeltaVdrift[1]->SetLineColor(2);
4418 fDeltaVdrift[0]->SetMarkerColor(6);
4419 fDeltaVdrift[0]->SetMarkerStyle(26);
4420 fDeltaVdrift[0]->SetLineColor(6);
4422 fErrorVdrift[1]->SetLineColor(2);
4423 fErrorVdrift[1]->SetLineStyle(2);
4424 fErrorVdrift[0]->SetLineColor(6);
4425 fErrorVdrift[0]->SetLineStyle(1);
4429 //_____________________________________________________________________________
4430 void AliTRDCalibra::CreateFitHistoT0(Int_t nbins, Double_t low, Double_t high)
4433 // Create the histos for fDebug = 1 and fDebug = 4 (Fit functions)
4436 // Histograms to store the coef
4437 fCoefT0[0] = new TH1F("coefT00" ,"",nbins,low ,high);
4438 fCoefT0[1] = new TH1F("coefT01" ,"",nbins,low ,high);
4439 fCoefT0[2] = new TH1F("coefT02" ,"",nbins,low ,high);
4441 // Histograms for Debug
4442 fDeltaT0[0] = new TH1F("deltaT00","",nbins,low ,high);
4443 fDeltaT0[1] = new TH1F("deltaT01","",nbins,low ,high);
4444 fErrorT0[0] = new TH1I("errorT00","",100,-0.2,0.2);
4445 fErrorT0[1] = new TH1I("errorT01","",100,-0.2,0.2);
4447 fCoefT0[0]->SetXTitle("Det/pad groups");
4448 fCoefT0[0]->SetYTitle("t0 [timebin]");
4449 fCoefT0[1]->SetXTitle("Det/pad groups");
4450 fCoefT0[1]->SetYTitle("t0 [timebin]");
4451 fCoefT0[2]->SetXTitle("Det/pad groups");
4452 fCoefT0[2]->SetYTitle("t0 [timebin]");
4454 fDeltaT0[0]->SetXTitle("Det/pad groups");
4455 fDeltaT0[0]->SetYTitle("#Deltat0 [timebin]");
4456 fDeltaT0[1]->SetXTitle("Det/pad groups");
4457 fDeltaT0[1]->SetYTitle("#Deltat0 [timebin]");
4459 fErrorT0[0]->SetXTitle("#Deltat0 [timebin]");
4460 fErrorT0[0]->SetYTitle("counts");
4461 fErrorT0[1]->SetXTitle("#Deltat0 [timebin]");
4462 fErrorT0[1]->SetYTitle("counts");
4464 fCoefT0[0]->SetStats(0);
4465 fCoefT0[1]->SetStats(0);
4466 fCoefT0[2]->SetStats(0);
4467 fDeltaT0[0]->SetStats(0);
4468 fDeltaT0[1]->SetStats(0);
4469 fErrorT0[0]->SetStats(0);
4470 fErrorT0[1]->SetStats(0);
4472 fCoefT0[0]->SetMarkerColor(6);
4473 fCoefT0[0]->SetMarkerStyle(26);
4474 fCoefT0[0]->SetLineColor(6);
4475 fCoefT0[1]->SetMarkerColor(2);
4476 fCoefT0[1]->SetMarkerStyle(24);
4477 fCoefT0[1]->SetLineColor(2);
4478 fCoefT0[2]->SetLineColor(4);
4480 fDeltaT0[1]->SetMarkerColor(2);
4481 fDeltaT0[1]->SetMarkerStyle(24);
4482 fDeltaT0[1]->SetLineColor(2);
4483 fDeltaT0[0]->SetMarkerColor(6);
4484 fDeltaT0[0]->SetMarkerStyle(26);
4485 fDeltaT0[0]->SetLineColor(6);
4487 fErrorT0[1]->SetLineColor(2);
4488 fErrorT0[1]->SetLineStyle(2);
4489 fErrorT0[0]->SetLineColor(6);
4490 fErrorT0[0]->SetLineStyle(1);
4494 //_____________________________________________________________________________
4495 void AliTRDCalibra::CreateFitHistoCH(Int_t nbins, Double_t low, Double_t high)
4498 // Create the histos for fDebug = 1 and fDebug = 4 (Fit functions)
4501 // Histograms to store the coef
4502 fCoefCharge[0] = new TH1F("coefcharge0" ,"",nbins,low ,high);
4503 fCoefCharge[1] = new TH1F("coefcharge1" ,"",nbins,low ,high);
4504 fCoefCharge[2] = new TH1F("coefcharge2" ,"",nbins,low ,high);
4505 fCoefCharge[3] = new TH1F("coefcharge3" ,"",nbins,low ,high);
4507 // Histograms for Debug
4508 fDeltaCharge[0] = new TH1F("deltacharge0","",nbins,low ,high);
4509 fDeltaCharge[1] = new TH1F("deltacharge1","",nbins,low ,high);
4510 fDeltaCharge[2] = new TH1F("deltacharge2","",nbins,low ,high);
4512 fErrorCharge[0] = new TH1I("errorcharge0","",100 ,-0.5,0.5);
4513 fErrorCharge[1] = new TH1I("errorcharge1","",100 ,-0.5,0.5);
4514 fErrorCharge[2] = new TH1I("errorcharge2","",100 ,-0.5,0.5);
4516 fCoefCharge[0]->SetXTitle("Det/Pad groups");
4517 fCoefCharge[0]->SetYTitle("gain factor");
4518 fCoefCharge[1]->SetXTitle("Det/Pad groups");
4519 fCoefCharge[1]->SetYTitle("gain factor");
4520 fCoefCharge[2]->SetXTitle("Det/Pad groups");
4521 fCoefCharge[2]->SetYTitle("gain factor");
4522 fCoefCharge[3]->SetXTitle("Det/Pad groups");
4523 fCoefCharge[3]->SetYTitle("gain factor");
4525 fDeltaCharge[0]->SetXTitle("Det/Pad groups");
4526 fDeltaCharge[0]->SetYTitle("#Deltag/g_{sim}");
4527 fDeltaCharge[1]->SetXTitle("Det/Pad groups");
4528 fDeltaCharge[1]->SetYTitle("#Deltag/g_{sim}");
4529 fDeltaCharge[2]->SetXTitle("Det/Pad groups");
4530 fDeltaCharge[2]->SetYTitle("#Deltag/g_{sim}");
4531 fDeltaCharge[0]->SetAxisRange(-0.5,0.5,"Y");
4532 fDeltaCharge[1]->SetAxisRange(-0.5,0.5,"Y");
4533 fDeltaCharge[2]->SetAxisRange(-0.5,0.5,"Y");
4535 fErrorCharge[0]->SetXTitle("#Deltag/g_{sim}");
4536 fErrorCharge[0]->SetYTitle("counts");
4537 fErrorCharge[1]->SetXTitle("#Deltag/g_{sim}");
4538 fErrorCharge[1]->SetYTitle("counts");
4539 fErrorCharge[2]->SetXTitle("#Deltag/g_{sim}");
4540 fErrorCharge[2]->SetYTitle("counts");
4542 fDeltaCharge[1]->SetMarkerColor(2);
4543 fDeltaCharge[1]->SetMarkerStyle(24);
4544 fDeltaCharge[1]->SetLineColor(2);
4545 fErrorCharge[1]->SetLineColor(2);
4546 fErrorCharge[1]->SetLineStyle(2);
4547 fDeltaCharge[2]->SetMarkerColor(8);
4548 fDeltaCharge[2]->SetLineColor(8);
4549 fDeltaCharge[2]->SetMarkerStyle(9);
4550 fErrorCharge[2]->SetLineColor(8);
4551 fErrorCharge[2]->SetLineStyle(5);
4552 fDeltaCharge[0]->SetMarkerColor(6);
4553 fDeltaCharge[0]->SetLineColor(6);
4554 fDeltaCharge[0]->SetMarkerStyle(26);
4555 fErrorCharge[0]->SetLineColor(6);
4556 fErrorCharge[0]->SetLineStyle(1);
4558 fCoefCharge[3]->SetLineColor(4);
4559 fCoefCharge[1]->SetMarkerColor(2);
4560 fCoefCharge[1]->SetLineColor(2);
4561 fCoefCharge[1]->SetMarkerStyle(24);
4562 fCoefCharge[2]->SetMarkerColor(8);
4563 fCoefCharge[2]->SetLineColor(8);
4564 fCoefCharge[2]->SetMarkerStyle(9);
4565 fCoefCharge[0]->SetMarkerColor(6);
4566 fCoefCharge[0]->SetLineColor(6);
4567 fCoefCharge[0]->SetMarkerStyle(26);
4569 fErrorCharge[2]->SetLineWidth(3);
4571 fDeltaCharge[1]->SetStats(0);
4572 fDeltaCharge[2]->SetStats(0);
4573 fDeltaCharge[0]->SetStats(0);
4574 fErrorCharge[1]->SetStats(0);
4575 fErrorCharge[2]->SetStats(0);
4576 fErrorCharge[0]->SetStats(0);
4577 fCoefCharge[1]->SetStats(0);
4578 fCoefCharge[0]->SetStats(0);
4579 fCoefCharge[3]->SetStats(0);
4580 fCoefCharge[2]->SetStats(0);
4584 //_____________________________________________________________________________
4585 void AliTRDCalibra::CreateFitHistoPRF(Int_t nbins, Double_t low, Double_t high)
4588 // Create the histos for fDebug = 1 and fDebug = 4 (Fit functions)
4591 // Histograms to store the coef
4592 fCoefPRF[0] = new TH1F("coefPRF0","",nbins,low ,high);
4593 fCoefPRF[1] = new TH1F("coefPRF1","",nbins,low ,high);
4595 // Histograms for Debug
4596 fDeltaPRF = new TH1F("deltaPRF","",nbins,low ,high);
4597 fErrorPRF = new TH1I("errorPRF","",300,-0.5,0.5);
4599 fDeltaPRF->SetMarkerColor(6);
4600 fDeltaPRF->SetMarkerStyle(26);
4601 fDeltaPRF->SetLineColor(6);
4602 fErrorPRF->SetLineColor(6);
4603 fErrorPRF->SetLineStyle(2);
4605 fCoefPRF[1]->SetLineColor(4);
4606 fCoefPRF[0]->SetMarkerColor(6);
4607 fCoefPRF[0]->SetMarkerStyle(26);
4608 fCoefPRF[0]->SetLineColor(6);
4610 fCoefPRF[0]->SetXTitle("Det/Pad groups");
4611 fCoefPRF[0]->SetYTitle("#sigma_{PRF}");
4612 fCoefPRF[1]->SetXTitle("Det/Pad groups");
4613 fCoefPRF[1]->SetYTitle("#sigma_{PRF}");
4615 fDeltaPRF->SetXTitle("Det/Pad groups");
4616 fDeltaPRF->SetYTitle("#Delta#sigma/#sigma_{sim}");
4618 fErrorPRF->SetXTitle("#Delta#sigma/#sigma_{sim}");
4619 fErrorPRF->SetYTitle("counts");
4621 fDeltaPRF->SetStats(0);
4622 fErrorPRF->SetStats(0);
4623 fCoefPRF[1]->SetStats(0);
4624 fCoefPRF[0]->SetStats(0);
4628 //_____________________________________________________________________________
4629 void AliTRDCalibra::CreateFitHistoPRFDB(Int_t rowMax, Int_t colMax)
4632 // Create the histos for fDebug = 3 and fDebug = 4 (Fit functions)
4635 fCoefPRFDB = new TH2F("coefPRF","",rowMax,0,rowMax,colMax,0,colMax);
4637 fCoefPRFDB->SetStats(0);
4638 fCoefPRFDB->SetXTitle("row Number");
4639 fCoefPRFDB->SetYTitle("col Number");
4640 fCoefPRFDB->SetZTitle("PRF width [pad width units]");
4642 fCoefPRFDB->SetFillColor(6);
4643 fCoefPRFDB->SetLineColor(6);
4647 //_____________________________________________________________________________
4648 void AliTRDCalibra::CreateFitHistoCHDB(Int_t rowMax, Int_t colMax)
4651 // Create the histos for fDebug = 3 and fDebug = 4 (Fit functions)
4654 fCoefChargeDB[0] = new TH2F("coefchargedb0","",rowMax,0,rowMax,colMax,0,colMax);
4655 fCoefChargeDB[1] = new TH2F("coefchargedb1","",rowMax,0,rowMax,colMax,0,colMax);
4656 fCoefChargeDB[2] = new TH2F("coefchargedb2","",rowMax,0,rowMax,colMax,0,colMax);
4658 fCoefChargeDB[0]->SetStats(0);
4659 fCoefChargeDB[1]->SetStats(0);
4660 fCoefChargeDB[2]->SetStats(0);
4661 fCoefChargeDB[0]->SetXTitle("row Number");
4662 fCoefChargeDB[0]->SetYTitle("col Number");
4663 fCoefChargeDB[1]->SetXTitle("row Number");
4664 fCoefChargeDB[1]->SetYTitle("col Number");
4665 fCoefChargeDB[2]->SetXTitle("row Number");
4666 fCoefChargeDB[2]->SetYTitle("col Number");
4667 fCoefChargeDB[0]->SetZTitle("f_{g} Fit method");
4668 fCoefChargeDB[1]->SetZTitle("f_{g} Mean method");
4669 fCoefChargeDB[2]->SetZTitle("f_{g} Fitbis method");
4671 fCoefChargeDB[0]->SetFillColor(6);
4672 fCoefChargeDB[0]->SetLineColor(6);
4673 fCoefChargeDB[0]->SetLineColor(6);
4674 fCoefChargeDB[1]->SetFillColor(2);
4675 fCoefChargeDB[1]->SetLineColor(2);
4676 fCoefChargeDB[1]->SetLineColor(2);
4677 fCoefChargeDB[2]->SetFillColor(8);
4678 fCoefChargeDB[2]->SetLineColor(8);
4679 fCoefChargeDB[2]->SetLineColor(8);
4683 //_____________________________________________________________________________
4684 void AliTRDCalibra::CreateFitHistoPHDB(Int_t rowMax, Int_t colMax)
4687 // Create the histos for fDebug = 3 and fDebug = 4 (Fit functions)
4690 fCoefVdriftDB[0] = new TH2F("coefvdriftdb0","",rowMax,0,rowMax,colMax,0,colMax);
4691 fCoefVdriftDB[1] = new TH2F("coefvdriftdb1","",rowMax,0,rowMax,colMax,0,colMax);
4693 fCoefVdriftDB[0]->SetStats(0);
4694 fCoefVdriftDB[1]->SetStats(0);
4695 fCoefVdriftDB[0]->SetXTitle("row Number");
4696 fCoefVdriftDB[0]->SetYTitle("col Number");
4697 fCoefVdriftDB[1]->SetXTitle("row Number");
4698 fCoefVdriftDB[1]->SetYTitle("col Number");
4699 fCoefVdriftDB[0]->SetZTitle("v_{drift} Fit method");
4700 fCoefVdriftDB[1]->SetZTitle("v_{drift} slope method");
4702 fCoefVdriftDB[0]->SetFillColor(6);
4703 fCoefVdriftDB[0]->SetLineColor(6);
4704 fCoefVdriftDB[0]->SetLineColor(6);
4705 fCoefVdriftDB[1]->SetFillColor(2);
4706 fCoefVdriftDB[1]->SetLineColor(2);
4707 fCoefVdriftDB[1]->SetLineColor(2);
4711 //_____________________________________________________________________________
4712 void AliTRDCalibra::CreateFitHistoT0DB(Int_t rowMax, Int_t colMax)
4715 // Create the histos for fDebug = 3 and fDebug = 4 (Fit functions)
4718 fCoefT0DB[0] = new TH2F("coefT0db0","",rowMax,0,rowMax,colMax,0,colMax);
4719 fCoefT0DB[1] = new TH2F("coefT0db1","",rowMax,0,rowMax,colMax,0,colMax);
4721 fCoefT0DB[0]->SetStats(0);
4722 fCoefT0DB[1]->SetStats(0);
4723 fCoefT0DB[0]->SetXTitle("row Number");
4724 fCoefT0DB[0]->SetYTitle("col Number");
4725 fCoefT0DB[1]->SetXTitle("row Number");
4726 fCoefT0DB[1]->SetYTitle("col Number");
4727 fCoefT0DB[0]->SetZTitle("t0 Fit method");
4728 fCoefT0DB[1]->SetZTitle("t0 slope method");
4730 fCoefT0DB[0]->SetFillColor(6);
4731 fCoefT0DB[0]->SetLineColor(6);
4732 fCoefT0DB[0]->SetLineColor(6);
4733 fCoefT0DB[1]->SetFillColor(2);
4734 fCoefT0DB[1]->SetLineColor(2);
4735 fCoefT0DB[1]->SetLineColor(2);
4739 //_____________________________________________________________________________
4740 Bool_t AliTRDCalibra::FillVectorFitCH(Int_t countdet)
4743 // For the Fit functions fill the vector FitCH special for the gain calibration
4746 AliTRDFitCHInfo *fitCHInfo = new AliTRDFitCHInfo();
4749 if (GetChamber(countdet) == 2) {
4756 Float_t *coef = new Float_t[ntotal];
4757 for (Int_t i = 0; i < ntotal; i++) {
4758 coef[i] = fCoefCH[i];
4761 Int_t detector = countdet;
4763 fitCHInfo->SetCoef(coef);
4764 fitCHInfo->SetDetector(detector);
4765 fVectorFitCH->Add((TObject *) fitCHInfo);
4771 //____________Functions for initialising the AliTRDCalibra in the code_________
4772 Bool_t AliTRDCalibra::InitFit(Int_t nbins, Double_t lowedge
4773 , Double_t upedge, Int_t i)
4776 // Init the calibration mode (Nz, Nrphi), the histograms for
4777 // debugging the fit methods if fDebug > 0,
4780 gStyle->SetPalette(1);
4781 gStyle->SetOptStat(1111);
4782 gStyle->SetPadBorderMode(0);
4783 gStyle->SetCanvasColor(10);
4784 gStyle->SetPadLeftMargin(0.13);
4785 gStyle->SetPadRightMargin(0.01);
4787 // Get the parameter object
4788 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
4790 AliInfo("Could not get CommonParam");
4794 // Mode groups of pads: the total number of bins!
4795 Int_t numberofbinsexpected = 0;
4796 ModePadCalibration(2,i);
4797 ModePadFragmentation(0,2,0,i);
4798 fDetChamb2[i] = fNfragZ[i] * fNfragRphi[i];
4800 AliInfo(Form("For the chamber 2: %d",fDetChamb2[i]));
4802 numberofbinsexpected += 6 * 18 * fDetChamb2[i];
4803 ModePadCalibration(0,i);
4804 ModePadFragmentation(0,0,0,i);
4805 fDetChamb0[i] = fNfragZ[i] * fNfragRphi[i];
4807 AliInfo(Form("For the other chamber 0: %d",fDetChamb0[i]));
4809 numberofbinsexpected += 6 * 4 * 18 * fDetChamb0[i];
4811 // Quick verification that we have the good pad calibration mode if 2D histos!
4813 if (numberofbinsexpected != nbins) {
4814 AliInfo("It doesn't correspond to the mode of pad group calibration!");
4819 // Security for fDebug 3 and 4
4820 if ((fDebug >= 3) &&
4824 AliInfo("This detector doesn't exit!");
4828 // Determine fDet1 and fDet2
4832 fDect1[i] = fFitVoir;
4833 fDect2[i] = fDect1[i] +1;
4837 fDect2[i] = numberofbinsexpected;
4840 CalculXBins(AliTRDgeometry::GetDetector(fDet[0],fDet[1],fDet[2]),i);
4841 fDect1[i] = fXbins[i];
4842 CalculXBins((AliTRDgeometry::GetDetector(fDet[0],fDet[1],fDet[2])+1),i);
4843 fDect2[i] = fXbins[i];
4846 // Create the histos for debugging
4851 // Init the VectorFitCH
4852 fVectorFitCH = new TObjArray();
4853 fCoefCH = new Float_t[2304];
4854 for (Int_t k = 0; k < 2304; k++) {
4857 fScaleFitFactor = 0.0;
4859 // Number of Xbins(detectors or groups of pads) if Vector2d
4860 // Quick verification that we are not out of range!
4861 if (fVectorCH && fPlaCH) {
4863 (fVectorCH->GetEntriesFast() > 0) &&
4864 ((Int_t) fPlaCH->GetEntriesFast() > 0)) {
4865 if ((Int_t) fVectorCH->GetEntriesFast() > numberofbinsexpected) {
4866 AliInfo("ch doesn't correspond to the mode of pad group calibration!");
4869 if ((Int_t) fVectorCH->GetEntriesFast() != (Int_t) fPlaCH->GetEntriesFast()) {
4870 AliInfo("VectorCH doesn't correspond to PlaCH!");
4877 // Debugging: Create the histos
4880 // fDebug == 0 nothing
4885 // Create the histos replique de ch if histos2D
4886 CreateFitHistoCH(nbins,lowedge,upedge);
4889 // Ccreate the histos replique de ch vector2d
4890 CreateFitHistoCH(numberofbinsexpected,0,numberofbinsexpected);
4894 // fDebug == 2 and fFitVoir no histo
4896 if (fFitVoir < numberofbinsexpected) {
4897 AliInfo(Form("We will see the fit of the object %d",fFitVoir));
4900 AliInfo("fFitVoir is out of range of the histo!");
4905 // fDebug == 3 or 4 and fDet
4907 if ((fNz[0] == 0) && (fNrphi[0] == 0)) {
4908 AliInfo("Do you really want to see one detector without pad groups?");
4912 AliInfo(Form("You will see the detector: iPlane %d, iChamb %d, iSect %d"
4913 ,fDet[0],fDet[1],fDet[2]));
4914 // A little geometry:
4915 Int_t rowMax = parCom->GetRowMax(fDet[0],fDet[1],fDet[2]);
4916 Int_t colMax = parCom->GetColMax(fDet[0]);
4917 // Create the histos to visualise
4918 CreateFitHistoCHDB(rowMax,colMax);
4920 CreateFitHistoCH((Int_t) (fDect2[0]-fDect1[0]),fDect1[0],fDect2[0]);
4930 // Number of Xbins (detectors or groups of pads) if vector2d
4931 // Quick verification that we are not out of range!
4932 if (fVectorPH && fPlaPH) {
4934 (fVectorPH->GetEntriesFast() > 0) &&
4935 ((Int_t) fPlaPH->GetEntriesFast() > 0)) {
4936 if ((Int_t) fVectorPH->GetEntriesFast() > numberofbinsexpected) {
4937 AliInfo("ph doesn't correspond to the mode of pad group calibration!");
4940 if ((Int_t) fVectorPH->GetEntriesFast() != (Int_t) fPlaPH->GetEntriesFast()) {
4941 AliInfo("VectorPH doesn't correspond to PlaPH!");
4952 // Debugging: Create the histos
4955 // fDebug == 0 nothing
4960 // Create the histos replique de ch
4961 CreateFitHistoPH(nbins,lowedge,upedge);
4962 CreateFitHistoT0(nbins,lowedge,upedge);
4965 // Create the histos replique de ch if vector2d
4966 CreateFitHistoPH(numberofbinsexpected,0,numberofbinsexpected);
4967 CreateFitHistoT0(numberofbinsexpected,0,numberofbinsexpected);
4971 // fDebug == 2 and fFitVoir no histo
4973 if (fFitVoir < numberofbinsexpected) {
4974 AliInfo(Form("We will see the fit of the object %d",fFitVoir));
4977 AliInfo("fFitVoir is out of range of the histo!");
4982 // fDebug == 3 or 4 and fDet
4984 if ((fNz[1] == 0) &&
4986 AliInfo("Do you really want to see one detector without pad groups?");
4990 AliInfo(Form("You will see the detector: iPlane %d, iChamb %d, iSect %d"
4991 ,fDet[0],fDet[1],fDet[2]));
4992 // A little geometry:
4993 Int_t rowMax = parCom->GetRowMax(fDet[0],fDet[1],fDet[2]);
4994 Int_t colMax = parCom->GetColMax(fDet[0]);
4995 // Create the histos to visualise
4996 CreateFitHistoPHDB(rowMax,colMax);
4997 CreateFitHistoT0DB(rowMax,colMax);
4999 CreateFitHistoPH((Int_t) (fDect2[1]-fDect1[1]),fDect1[1],fDect2[1]);
5000 CreateFitHistoT0((Int_t) (fDect2[1]-fDect1[1]),fDect1[1],fDect2[1]);
5010 // Number of Xbins(detectors or groups of pads) if vector2d
5011 if (fVectorPRF && fPlaPRF){
5013 (fVectorPRF->GetEntriesFast() > 0) &&
5014 (fPlaPRF->GetEntriesFast() > 0)) {
5015 // Quick verification that we are not out of range!
5016 if ((Int_t) fVectorPRF->GetEntriesFast() > numberofbinsexpected) {
5017 AliInfo("ch doesn't correspond to the mode of pad group calibration!");
5020 if ((Int_t) fVectorPRF->GetEntriesFast() != (Int_t) fPlaPRF->GetEntriesFast()) {
5021 AliInfo("VectorPRF doesn't correspond to PlaCH!");
5031 // Debugging: Create the histos
5034 // fDebug == 0 nothing
5039 // Create the histos replique de ch
5040 CreateFitHistoPRF(nbins,lowedge,upedge);
5043 // Create the histos replique de ch
5044 CreateFitHistoPRF(numberofbinsexpected,0,numberofbinsexpected);
5048 // fDebug == 2 and fFitVoir no histo
5050 if (fFitVoir < numberofbinsexpected) {
5051 AliInfo(Form("We will see the fit of the object %d",fFitVoir));
5054 AliInfo("fFitVoir is out of range of the histo!");
5059 // fDebug == 3 or 4 and fDet
5061 if ((fNz[2] == 0) &&
5063 AliInfo("Do you really want to see one detector without pad groups?");
5067 AliInfo(Form("You will see the detector: iPlane %d, iChamb %d, iSect %d"
5068 ,fDet[0],fDet[1],fDet[2]));
5069 // A little geometry:
5070 Int_t rowMax = parCom->GetRowMax(fDet[0],fDet[1],fDet[2]);
5071 Int_t colMax = parCom->GetColMax(fDet[0]);
5072 // Create the histos to visualise
5073 CreateFitHistoPRFDB(rowMax,colMax);
5075 CreateFitHistoPRF((Int_t) (fDect2[2]-fDect1[2]),fDect1[2],fDect2[2]);
5086 //____________Functions for initialising the AliTRDCalibra in the code_________
5087 void AliTRDCalibra::InitfCountDetAndfCount(Int_t i)
5090 // Init the current detector where we are fCountDet and the
5091 // next fCount for the functions Fit...
5094 // Loop on the Xbins of ch!!
5095 fCountDet[i] = -1; // Current detector
5096 fCount[i] = 0; // To find the next detector
5101 // Set countdet to the detector
5102 fCountDet[i] = AliTRDgeometry::GetDetector(fDet[0],fDet[1],fDet[2]);
5104 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
5105 ModePadCalibration(fDet[1],i);
5106 ModePadFragmentation(fDet[0],fDet[1],fDet[2],i);
5108 // Set counter to write at the end of the detector
5109 fCount[i] = fDect1[i] + fNfragZ[i]*fNfragRphi[i];
5115 //____________Functions for initialising the AliTRDCalibra in the code_________
5116 void AliTRDCalibra::UpdatefCountDetAndfCount(Int_t idect, Int_t i)
5119 // See if we are in a new detector and update the
5120 // variables fNfragZ and fNfragRphi if yes
5123 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
5124 // If fDebug == 1 or 0
5125 if ((fDebug == 0) ||
5128 if (fCount[i] == idect) {
5130 // On en est au detector
5133 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
5134 ModePadCalibration((Int_t) GetChamber(fCountDet[i]),i);
5135 ModePadFragmentation((Int_t) GetPlane(fCountDet[i])
5136 ,(Int_t) GetChamber(fCountDet[i])
5137 ,(Int_t) GetSector(fCountDet[i]),i);
5139 // Set for the next detector
5140 fCount[i] += fNfragZ[i]*fNfragRphi[i];
5148 //____________Functions for initialising the AliTRDCalibra in the code_________
5149 void AliTRDCalibra::ReconstructFitRowMinRowMax(Int_t idect, Int_t i)
5152 // Reconstruct the min pad row, max pad row, min pad col and
5153 // max pad col of the calibration group for the Fit functions
5157 ReconstructionRowPadGroup((Int_t) (idect-(fCount[i]-(fNfragZ[i]*fNfragRphi[i]))),i);
5160 ReconstructionRowPadGroup((Int_t) (idect-fDect1[i]),i);
5165 //____________Functions for initialising the AliTRDCalibra in the code_________
5166 Bool_t AliTRDCalibra::NotEnoughStatistic(Int_t idect, Int_t i)
5169 // For the case where there are not enough entries in the histograms
5170 // of the calibration group, the value present in the choosen database
5171 // will be put. A negativ sign enables to know that a fit was not possible.
5174 // Get the parameter object
5175 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
5177 AliInfo("Could not get CommonParam Manager");
5182 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
5184 AliInfo("Could not get calibDB");
5189 AliInfo(Form("The element %d in this detector %d has not enough statistic to be fitted"
5190 ,idect-(fCount[i]-(fNfragZ[i]*fNfragRphi[i])),fCountDet[i]));
5193 AliInfo("The element has not enough statistic to be fitted");
5196 if ((i == 0) && (fDebug != 2)) {
5198 // Calcul the coef from the database choosen
5199 CalculChargeCoefMean(fCountDet[0],(Int_t) (idect-fDect1[0]),kFALSE);
5201 // Fill the coefCH[2304] with negative value to say: not fitted
5202 for (Int_t k = fRowMin[0]; k < fRowMax[0]; k++) {
5203 for (Int_t j = fColMin[0]; j < fColMax[0]; j++) {
5204 if (GetChamber(fCountDet[0]) == 2) {
5205 fCoefCH[(Int_t)(j*12+k)] = -TMath::Abs(fChargeCoef[3]);
5207 if (GetChamber(fCountDet[0]) != 2) {
5208 fCoefCH[(Int_t)(j*16+k)] = -TMath::Abs(fChargeCoef[3]);
5213 // End of one detector
5214 if ((idect == (fCount[0]-1))) {
5215 FillVectorFitCH((Int_t) fCountDet[0]);
5217 for (Int_t k = 0; k < 2304; k++) {
5224 if ((i == 1) && (fDebug != 2)) {
5226 CalculVdriftCoefMean(fCountDet[1],(Int_t) (idect-fDect1[1]));
5227 CalculT0CoefMean(fCountDet[1],(Int_t) (idect-fDect1[1]));
5229 // Put the default value
5230 if ((fDebug == 1) ||
5234 fCoefVdrift[0]->SetBinContent(idect-fDect1[1]+1,-fVdriftCoef[2]);
5235 fCoefT0[0]->SetBinContent(idect-fDect1[1]+1,-fT0Coef[2]);
5238 fCoefVdrift[1]->SetBinContent(idect-fDect1[1]+1,-fVdriftCoef[2]);
5239 fCoefT0[1]->SetBinContent(idect-fDect1[1]+1,-fT0Coef[2]);
5243 // Put the default value
5245 fVdriftCoef[0] = fVdriftCoef[2];
5246 fVdriftCoef[1] = fVdriftCoef[2];
5248 fT0Coef[0] = fT0Coef[2];
5249 fT0Coef[1] = fT0Coef[2];
5253 // Fill the tree if end of a detector.
5254 // The pointer to the branch stays with the default value 1.5!!!
5256 // Pointer to the branch
5257 for (Int_t k = fRowMin[1]; k < fRowMax[1]; k++) {
5258 for (Int_t j = fColMin[1]; j < fColMax[1]; j++) {
5259 if (GetChamber(fCountDet[1]) == 2) {
5260 fVdriftPad[(Int_t)(j*12+k)] = -TMath::Abs(fVdriftCoef[2]);
5262 if (GetChamber(fCountDet[1]) != 2) {
5263 fVdriftPad[(Int_t)(j*16+k)] = -TMath::Abs(fVdriftCoef[2]);
5268 // End of one detector
5269 if ((idect == (fCount[1]-1)) && (fDebug != 2)) {
5270 FillTreeVdrift((Int_t) fCountDet[1]);
5274 // Fill the tree if end of a detector.
5275 // The pointer to the branch stays with the default value 1.5!!!
5276 // Pointer to the branch
5277 for (Int_t k = fRowMin[1]; k < fRowMax[1]; k++) {
5278 for (Int_t j = fColMin[1]; j < fColMax[1]; j++) {
5279 if (GetChamber(fCountDet[1]) == 2) {
5280 fT0Pad[(Int_t)(j*12+k)] = -TMath::Abs(fT0Coef[2]);
5282 if (GetChamber(fCountDet[1]) != 2) {
5283 fT0Pad[(Int_t)(j*16+k)] = -TMath::Abs(fT0Coef[2]);
5288 // End of one detector
5289 if ((idect == (fCount[1]-1)) && (fDebug != 2)) {
5290 FillTreeT0((Int_t) fCountDet[1]);
5295 if ((i == 2) && (fDebug != 2)) {
5297 CalculPRFCoefMean(fCountDet[2],(Int_t) (idect-fDect1[2]));
5299 if ((fDebug == 1) ||
5301 fCoefPRF[0]->SetBinContent(idect-fDect1[2]+1,fPRFCoef[1]);
5305 fPRFCoef[0] = fPRFCoef[1];
5309 // Fill the tree if end of a detector.
5310 // The pointer to the branch stays with the default value 1.5!!!
5311 // Pointer to the branch
5312 for (Int_t k = fRowMin[2]; k < fRowMax[2]; k++) {
5313 for (Int_t j = fColMin[2]; j < fColMax[2]; j++) {
5314 if((parCom->GetColMax(GetPlane(fCountDet[2])) != (j+1)) && (j != 0)){
5315 if (GetChamber(fCountDet[2]) == 2) {
5316 fPRFPad[(Int_t)(j*12+k)] = -fPRFCoef[1];
5318 if (GetChamber(fCountDet[2]) != 2) {
5319 fPRFPad[(Int_t)(j*16+k)] = -fPRFCoef[1];
5324 if (GetChamber(fCountDet[2]) == 2) {
5325 fPRFPad[(Int_t)(j*12+k)] = -((Float_t) cal->GetPRFWidth(fCountDet[2],j,k));
5327 if (GetChamber(fCountDet[2]) != 2) {
5328 fPRFPad[(Int_t)(j*16+k)] = -((Float_t) cal->GetPRFWidth(fCountDet[2],j,k));
5332 if (GetChamber(fCountDet[2]) == 2) {
5333 fPRFPad[(Int_t)(j*12+k)] = -((Float_t) GetPRFDefault(GetPlane(fCountDet[2])));
5335 if (GetChamber(fCountDet[2]) != 2) {
5336 fPRFPad[(Int_t)(j*16+k)] = -((Float_t) GetPRFDefault(GetPlane(fCountDet[2])));
5343 // End of one detector
5344 if ((idect == (fCount[2]-1)) && (fDebug != 2)) {
5345 FillTreePRF((Int_t) fCountDet[2]);
5354 //____________Functions for initialising the AliTRDCalibra in the code_________
5355 Bool_t AliTRDCalibra::FillInfosFit(Int_t idect, Int_t i)
5358 // Fill the coefficients found with the fits or other
5359 // methods from the Fit functions
5362 // Get the parameter object
5363 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
5365 AliInfo("Could not get CommonParam Manager");
5370 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
5372 AliInfo("Could not get calibDB");
5376 if ((i == 0) && (fDebug != 2)) {
5377 // Fill the coefCH[2304] with fChargeCoef[0]
5378 // that would be negativ only if the fit failed totally
5379 for (Int_t k = fRowMin[0]; k < fRowMax[0]; k++) {
5380 for (Int_t j = fColMin[0]; j < fColMax[0]; j++) {
5381 if (GetChamber(fCountDet[0]) == 2) {
5382 fCoefCH[(Int_t)(j*12+k)] = fChargeCoef[0];
5384 if (GetChamber(fCountDet[0]) != 2) {
5385 fCoefCH[(Int_t)(j*16+k)] = fChargeCoef[0];
5389 // End of one detector
5390 if ((idect == (fCount[0]-1))) {
5391 FillVectorFitCH((Int_t) fCountDet[0]);
5393 for (Int_t k = 0; k < 2304; k++) {
5399 if ((i == 1) && (fDebug != 2)) {
5402 // Pointer to the branch: fVdriftCoef[1] will ne negativ only if the fit failed totally
5403 for (Int_t k = fRowMin[1]; k < fRowMax[1]; k++) {
5404 for (Int_t j = fColMin[1]; j < fColMax[1]; j++) {
5405 if (GetChamber(fCountDet[1]) == 2) {
5406 fVdriftPad[(Int_t)(j*12+k)]=fVdriftCoef[1];
5408 if (GetChamber(fCountDet[1]) != 2) {
5409 fVdriftPad[(Int_t)(j*16+k)]=fVdriftCoef[1];
5413 // End of one detector
5414 if ((idect == (fCount[1]-1)) && (fDebug != 2)) {
5415 FillTreeVdrift((Int_t) fCountDet[1]);
5419 // Pointer to the branch: fT0Coef[1] will ne negativ only if the fit failed totally
5420 for (Int_t k = fRowMin[1]; k < fRowMax[1]; k++) {
5421 for (Int_t j = fColMin[1]; j < fColMax[1]; j++) {
5422 if (GetChamber(fCountDet[1]) == 2) {
5423 fT0Pad[(Int_t)(j*12+k)]=fT0Coef[1];
5425 if (GetChamber(fCountDet[1]) != 2) {
5426 fT0Pad[(Int_t)(j*16+k)]=fT0Coef[1];
5430 // End of one detector
5431 if ((idect == (fCount[1]-1)) && (fDebug != 2)) {
5432 FillTreeT0((Int_t) fCountDet[1]);
5437 if ((i == 2) && (fDebug != 2)) {
5438 // Pointer to the branch
5439 for (Int_t k = fRowMin[2]; k < fRowMax[2]; k++) {
5440 for (Int_t j = fColMin[2]; j < fColMax[2]; j++) {
5441 if ((parCom->GetColMax(GetPlane(fCountDet[2])) != (j+1)) && (j != 0)) {
5442 if (GetChamber(fCountDet[2]) == 2) {
5443 fPRFPad[(Int_t)(j*12+k)] = fPRFCoef[0];
5445 if (GetChamber(fCountDet[2]) != 2) {
5446 fPRFPad[(Int_t)(j*16+k)] = fPRFCoef[0];
5451 if (GetChamber(fCountDet[2]) == 2) {
5452 fPRFPad[(Int_t)(j*12+k)] = (Float_t) cal->GetPRFWidth(fCountDet[2],j,k);
5454 if (GetChamber(fCountDet[2]) != 2) {
5455 fPRFPad[(Int_t)(j*16+k)] = (Float_t) cal->GetPRFWidth(fCountDet[2],j,k);
5459 if (GetChamber(fCountDet[2]) == 2) {
5460 fPRFPad[(Int_t)(j*12+k)] = (Float_t) GetPRFDefault(GetPlane(fCountDet[2]));
5462 if (GetChamber(fCountDet[2]) != 2) {
5463 fPRFPad[(Int_t)(j*16+k)] = (Float_t) GetPRFDefault(GetPlane(fCountDet[2]));
5469 // End of one detector
5470 if ((idect == (fCount[2]-1)) && (fDebug != 2)) {
5471 FillTreePRF((Int_t) fCountDet[2]);
5479 //____________Functions for initialising the AliTRDCalibra in the code_________
5480 Bool_t AliTRDCalibra::WriteFitInfos(Int_t i)
5483 // In the case the user wants to write a file with a tree of the found
5484 // coefficients for the calibration before putting them in the database
5487 TFile *fout = TFile::Open(fWriteNameCoef,"UPDATE");
5488 // Check if the file could be opened
5489 if (!fout || !fout->IsOpen()) {
5490 AliInfo("No File found!");
5494 if ((i == 0) && (fDebug != 2)) {
5496 if ((fDebug == 1) ||
5501 if ((fDebug == 4) ||
5506 fout->WriteTObject(fGain,fGain->GetName(),(Option_t *) "writedelete");
5509 if ((i == 1) && (fDebug != 2)) {
5512 if ((fDebug == 1) ||
5517 if ((fDebug == 4) ||
5522 fout->WriteTObject(fVdrift,fVdrift->GetName(),(Option_t *) "writedelete");
5525 if ((fDebug == 1) ||
5530 if ((fDebug == 4) ||
5535 fout->WriteTObject(fT0,fT0->GetName(),(Option_t *) "writedelete");
5538 if ((i == 2) && (fDebug != 2)) {
5540 if ((fDebug == 1) ||
5545 if ((fDebug == 4) ||
5550 fout->WriteTObject(fPRF,fPRF->GetName(),(Option_t *) "writedelete");
5560 //____________Fill the Error histos in case of fDebug == 1_____________________
5563 //_____________________________________________________________________________
5564 void AliTRDCalibra::ErrorPRF()
5567 // Fill the error histos for fDebug = 1 and fDebug = 4 from the delta histos
5570 for (Int_t k= 0; k < fDeltaPRF->GetNbinsX(); k++) {
5571 if (fDeltaPRF->GetBinContent(k+1) != 0.0) {
5572 fErrorPRF->Fill(fDeltaPRF->GetBinContent(k+1));
5578 //_____________________________________________________________________________
5579 void AliTRDCalibra::ErrorCH()
5582 // Fill the error histos for fDebug = 1 and fDebug = 4 from the delta histos
5585 for (Int_t k= 0; k < fDeltaCharge[0]->GetNbinsX(); k++) {
5586 if (fDeltaCharge[0]->GetBinContent(k+1) != 0.0) {
5587 fErrorCharge[0]->Fill(fDeltaCharge[0]->GetBinContent(k+1));
5590 if (fMeanChargeOn) {
5591 for (Int_t k= 0; k < fDeltaCharge[1]->GetNbinsX(); k++) {
5592 if (fDeltaCharge[1]->GetBinContent(k+1) != 0.0) {
5593 fErrorCharge[1]->Fill(fDeltaCharge[1]->GetBinContent(k+1));
5597 if (fFitChargeBisOn ) {
5598 for (Int_t k= 0; k < fDeltaCharge[2]->GetNbinsX(); k++) {
5599 if (fDeltaCharge[2]->GetBinContent(k+1) != 0.0) {
5600 fErrorCharge[2]->Fill(fDeltaCharge[2]->GetBinContent(k+1));
5607 //_____________________________________________________________________________
5608 void AliTRDCalibra::ErrorPH()
5611 // Fill the error histos for fDebug = 1 and fDebug = 4 from the delta histos
5615 for (Int_t k = 0; k < fDeltaVdrift[0]->GetNbinsX(); k++) {
5616 if (fDeltaVdrift[0]->GetBinContent(k+1) != 0.0) {
5617 fErrorVdrift[0]->Fill(fDeltaVdrift[0]->GetBinContent(k+1));
5621 for (Int_t k = 0; k < fDeltaVdrift[1]->GetNbinsX(); k++) {
5622 if (fDeltaVdrift[1]->GetBinContent(k+1) != 0.0) {
5623 fErrorVdrift[1]->Fill(fDeltaVdrift[1]->GetBinContent(k+1));
5629 //_____________________________________________________________________________
5630 void AliTRDCalibra::ErrorT0()
5633 // Fill the error histos for fDebug = 1 and fDebug = 4 from the delta histos
5637 for (Int_t k = 0; k < fDeltaT0[0]->GetNbinsX(); k++) {
5638 if (fDeltaT0[0]->GetBinContent(k+1) != 0.0) {
5639 fErrorT0[0]->Fill(fDeltaT0[0]->GetBinContent(k+1));
5643 for (Int_t k = 0; k < fDeltaT0[1]->GetNbinsX(); k++) {
5644 if (fDeltaT0[1]->GetBinContent(k+1) != 0.0) {
5645 fErrorT0[1]->Fill(fDeltaT0[1]->GetBinContent(k+1));
5652 //____________Fill Coef DB in case of visualisation of one detector____________
5655 //_____________________________________________________________________________
5656 void AliTRDCalibra::FillCoefVdriftDB()
5659 // Fill the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5662 for (Int_t row = fRowMin[1]; row < fRowMax[1]; row++) {
5663 for (Int_t col = fColMin[1]; col < fColMax[1]; col++) {
5664 fCoefVdriftDB[1]->SetBinContent(row+1,col+1,TMath::Abs(fVdriftCoef[1]));
5666 fCoefVdriftDB[0]->SetBinContent(row+1,col+1,TMath::Abs(fVdriftCoef[0]));
5673 //_____________________________________________________________________________
5674 void AliTRDCalibra::FillCoefT0DB()
5677 // Fill the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5680 for (Int_t row = fRowMin[1]; row < fRowMax[1]; row++) {
5681 for (Int_t col = fColMin[1]; col < fColMax[1]; col++) {
5682 fCoefT0DB[1]->SetBinContent(row+1,col+1,TMath::Abs(fT0Coef[1]));
5684 fCoefT0DB[0]->SetBinContent(row+1,col+1,TMath::Abs(fT0Coef[0]));
5691 //_____________________________________________________________________________
5692 void AliTRDCalibra::FillCoefChargeDB()
5695 // Fill the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5698 for (Int_t row = fRowMin[0]; row < fRowMax[0]; row++) {
5699 for (Int_t col = fColMin[0]; col < fColMax[0]; col++) {
5700 if (fMeanChargeOn) {
5701 fCoefChargeDB[1]->SetBinContent(row+1,col+1,TMath::Abs(fChargeCoef[1]));
5703 if (fFitChargeBisOn) {
5704 fCoefChargeDB[2]->SetBinContent(row+1,col+1,TMath::Abs(fChargeCoef[2]));
5706 fCoefChargeDB[0]->SetBinContent(row+1,col+1,TMath::Abs(fChargeCoef[0]));
5712 //_____________________________________________________________________________
5713 void AliTRDCalibra::FillCoefPRFDB()
5716 // Fill the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5719 for (Int_t row = fRowMin[2]; row < fRowMax[2]; row++) {
5720 for (Int_t col = fColMin[2]; col < fColMax[2]; col++) {
5721 fCoefPRFDB->SetBinContent(row+1,col+1,fPRFCoef[0]);
5728 //____________Plot histos CoefPRF....__________________________________________
5731 //_____________________________________________________________________________
5732 void AliTRDCalibra::PlotCH()
5735 // Plot the histos for fDebug = 1 and fDebug = 4 for the errors
5738 TCanvas *cch1 = new TCanvas("cch1","",50,50,600,800);
5740 TLegend *legch1 = new TLegend(0.4,0.6,0.89,0.89);
5741 legch1->AddEntry(fCoefCharge[3],"f_{g} simulated","l");
5742 if (fMeanChargeOn) {
5743 legch1->AddEntry(fCoefCharge[1],"f_{g} mean","p");
5745 legch1->AddEntry(fCoefCharge[0],"f_{g} fit","p");
5746 if (fFitChargeBisOn ) {
5747 legch1->AddEntry(fCoefCharge[2],"f_{g} fitbis","p");
5750 fCoefCharge[0]->Draw("E2");
5751 if (fMeanChargeOn) {
5752 fCoefCharge[1]->Draw("E2 same");
5754 if (fFitChargeBisOn ) {
5755 fCoefCharge[2]->Draw("E2 same");
5757 fCoefCharge[3]->Draw("same");
5758 legch1->Draw("same");
5760 TCanvas *cch2 = new TCanvas("cch2","",50,50,600,800);
5763 TLegend *legch2 = new TLegend(0.4,0.6,0.89,0.89);
5764 if (fMeanChargeOn) {
5765 legch2->AddEntry(fErrorCharge[1],"f_{g} mean","l");
5767 legch2->AddEntry(fErrorCharge[0],"f_{g} fit","l");
5768 if (fFitChargeBisOn) {
5769 legch2->AddEntry(fErrorCharge[2],"f_{g} fitbis","l");
5771 fErrorCharge[0]->Draw();
5772 if (fMeanChargeOn) {
5773 fErrorCharge[1]->Draw("same");
5775 if (fFitChargeBisOn) {
5776 fErrorCharge[2]->Draw("same");
5778 legch2->Draw("same");
5780 TLegend *legch3 = new TLegend(0.4,0.6,0.89,0.89);
5781 if (fMeanChargeOn) {
5782 legch3->AddEntry(fDeltaCharge[1],"mean","p");
5784 legch3->AddEntry(fDeltaCharge[0],"fit","p");
5785 if (fFitChargeBisOn) {
5786 legch3->AddEntry(fDeltaCharge[2],"fit","p");
5788 fDeltaCharge[0]->Draw("E2");
5789 if (fMeanChargeOn) {
5790 fDeltaCharge[1]->Draw("E2 same");
5792 if (fFitChargeBisOn) {
5793 fDeltaCharge[2]->Draw("E2 same");
5795 legch3->Draw("same");
5799 //_____________________________________________________________________________
5800 void AliTRDCalibra::PlotPH()
5803 // Plot the histos for fDebug = 1 and fDebug = 4 for the errors
5806 TCanvas *cph1 = new TCanvas("cph1","",50,50,600,800);
5808 TLegend *legph1 = new TLegend(0.4,0.6,0.89,0.89);
5809 legph1->AddEntry(fCoefVdrift[2],"v_{real} simulated","l");
5810 legph1->AddEntry(fCoefVdrift[1],"v_{sm} slope 1 method","p");
5813 legph1->AddEntry(fCoefVdrift[0],"v_{fit} fit","p");
5815 fCoefVdrift[1]->Draw("E2");
5816 fCoefVdrift[2]->Draw("same");
5818 fCoefVdrift[0]->Draw("E2 same");
5820 legph1->Draw("same");
5822 TCanvas *cph2 = new TCanvas("cph2","",50,50,600,800);
5825 TLegend *legph2 = new TLegend(0.4,0.6,0.89,0.89);
5826 legph2->AddEntry(fErrorVdrift[1],"v_{sm} slope 1 method","l");
5828 legph2->AddEntry(fErrorVdrift[0],"v_{fit} fit","l");
5830 fErrorVdrift[1]->Draw();
5832 fErrorVdrift[0]->Draw("l,same");
5834 legph2->Draw("same");
5836 TLegend *legph3 = new TLegend(0.4,0.6,0.89,0.89);
5837 legph3->AddEntry(fDeltaVdrift[1],"v_{sm} slope 1 method","p");
5839 legph3->AddEntry(fDeltaVdrift[0],"v_{fit} fit","p");
5841 fDeltaVdrift[1]->Draw("E2");
5843 fDeltaVdrift[0]->Draw("E2 same");
5845 legph3->Draw("same");
5849 //_____________________________________________________________________________
5850 void AliTRDCalibra::PlotT0()
5853 // Plot the histos for fDebug = 1 and fDebug = 4 for the errors
5856 TCanvas *ct01 = new TCanvas("ct01","",50,50,600,800);
5858 TLegend *legt01 = new TLegend(0.4,0.6,0.89,0.89);
5859 legt01->AddEntry(fCoefT0[2],"t0 simulated","l");
5860 legt01->AddEntry(fCoefT0[1],"t0 slope 1 method","p");
5863 legt01->AddEntry(fCoefT0[0],"t0 fit","p");
5865 fCoefT0[1]->Draw("E2");
5866 fCoefT0[2]->Draw("same");
5868 fCoefT0[0]->Draw("E2 same");
5870 legt01->Draw("same");
5872 TCanvas *ct02 = new TCanvas("ct02","",50,50,600,800);
5875 TLegend *legt02 = new TLegend(0.4,0.6,0.89,0.89);
5876 legt02->AddEntry(fErrorT0[1],"t0 slope 1 method","l");
5878 legt02->AddEntry(fErrorT0[0],"t0 fit","l");
5880 fErrorT0[1]->Draw();
5882 fErrorT0[0]->Draw("l,same");
5884 legt02->Draw("same");
5886 TLegend *legt03 = new TLegend(0.4,0.6,0.89,0.89);
5887 legt03->AddEntry(fDeltaT0[1],"t0 slope 1 method","p");
5889 legt03->AddEntry(fDeltaT0[0],"t0 fit","p");
5891 fDeltaT0[1]->Draw("E2");
5893 fDeltaT0[0]->Draw("E2 same");
5895 legt03->Draw("same");
5899 //_____________________________________________________________________________
5900 void AliTRDCalibra::PlotPRF()
5903 // Plot the histos for fDebug = 1 and fDebug = 4 for the errors
5906 TCanvas *cprf1 = new TCanvas("cprf1","",50,50,600,800);
5908 TLegend *legprf1 = new TLegend(0.4,0.6,0.89,0.89);
5909 legprf1->AddEntry(fCoefPRF[1],"#sigma_{real} simulated","l");
5910 legprf1->AddEntry(fCoefPRF[0],"#sigma_{fit} reconstructed","p");
5912 fCoefPRF[0]->Draw("E2");
5913 fCoefPRF[1]->Draw("same");
5914 legprf1->Draw("same");
5916 TCanvas *cprf2 = new TCanvas("cprf2","",50,50,600,800);
5919 TLegend *legprf2 = new TLegend(0.4,0.6,0.89,0.89);
5920 legprf2->AddEntry(fErrorPRF,"#sigma_{fit} reconstructed","l");
5921 fErrorPRF->Draw("");
5922 legprf2->Draw("same");
5924 TLegend *legprf3 = new TLegend(0.4,0.6,0.89,0.89);
5925 legprf3->AddEntry(fDeltaPRF,"#sigma_{fit} reconstructed","p");
5926 fDeltaPRF->Draw("E2");
5927 legprf3->Draw("same");
5932 //____________Plot histos DB___________________________________________________
5935 //_____________________________________________________________________________
5936 void AliTRDCalibra::PlotCHDB()
5939 // Plot the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5942 TCanvas *cchdb = new TCanvas("cchdb","",50,50,600,800);
5943 if ((fFitChargeBisOn) && (fMeanChargeOn)) {
5946 fCoefChargeDB[0]->Draw("LEGO");
5948 fCoefChargeDB[1]->Draw("LEGO");
5950 fCoefChargeDB[2]->Draw("LEGO");
5952 if ((!fFitChargeBisOn) && (fMeanChargeOn)) {
5955 fCoefChargeDB[0]->Draw("LEGO");
5957 fCoefChargeDB[1]->Draw("LEGO");
5961 fCoefChargeDB[0]->Draw("LEGO");
5966 //_____________________________________________________________________________
5967 void AliTRDCalibra::PlotPHDB()
5970 // Plot the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5973 TCanvas *cphdb = new TCanvas("cphdb","",50,50,600,800);
5977 fCoefVdriftDB[0]->Draw("LEGO");
5979 fCoefVdriftDB[1]->Draw("LEGO");
5983 fCoefVdriftDB[1]->Draw("LEGO");
5988 //_____________________________________________________________________________
5989 void AliTRDCalibra::PlotT0DB()
5992 // Plot the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5995 TCanvas *ct0db = new TCanvas("ct0db","",50,50,600,800);
5999 fCoefT0DB[0]->Draw("LEGO");
6001 fCoefT0DB[1]->Draw("LEGO");
6005 fCoefT0DB[1]->Draw("LEGO");
6010 //_____________________________________________________________________________
6011 void AliTRDCalibra::PlotPRFDB()
6014 // Plot the histos for fDebug = 3 and fDebug = 4 to visualise the detector
6017 TCanvas *cprfdb = new TCanvas("cprfdb","",50,50,600,800);
6019 fCoefPRFDB->Draw("LEGO");
6024 //____________Write histos Coef________________________________________________
6027 //_____________________________________________________________________________
6028 void AliTRDCalibra::WriteCH(TFile *fout)
6031 // If wanted, write the debug histos for fDebug = 1 and fDebug = 4
6034 fout->WriteTObject(fCoefCharge[0],fCoefCharge[0]->GetName(),(Option_t *) "OverWrite");
6035 if (fMeanChargeOn) {
6036 fout->WriteTObject(fCoefCharge[1],fCoefCharge[1]->GetName(),(Option_t *) "OverWrite");
6038 if (fFitChargeBisOn) {
6039 fout->WriteTObject(fCoefCharge[2],fCoefCharge[2]->GetName(),(Option_t *) "OverWrite");
6042 fout->WriteTObject(fCoefCharge[3],fCoefCharge[3]->GetName(),(Option_t *) "OverWrite");
6044 fout->WriteTObject(fDeltaCharge[0],fDeltaCharge[0]->GetName(),(Option_t *) "OverWrite");
6045 if (fMeanChargeOn) {
6046 fout->WriteTObject(fDeltaCharge[1],fDeltaCharge[1]->GetName(),(Option_t *) "OverWrite");
6048 if (fFitChargeBisOn) {
6049 fout->WriteTObject(fDeltaCharge[2],fDeltaCharge[2]->GetName(),(Option_t *) "OverWrite");
6052 fout->WriteTObject(fErrorCharge[0],fErrorCharge[0]->GetName(),(Option_t *) "OverWrite");
6053 if (fMeanChargeOn) {
6054 fout->WriteTObject(fErrorCharge[1],fErrorCharge[1]->GetName(),(Option_t *) "OverWrite");
6056 if (fFitChargeBisOn) {
6057 fout->WriteTObject(fErrorCharge[2],fErrorCharge[2]->GetName(),(Option_t *) "OverWrite");
6062 //_____________________________________________________________________________
6063 void AliTRDCalibra::WritePH(TFile *fout)
6066 // If wanted, write the debug histos for fDebug = 1 and fDebug = 4
6070 fout->WriteTObject(fCoefVdrift[0],fCoefVdrift[0]->GetName(),(Option_t *) "OverWrite");
6072 fout->WriteTObject(fCoefVdrift[1],fCoefVdrift[1]->GetName(),(Option_t *) "OverWrite");
6073 fout->WriteTObject(fCoefVdrift[2],fCoefVdrift[2]->GetName(),(Option_t *) "OverWrite");
6076 fout->WriteTObject(fDeltaVdrift[0],fDeltaVdrift[0]->GetName(),(Option_t *) "OverWrite");
6078 fout->WriteTObject(fDeltaVdrift[1],fDeltaVdrift[1]->GetName(),(Option_t *) "OverWrite");
6081 fout->WriteTObject(fErrorVdrift[0],fErrorVdrift[0]->GetName(),(Option_t *) "OverWrite");
6083 fout->WriteTObject(fErrorVdrift[1],fErrorVdrift[1]->GetName(),(Option_t *) "OverWrite");
6087 //_____________________________________________________________________________
6088 void AliTRDCalibra::WriteT0(TFile *fout)
6091 // If wanted, write the debug histos for fDebug = 1 and fDebug = 4
6095 fout->WriteTObject(fCoefT0[0],fCoefT0[0]->GetName(),(Option_t *) "OverWrite");
6097 fout->WriteTObject(fCoefT0[1],fCoefT0[1]->GetName(),(Option_t *) "OverWrite");
6098 fout->WriteTObject(fCoefT0[2],fCoefT0[2]->GetName(),(Option_t *) "OverWrite");
6101 fout->WriteTObject(fDeltaT0[0],fDeltaT0[0]->GetName(),(Option_t *) "OverWrite");
6103 fout->WriteTObject(fDeltaT0[1],fDeltaT0[1]->GetName(),(Option_t *) "OverWrite");
6106 fout->WriteTObject(fErrorT0[0],fErrorT0[0]->GetName(),(Option_t *) "OverWrite");
6108 fout->WriteTObject(fErrorT0[1],fErrorT0[1]->GetName(),(Option_t *) "OverWrite");
6112 //________________________________________________________________________________
6113 void AliTRDCalibra::WritePRF(TFile *fout)
6116 // If wanted, write the debug histos for fDebug = 1 and fDebug = 4
6119 fout->WriteTObject(fCoefPRF[0],fCoefPRF[0]->GetName(),(Option_t *) "OverWrite");
6120 fout->WriteTObject(fCoefPRF[1],fCoefPRF[1]->GetName(),(Option_t *) "OverWrite");
6122 fout->WriteTObject(fDeltaPRF,fDeltaPRF->GetName(), (Option_t *)"OverWrite");
6123 fout->WriteTObject(fErrorPRF,fErrorPRF->GetName(), (Option_t *)"OverWrite");
6128 //____________Write DB Histos__________________________________________________
6131 //_____________________________________________________________________________
6132 void AliTRDCalibra::WriteCHDB(TFile *fout)
6135 // If wanted, write the debug histos for fDebug = 3 and fDebug = 4
6138 fout->WriteTObject(fCoefChargeDB[0],fCoefChargeDB[0]->GetName(),(Option_t *) "OverWrite");
6139 if (fMeanChargeOn) {
6140 fout->WriteTObject(fCoefChargeDB[1],fCoefChargeDB[1]->GetName(),(Option_t *) "OverWrite");
6142 if (fFitChargeBisOn ) {
6143 fout->WriteTObject(fCoefChargeDB[2],fCoefChargeDB[2]->GetName(),(Option_t *) "OverWrite");
6148 //_____________________________________________________________________________
6149 void AliTRDCalibra::WritePHDB(TFile *fout)
6152 // If wanted, write the debug histos for fDebug = 3 and fDebug = 4
6156 fout->WriteTObject(fCoefVdriftDB[0],fCoefVdriftDB[0]->GetName(),(Option_t *) "OverWrite");
6158 fout->WriteTObject(fCoefVdriftDB[1],fCoefVdriftDB[1]->GetName(),(Option_t *) "OverWrite");
6162 //_____________________________________________________________________________
6163 void AliTRDCalibra::WriteT0DB(TFile *fout)
6166 // If wanted, write the debug histos for fDebug = 3 and fDebug = 4
6170 fout->WriteTObject(fCoefT0DB[0],fCoefT0DB[0]->GetName(),(Option_t *) "OverWrite");
6172 fout->WriteTObject(fCoefT0DB[1],fCoefT0DB[1]->GetName(),(Option_t *) "OverWrite");
6176 //_____________________________________________________________________________
6177 void AliTRDCalibra::WritePRFDB(TFile *fout)
6180 // If wanted, write the debug histos for fDebug = 3 and fDebug = 4
6183 fout->WriteTObject(fCoefPRFDB,fCoefPRFDB->GetName(),(Option_t *) "OverWrite");
6188 //____________Calcul Coef Mean_________________________________________________
6191 //_____________________________________________________________________________
6192 Bool_t AliTRDCalibra::CalculT0CoefMean(Int_t dect, Int_t idect)
6195 // For the detector Dect calcul the mean time 0
6196 // for the calibration group idect from the choosen database
6199 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
6201 AliInfo("Could not get calibDB Manager");
6207 if ((fDebug != 2) && fAccCDB) {
6209 for (Int_t row = fRowMin[1]; row < fRowMax[1]; row++) {
6210 for (Int_t col = fColMin[1]; col < fColMax[1]; col++) {
6214 fT0Coef[2] += (Float_t) cal->GetT0(dect,col,row);
6218 fT0Coef[2] += (Float_t) cal->GetT0Average(dect);
6223 fT0Coef[2] = fT0Coef[2] / ((fColMax[1]-fColMin[1])*(fRowMax[1]-fRowMin[1]));
6224 if ((fDebug == 1) ||
6226 fCoefT0[2]->SetBinContent(idect+1,fT0Coef[2]);
6235 //_____________________________________________________________________________
6236 Bool_t AliTRDCalibra::CalculChargeCoefMean(Int_t dect, Int_t idect, Bool_t vrai)
6239 // For the detector Dect calcul the mean gain factor
6240 // for the calibration group idect from the choosen database
6243 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
6245 AliInfo("Could not get calibDB Manager");
6248 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
6250 AliInfo("Could not get CommonParam Manager");
6254 fChargeCoef[3] = 0.0;
6258 for (Int_t row = fRowMin[0]; row < fRowMax[0]; row++) {
6259 for (Int_t col = fColMin[0]; col < fColMax[0]; col++) {
6264 fChargeCoef[3] += (Float_t) cal->GetGainFactor(dect,col,row);
6266 if (vrai && fAccCDB) {
6267 fScaleFitFactor += (Float_t) cal->GetGainFactor(dect,col,row);
6270 fChargeCoef[3] += 1.0;
6272 if (vrai && (!fAccCDB)) {
6273 fScaleFitFactor += 1.0;
6279 fChargeCoef[3] += (Float_t) cal->GetGainFactorAverage(dect);
6281 if (vrai && fAccCDB) {
6282 fScaleFitFactor += ((Float_t) cal->GetGainFactorAverage(dect));
6285 fChargeCoef[3] += 1.0;
6287 if (vrai && (!fAccCDB)) {
6288 fScaleFitFactor += 1.0;
6294 fChargeCoef[3] = fChargeCoef[3] / ((fColMax[0]-fColMin[0])*(fRowMax[0]-fRowMin[0]));
6295 if ((fDebug == 1) ||
6297 fCoefCharge[3]->SetBinContent(idect+1,fChargeCoef[3]);
6306 //_____________________________________________________________________________
6307 Bool_t AliTRDCalibra::CalculPRFCoefMean(Int_t dect, Int_t idect)
6310 // For the detector Dect calcul the mean sigma of pad response
6311 // function for the calibration group idect from the choosen database
6314 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
6316 AliInfo("Could not get calibDB Manager");
6320 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
6322 AliInfo("Could not get CommonParam Manager");
6331 for (Int_t row = fRowMin[2]; row < fRowMax[2]; row++) {
6332 for (Int_t col = fColMin[2]; col < fColMax[2]; col++) {
6333 if ((parCom->GetColMax(GetPlane(dect)) != (col+1)) && (col != 0)) {
6336 fPRFCoef[1] += (Float_t) cal->GetPRFWidth(dect,col,row);
6339 fPRFCoef[1] += GetPRFDefault(GetPlane(dect));
6346 fPRFCoef[1] = fPRFCoef[1]/cot;
6347 if ((fDebug == 1) ||
6349 fCoefPRF[1]->SetBinContent(idect+1,fPRFCoef[1]);
6353 if ((fDebug == 1) ||
6356 fCoefPRF[1]->SetBinContent(idect+1,cal->GetPRFWidth(dect,fColMin[2],fRowMin[2]));
6359 fCoefPRF[1]->SetBinContent(idect+1,GetPRFDefault(GetPlane(dect)));
6370 //_____________________________________________________________________________
6371 Bool_t AliTRDCalibra::CalculVdriftCoefMean(Int_t dect, Int_t idect)
6374 // For the detector dect calcul the mean drift velocity for the
6375 // calibration group idect from the choosen database
6378 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
6380 AliInfo("Could not get calibDB Manager");
6384 fVdriftCoef[2] = 0.0;
6387 for (Int_t row = fRowMin[1]; row < fRowMax[1]; row++) {
6388 for (Int_t col = fColMin[1]; col < fColMax[1]; col++) {
6393 fVdriftCoef[2] += (Float_t) cal->GetVdrift(dect,col,row);
6396 fVdriftCoef[2] += 1.5;
6402 fVdriftCoef[2] += (Float_t) cal->GetVdriftAverage(dect);
6405 fVdriftCoef[2] += 1.5;
6410 fVdriftCoef[2] = fVdriftCoef[2] / ((fColMax[1]-fColMin[1])*(fRowMax[1]-fRowMin[1]));
6411 if ((fDebug == 1) ||
6413 fCoefVdrift[2]->SetBinContent(idect+1,fVdriftCoef[2]);
6421 //_____________________________________________________________________________
6422 Float_t AliTRDCalibra::GetPRFDefault(Int_t plane) const
6425 // Default width of the PRF if there is no database as reference
6453 //____________Pad group calibration mode_______________________________________
6456 //_____________________________________________________________________________
6457 void AliTRDCalibra::ReconstructionRowPadGroup(Int_t idect, Int_t i)
6460 // For the calibration group idect in a detector calculate the
6461 // first and last row pad and col pad.
6462 // The pads in the interval will have the same calibrated coefficients
6472 if (fNfragZ[i] != 0) {
6473 posc = (Int_t) idect / fNfragZ[i];
6475 if (fNfragRphi[i] != 0) {
6476 posr = (Int_t) idect % fNfragZ[i];
6478 fRowMin[i] = posr * fNnZ[i];
6479 fRowMax[i] = (posr+1) * fNnZ[i];
6480 fColMin[i] = posc * fNnRphi[i];
6481 fColMax[i] = (posc+1) * fNnRphi[i];
6485 //_____________________________________________________________________________
6486 void AliTRDCalibra::CalculXBins(Int_t idect, Int_t i)
6489 // For the detector idect calcul the first Xbins
6494 AliInfo(Form("detector: %d", idect));
6498 Int_t sector = GetSector(idect);
6499 fXbins[i] += sector*(6*fDetChamb2[i]+6*4*fDetChamb0[i]);
6501 // In which chamber?
6502 Int_t chamber = GetChamber(idect);
6504 while (kc < chamber) {
6506 fXbins[i] += 6 * fDetChamb2[i];
6509 fXbins[i] += 6 * fDetChamb0[i];
6515 Int_t plane = GetPlane(idect);
6517 fXbins[i] += plane*fDetChamb2[i];
6520 fXbins[i] += plane*fDetChamb0[i];
6525 //_____________________________________________________________________________
6526 Int_t AliTRDCalibra::SearchInVector(Int_t group, Int_t i) const
6529 // Search if the calibration group "group" has already been
6530 // initialised by a previous track in the vector
6534 for (Int_t k = 0; k < (Int_t) fPlaCH->GetEntriesFast(); k++) {
6535 if (((AliTRDPlace *) fPlaCH->At(k))->GetPlace() == group) {
6543 for (Int_t k = 0; k < (Int_t) fPlaPH->GetEntriesFast(); k++) {
6544 if (((AliTRDPlace *) fPlaPH->At(k))->GetPlace() == group) {
6552 for (Int_t k = 0; k < (Int_t) fPlaPRF->GetEntriesFast(); k++) {
6553 if (((AliTRDPlace *) fPlaPRF->At(k))->GetPlace() == group) {
6564 //_____________________________________________________________________________
6565 Int_t AliTRDCalibra::SearchInTreeVector(TObjArray *vectorplace, Int_t group) const
6568 // Search if the calibration group "group" is present in the tree
6571 for (Int_t k = 0; k < (Int_t) vectorplace->GetEntriesFast(); k++) {
6572 if (((AliTRDPlace *) vectorplace->At(k))->GetPlace() == group) {
6581 //_____________________________________________________________________________
6582 Int_t AliTRDCalibra::SearchBin(Float_t value, Int_t i) const
6590 Int_t fbinmax = (Int_t) value;
6591 Int_t fNumberOfBin = -1;
6597 fNumberOfBin = fNumberBinCharge;
6604 fNumberOfBin = fNumberBinPRF;
6608 if ((value >= fbinmax) ||
6609 (value < fbinmin)) {
6614 reponse = (Int_t) ((fNumberOfBin*(value-fbinmin)) / (fbinmax-fbinmin));
6621 //_____________________________________________________________________________
6622 Bool_t AliTRDCalibra::UpdateVectorCH(Int_t group, Float_t value)
6625 // Fill the vector if a new calibration group "group" or update the
6626 // values of the calibration group "group" if already here
6630 Int_t bin = SearchBin(value,0);
6632 if ((bin < 0) || (bin >= fNumberBinCharge)) {
6637 Int_t place = SearchInVector(group,0);
6641 AliTRDPlace *placegroup = new AliTRDPlace();
6642 placegroup->SetPlace(group);
6643 fPlaCH->Add((TObject *) placegroup);
6645 AliTRDCTInfo *fCHInfo = new AliTRDCTInfo();
6646 UShort_t *entries = new UShort_t[fNumberBinCharge];
6648 for(Int_t k = 0; k < fNumberBinCharge; k++) {
6654 fCHInfo->SetEntries(entries);
6655 // Set in the vector
6656 fVectorCH->Add((TObject *) fCHInfo);
6658 // Group already exits
6661 AliTRDCTInfo *fCHInfo = new AliTRDCTInfo();
6663 fCHInfo = ((AliTRDCTInfo *) fVectorCH->At(place));
6664 UShort_t *entries = fCHInfo->GetEntries();
6668 fCHInfo->SetEntries(entries);
6669 // Update the vector
6670 fVectorCH->AddAt((TObject *) fCHInfo,place);
6677 //_____________________________________________________________________________
6678 Bool_t AliTRDCalibra::UpdateVectorPRF(Int_t group, Float_t x, Float_t y)
6681 // Fill the vector if a new calibration group "group" or update the
6682 // values of the calibration group "group" if already here
6686 Int_t bin = SearchBin(x,2);
6688 if ((bin < 0) || (bin >= fNumberBinPRF)) {
6693 Int_t place = SearchInVector(group,2);
6698 AliTRDPlace *placegroup = new AliTRDPlace();
6699 placegroup->SetPlace(group);
6700 fPlaPRF->Add((TObject *) placegroup);
6701 AliTRDPInfo *fPRFInfo = new AliTRDPInfo();
6703 Float_t *sum = new Float_t[fNumberBinPRF];
6704 Float_t *sumsquare = new Float_t[fNumberBinPRF];
6705 UShort_t *entries = new UShort_t[fNumberBinPRF];
6708 for (Int_t k = 0; k < fNumberBinPRF; k++) {
6716 sumsquare[bin] += y*y;
6720 fPRFInfo->SetSum(sum);
6721 fPRFInfo->SetSumSquare(sumsquare);
6722 fPRFInfo->SetEntries(entries);
6724 // Set in the vector
6725 fVectorPRF->Add((TObject *) fPRFInfo);
6728 // Group already exits
6731 AliTRDPInfo *fPRFInfo = new AliTRDPInfo();
6733 fPRFInfo = (AliTRDPInfo *) fVectorPRF->At(place);
6735 Float_t *sum = fPRFInfo->GetSum();
6736 Float_t *sumsquare = fPRFInfo->GetSumSquare();
6737 UShort_t *entries = fPRFInfo->GetEntries();
6740 Double_t calcul = (((Double_t) fPRFInfo->GetEntries()[bin])
6741 * ((Double_t) fPRFInfo->GetSum()[bin]) + (Double_t) y)
6742 / (((Double_t) fPRFInfo->GetEntries()[bin]) + 1);
6743 sum[bin] = (Float_t) calcul;
6744 Double_t calculsquare = (((Double_t) fPRFInfo->GetSumSquare()[bin])
6745 * ((Double_t) fPRFInfo->GetEntries()[bin]) + ((Double_t) y)*((Double_t) y))
6746 / (((Double_t) fPRFInfo->GetEntries()[bin]) + 1);
6747 sumsquare[bin] = (Float_t) calculsquare;
6751 fPRFInfo->SetSum(sum);
6752 fPRFInfo->SetSumSquare(sumsquare);
6753 fPRFInfo->SetEntries(entries);
6755 // Update the vector
6756 fVectorPRF->AddAt((TObject *) fPRFInfo,place);
6764 //_____________________________________________________________________________
6765 Bool_t AliTRDCalibra::UpdateVectorPH(Int_t group, Int_t time, Float_t value)
6768 // Fill the vector if a new calibration group "group" or update
6769 // the values of the calibration group "group" if already here
6776 (bin >= fTimeMax)) {
6781 Int_t place = SearchInVector(group,1);
6786 AliTRDPlace *placegroup = new AliTRDPlace();
6787 placegroup->SetPlace(group);
6788 fPlaPH->Add((TObject *) placegroup);
6789 AliTRDPInfo *fPHInfo = new AliTRDPInfo();
6791 Float_t *sum = new Float_t[fTimeMax];
6792 Float_t *sumsquare = new Float_t[fTimeMax];
6793 UShort_t *entries = new UShort_t[fTimeMax];
6796 for (Int_t k = 0; k < fTimeMax; k++) {
6804 sumsquare[bin] += value*value;
6808 fPHInfo->SetSum(sum);
6809 fPHInfo->SetSumSquare(sumsquare);
6810 fPHInfo->SetEntries(entries);
6812 // Set in the vector
6813 fVectorPH->Add((TObject *) fPHInfo);
6816 // Group already exits
6819 AliTRDPInfo *fPHInfo = new AliTRDPInfo();
6821 fPHInfo = (AliTRDPInfo *) fVectorPH->At(place);
6823 Float_t *sum = fPHInfo->GetSum();
6824 Float_t *sumsquare = fPHInfo->GetSumSquare();
6825 UShort_t *entries = fPHInfo->GetEntries();
6828 Double_t calcul = (((Double_t) fPHInfo->GetEntries()[bin])
6829 * ((Double_t) fPHInfo->GetSum()[bin]) + (Double_t) value)
6830 / (((Double_t) fPHInfo->GetEntries()[bin]) + 1);
6831 sum[bin] = (Float_t) calcul;
6832 Double_t calculsquare = ((((Double_t) fPHInfo->GetSumSquare()[bin])
6833 * ((Double_t) fPHInfo->GetEntries()[bin]))
6834 + (((Double_t) value) * ((Double_t)value)))
6835 / (((Double_t) fPHInfo->GetEntries()[bin]) + 1);
6836 sumsquare[bin] = (Float_t) calculsquare;
6840 fPHInfo->SetSum(sum);
6841 fPHInfo->SetSumSquare(sumsquare);
6842 fPHInfo->SetEntries(entries);
6844 // Update the vector
6845 fVectorPH->AddAt((TObject *) fPHInfo,place);
6853 //_____________________________________________________________________________
6854 TGraphErrors *AliTRDCalibra::ConvertVectorPHisto(AliTRDPInfo *pInfo
6855 , const Char_t *name) const
6858 // Convert the PInfo in a 1D grapherror, name must contains "PRF"
6859 // if PRF calibration and not "PRF" for Vdrift calibration
6862 TGraphErrors *histo;
6863 const Char_t *pattern1 = "PRF";
6870 Double_t step = 0.0;
6875 if (strstr(name,pattern1)) {
6876 ntimes = fNumberBinPRF;
6881 x = new Double_t[ntimes]; // Xaxis
6882 y = new Double_t[ntimes]; // Mean
6883 ex = new Double_t[ntimes]; // Nentries
6884 ey = new Double_t[ntimes]; // Sum of square/nentries
6887 if (!strstr(name,pattern1)) {
6892 step = (1.0 - (-1.0)) / fNumberBinPRF;
6893 min = -1.0 + step / 2.0;
6897 for (Int_t k = 0; k < ntimes; k++) {
6898 x[k] = min + k*step;
6902 // Fill only if there is more than 0 something
6903 if (pInfo->GetEntries()[k] > 0) {
6904 ex[k] = pInfo->GetEntries()[k];
6905 y[k] = pInfo->GetSum()[k];
6906 ey[k] = pInfo->GetSumSquare()[k];
6910 // Define the TGraphErrors
6911 histo = new TGraphErrors(ntimes,x,y,ex,ey);
6912 histo->SetTitle(name);
6917 //_____________________________________________________________________________
6918 TH1F *AliTRDCalibra::ConvertVectorCTHisto(AliTRDCTInfo *cTInfo
6919 , const Char_t * name) const
6922 // Convert the CTInfo in a 1D histo
6927 Int_t ntimes = fNumberBinCharge;
6928 UShort_t *entries = cTInfo->GetEntries();
6931 histo = new TH1F(name,name,fNumberBinCharge,0,300);
6934 for (Int_t k = 0; k < ntimes; k++) {
6935 histo->SetBinContent(k+1,entries[k]);
6936 histo->SetBinError(k+1,TMath::Sqrt(TMath::Abs(entries[k])));
6943 //_____________________________________________________________________________
6944 TTree *AliTRDCalibra::ConvertVectorCTTreeHisto(TObjArray *vVectorCT
6946 , const Char_t *name
6947 , const Char_t *nametitle) const
6950 // Convert the vector in a tree with two branchs: the group number
6951 // and the TH1F histo reconstructed from the vector
6954 // Size of the things
6955 Int_t ntotal = (Int_t) pPlaCT->GetEntriesFast();
6957 AliInfo("nothing to write!");
6958 TTree *treeCT = new TTree(name,nametitle);
6962 // Variable of the tree
6963 Int_t groupnumber = -1; // Group calibration
6965 TObjArray vectorCT = *vVectorCT;
6966 TObjArray plaCT = *pPlaCT;
6969 TTree *treeCT = new TTree(name,nametitle);
6970 treeCT->Branch("groupnumber",&groupnumber,"groupnumber/I");
6971 treeCT->Branch("histo","TH1F",&histo,32000,0);
6975 while (k < ntotal) {
6977 groupnumber = ((AliTRDPlace *) plaCT.At(0))->GetPlace();
6978 nome += groupnumber;
6979 histo = ConvertVectorCTHisto(((AliTRDCTInfo *) vectorCT.At(0)),nome);
6981 vectorCT.RemoveAt(0);
6982 vectorCT.Compress();
6992 //_____________________________________________________________________________
6993 TTree *AliTRDCalibra::ConvertVectorPTreeHisto(TObjArray *vVectorP
6995 , const Char_t *name
6996 , const Char_t *nametitle) const
6999 // Convert the vector in a tree with two branchs: the group number
7000 // and the TGraphErrors histo reconstructed from the vector.
7001 // The name must contain "PRF" for PRF calibration and not "PRF"
7002 // for Vdrift calibration
7005 // Size of the things
7006 Int_t ntotal = (Int_t) pPlaP->GetEntriesFast();
7008 AliInfo("nothing to write!");
7009 TTree *treeP = new TTree(name,nametitle);
7013 // Variable of the tree
7014 Int_t groupnumber = -1; // Group calibration
7015 TGraphErrors *histo = 0x0;
7016 TObjArray vectorP = *vVectorP;
7017 TObjArray plaP = *pPlaP;
7020 TTree *treeP = new TTree(name,nametitle);
7021 treeP->Branch("groupnumber",&groupnumber,"groupnumber/I");
7022 treeP->Branch("histo","TGraphErrors",&histo,32000,0);
7026 while (k < ntotal) {
7028 groupnumber = ((AliTRDPlace *) plaP.At(0))->GetPlace();
7029 nome += groupnumber;
7030 histo = ConvertVectorPHisto((AliTRDPInfo *) vectorP.At(0),nome);
7032 vectorP.RemoveAt(0);
7043 //_____________________________________________________________________________
7044 TObjArray *AliTRDCalibra::ConvertTreeVector(TTree *tree) const
7047 // Convert the branch groupnumber of the tree taken from
7048 // TRD.calibration.root in case of vector method in a std::vector
7053 TObjArray *vectorplace = new TObjArray();
7055 // Variable of the tree
7056 Int_t groupnumber = -1; // Group calibration
7059 tree->SetBranchAddress("groupnumber",&groupnumber);
7062 Int_t ntotal = tree->GetEntries();
7063 for (Int_t k = 0; k < ntotal; k++) {
7065 AliTRDPlace *placegroupnumber = new AliTRDPlace();
7066 placegroupnumber->SetPlace(groupnumber);
7067 vectorplace->Add((TObject *) placegroupnumber);
7074 //_____________________________________________________________________________
7075 Bool_t AliTRDCalibra::MergeVectorCT(TObjArray *vVectorCT2, TObjArray *pPlaCT2)
7078 // Add the two vectors and place the result in the first
7081 if (((Int_t) pPlaCT2->GetEntriesFast()) != ((Int_t) vVectorCT2->GetEntriesFast())){
7082 AliInfo("VectorCT2 doesn't correspond to PlaCT2!");
7087 for (Int_t k = 0; k < (Int_t) fPlaCH->GetEntriesFast(); k++) {
7089 // Look if PlaCT1[k] it is also in the second vector
7091 for (Int_t j = 0; j < (Int_t) pPlaCT2->GetEntriesFast(); j++) {
7092 if (((AliTRDPlace *) pPlaCT2->At(j))->GetPlace() ==
7093 ((AliTRDPlace *) fPlaCH->At(k))->GetPlace()) {
7099 // If not in the second vector nothing to do
7101 // If in the second vector
7104 AliTRDCTInfo *fCTInfo = new AliTRDCTInfo();
7105 UShort_t *entries = new UShort_t[fNumberBinCharge];
7107 for (Int_t nu = 0; nu < fNumberBinCharge; nu++) {
7108 entries[nu] = ((AliTRDCTInfo *) fVectorCH->At(((AliTRDPlace *) fPlaCH->At(k))->GetPlace()))->GetEntries()[nu]
7109 + ((AliTRDCTInfo *) vVectorCT2->At(((AliTRDPlace *) fPlaCH->At(k))->GetPlace()))->GetEntries()[nu];
7113 fCTInfo->SetEntries(entries);
7115 // Nothing to do on PlaCT1
7117 // Update the vector
7118 fVectorCH->AddAt((TObject *) fCTInfo,((AliTRDPlace *) fPlaCH->At(k))->GetPlace());
7124 // And at the end the vector in CT2 but not in CH1
7125 for (Int_t k = 0; k < (Int_t) pPlaCT2->GetEntriesFast(); k++) {
7127 // Look if pPlaCT2[k] it is also in the second vector
7129 for (Int_t j = 0; j < (Int_t) fPlaCH->GetEntriesFast(); j++) {
7130 if (((AliTRDPlace *) fPlaCH->At(j))->GetPlace() == ((AliTRDPlace *) pPlaCT2->At(k))->GetPlace()) {
7136 // If not in the first vector
7139 AliTRDCTInfo *fCTInfo = new AliTRDCTInfo();
7140 fCTInfo = ((AliTRDCTInfo *) vVectorCT2->At(((AliTRDPlace *) pPlaCT2->At(k))->GetPlace()));
7143 fPlaCH->Add((TObject *) (pPlaCT2->At(k)));
7144 fVectorCH->Add((TObject *) fCTInfo);
7154 //_____________________________________________________________________________
7155 Bool_t AliTRDCalibra::MergeVectorP(TObjArray *vVectorP2
7160 // Add the two vectors and place the result in the first
7163 if (((Int_t) pPlaP2->GetEntriesFast()) != ((Int_t) vVectorP2->GetEntriesFast())) {
7164 AliInfo("VectorP2 doesn't correspond to PlaP2!");
7171 for (Int_t k = 0; k < (Int_t) fPlaPH->GetEntriesFast(); k++) {
7173 // Look if fPlaPH[k] it is also in the second vector
7175 for (Int_t j = 0; j < (Int_t) pPlaP2->GetEntriesFast(); j++) {
7176 if (((AliTRDPlace *) pPlaP2->At(j))->GetPlace() == ((AliTRDPlace *) fPlaPH->At(k))->GetPlace()) {
7182 // If not in the second vector nothing to do
7184 // If in the second vector
7187 AliTRDPInfo *fPInfo = new AliTRDPInfo();
7188 UShort_t *entries = new UShort_t[fTimeMax];
7189 Float_t *sum = new Float_t[fTimeMax];
7190 Float_t *sumsquare = new Float_t[fTimeMax];
7192 for (Int_t nu = 0; nu < fTimeMax; nu++) {
7194 entries[nu] = ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]
7195 + ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu];
7197 Double_t calcul = ((((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSum()[nu])
7198 * ((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]))
7199 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSum()[nu])
7200 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu])))
7201 / ((Double_t) fPInfo->GetEntries()[nu]);
7203 sum[nu] = (Float_t) calcul;
7205 Double_t calculsquare = ((((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSumSquare()[nu])
7206 * ((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]))
7207 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSumSquare()[nu])
7208 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu])))
7209 / ((Double_t) fPInfo->GetEntries()[nu]);
7212 sumsquare[nu] = calculsquare;
7217 fPInfo->SetSum(sum);
7218 fPInfo->SetSumSquare(sumsquare);
7219 fPInfo->SetEntries(entries);
7221 // Nothing to do on PlaCT1
7223 // Update the vector VectorCT1
7224 fVectorPH->AddAt((TObject *) fPInfo,((AliTRDPlace *) fPlaPH->At(k))->GetPlace());
7230 // And at the end the vector in P2 but not in CH1
7231 for (Int_t k = 0; k < (Int_t) pPlaP2->GetEntriesFast(); k++) {
7233 // Look if PlaCT2[k] it is also in the second vector
7235 for (Int_t j = 0; j < (Int_t) fPlaPH->GetEntriesFast(); j++) {
7236 if (((AliTRDPlace *) fPlaPH->At(j))->GetPlace() == ((AliTRDPlace *) pPlaP2->At(k))->GetPlace()) {
7242 // If not in the first vector
7245 AliTRDPInfo *fPInfo = new AliTRDPInfo();
7246 fPInfo = (AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) pPlaP2->At(k))->GetPlace());
7248 // Add at the end of CH1
7249 fPlaPH->Add(((TObject *) pPlaP2->At(k)));
7250 fVectorPH->Add((TObject *) fPInfo);
7262 for (Int_t k = 0; k < (Int_t) fPlaPRF->GetEntriesFast(); k++) {
7264 // Look if fPlaPRF[k] it is also in the second vector
7266 for (Int_t j = 0; j < (Int_t) pPlaP2->GetEntriesFast(); j++) {
7267 if (((AliTRDPlace *) pPlaP2->At(j))->GetPlace() == ((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()) {
7273 // If not in the second vector nothing to do
7275 // If in the second vector
7278 AliTRDPInfo *fPInfo = new AliTRDPInfo();
7279 UShort_t *entries = new UShort_t[fNumberBinPRF];
7280 Float_t *sum = new Float_t[fNumberBinPRF];
7281 Float_t *sumsquare = new Float_t[fNumberBinPRF];
7283 for (Int_t nu = 0; nu < fNumberBinPRF; nu++) {
7285 entries[nu] = ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]
7286 + ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu];
7288 Double_t calcul = ((((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSum()[nu])
7289 * ((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]))
7290 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSum()[nu])
7291 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu])))
7292 / ((Double_t) fPInfo->GetEntries()[nu]);
7294 sum[nu] = (Float_t) calcul;
7296 Double_t calculsquare = ((((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSumSquare()[nu])
7297 * ((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]))
7298 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSumSquare()[nu])
7299 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu])))
7300 / ((Double_t) fPInfo->GetEntries()[nu]);
7302 sumsquare[nu] = calculsquare;
7307 fPInfo->SetSum(sum);
7308 fPInfo->SetSumSquare(sumsquare);
7309 fPInfo->SetEntries(entries);
7311 // Nothing to do on PlaCT1
7313 // Update the vector VectorCT1
7314 fVectorPRF->AddAt((TObject *) fPInfo,((AliTRDPlace *) fPlaPRF->At(k))->GetPlace());
7320 // And at the end the vector in P2 but not in CH1
7321 for (Int_t k = 0; k < (Int_t) pPlaP2->GetEntriesFast(); k++) {
7323 // Look if PlaCT2[k] it is also in the second vector
7325 for (Int_t j = 0; j < (Int_t) fPlaPRF->GetEntriesFast(); j++) {
7326 if (((AliTRDPlace *) fPlaPRF->At(j))->GetPlace() == ((AliTRDPlace *) pPlaP2->At(k))->GetPlace()) {
7332 // If not in the first vector
7335 AliTRDPInfo *fPInfo = new AliTRDPInfo();
7336 fPInfo = (AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) pPlaP2->At(k))->GetPlace());
7338 // Add at the end of CH1
7339 fPlaPRF->Add(((TObject *) pPlaP2->At(k)));
7340 fVectorPRF->Add((TObject *) fPInfo);
7352 //____________Fit Methods______________________________________________________
7354 //_____________________________________________________________________________
7355 void AliTRDCalibra::FitPente(TH1* projPH, Int_t idect)
7358 // Slope methode for the drift velocity
7362 const Float_t kDrWidth = AliTRDgeometry::DrThick();
7369 fVdriftCoef[1] = 0.0;
7371 TLine *line = new TLine();
7374 TAxis *xpph = projPH->GetXaxis();
7375 Int_t nbins = xpph->GetNbins();
7376 Double_t lowedge = xpph->GetBinLowEdge(1);
7377 Double_t upedge = xpph->GetBinUpEdge(xpph->GetNbins());
7378 Double_t widbins = (upedge - lowedge) / nbins;
7379 Double_t limit = upedge + 0.5 * widbins;
7381 // Beginning of the signal
7382 TH1D *pentea = new TH1D("pentea","pentea",projPH->GetNbinsX(),0,(Float_t) limit);
7383 for (Int_t k = 1; k < projPH->GetNbinsX(); k++) {
7384 pentea->SetBinContent(k,(Double_t) (projPH->GetBinContent(k+1) - projPH->GetBinContent(k)));
7387 binmax = (Int_t) pentea->GetMaximumBin();
7390 AliInfo("Put the binmax from 1 to 2 to enable the fit");
7392 pentea->Fit("pol2","0MR","",TMath::Max(pentea->GetBinCenter(binmax-1),0.0),pentea->GetBinCenter(binmax+1));
7393 Float_t l3P1am = pentea->GetFunction("pol2")->GetParameter(1);
7394 Float_t l3P2am = pentea->GetFunction("pol2")->GetParameter(2);
7396 fPhd[0] = -(l3P1am / (2 * l3P2am));
7399 // Amplification region
7402 for (Int_t kbin = 1; kbin < projPH->GetNbinsX(); kbin ++) {
7403 if (((projPH->GetBinContent(kbin+1) - projPH->GetBinContent(kbin)) <= 0.0) && (ju == 0)) {
7410 AliInfo("Put the binmax from 1 to 2 to enable the fit");
7412 projPH->Fit("pol2","0MR","",TMath::Max(projPH->GetBinCenter(binmax-1),0.0),projPH->GetBinCenter(binmax+1));
7413 Float_t l3P1amf = projPH->GetFunction("pol2")->GetParameter(1);
7414 Float_t l3P2amf = projPH->GetFunction("pol2")->GetParameter(2);
7417 fPhd[1] = -(l3P1amf / (2 * l3P2amf));
7421 TH1D *pente = new TH1D("pente","pente",projPH->GetNbinsX(),0,(Float_t) limit);
7422 for (Int_t k = binmax+4; k < projPH->GetNbinsX(); k++) {
7423 pente->SetBinContent(k,(Double_t) (projPH->GetBinContent(k+1) - projPH->GetBinContent(k)));
7425 binmin = (Int_t) pente->GetMinimumBin();
7428 AliInfo("Put the binmax from 1 to 2 to enable the fit");
7433 ,TMath::Max(pente->GetBinCenter(binmin-1), 0.0)
7434 ,TMath::Min(pente->GetBinCenter(binmin+2),(Double_t) limit));
7435 Float_t l3P1dr = pente->GetFunction("pol2")->GetParameter(1);
7436 Float_t l3P2dr = pente->GetFunction("pol2")->GetParameter(2);
7438 fPhd[2] = -(l3P1dr / (2 * l3P2dr));
7441 if ((fPhd[2] > fPhd[0]) &&
7442 (fPhd[2] > fPhd[1]) &&
7443 (fPhd[1] > fPhd[0])) {
7444 fVdriftCoef[1] = (kDrWidth) / (fPhd[2]-fPhd[1]);
7445 if (fPhd[0] >= 0.0) {
7446 fT0Coef[1] = (fPhd[0] - fT0Shift) / widbins;
7447 if (fT0Coef[1] < 0.0) {
7448 fT0Coef[1] = - TMath::Abs(fT0Coef[2]);
7452 fT0Coef[1] = -TMath::Abs(fT0Coef[2]);
7456 fVdriftCoef[1] = -TMath::Abs(fVdriftCoef[2]);
7457 fT0Coef[1] = -TMath::Abs(fT0Coef[2]);
7460 if ((fDebug == 1) ||
7462 fCoefVdrift[1]->SetBinContent(idect+1,fVdriftCoef[1]);
7463 fCoefT0[1]->SetBinContent(idect+1,fT0Coef[1]);
7464 if (fVdriftCoef[1] > 0.0) {
7465 if (fVdriftCoef[2] != 0.0) {
7466 fDeltaVdrift[1]->SetBinContent(idect+1,(fVdriftCoef[1] - fVdriftCoef[2]) / fVdriftCoef[2]);
7469 if(fT0Coef[1] >= 0.0){
7470 fDeltaT0[1]->SetBinContent(idect+1,(fT0Coef[1] - fT0Coef[2]));
7475 TCanvas *cpentei = new TCanvas("cpentei","cpentei",50,50,600,800);
7478 line->SetLineColor(2);
7479 line->DrawLine(fPhd[0],0,fPhd[0],projPH->GetMaximum());
7480 line->DrawLine(fPhd[1],0,fPhd[1],projPH->GetMaximum());
7481 line->DrawLine(fPhd[2],0,fPhd[2],projPH->GetMaximum());
7482 AliInfo(Form("fPhd[0] (beginning of the signal): %f" ,(Float_t) fPhd[0]));
7483 AliInfo(Form("fPhd[1] (end of the amplification region): %f" ,(Float_t) fPhd[1]));
7484 AliInfo(Form("fPhd[2] (end of the drift region): %f" ,(Float_t) fPhd[2]));
7485 AliInfo(Form("fVriftCoef[1] (with only the drift region(default)): %f",(Float_t) fVdriftCoef[1]));
7497 //_____________________________________________________________________________
7498 void AliTRDCalibra::FitPH(TH1* projPH, Int_t idect)
7501 // Fit methode for the drift velocity
7505 const Float_t kDrWidth = AliTRDgeometry::DrThick();
7508 TAxis *xpph = projPH->GetXaxis();
7509 Double_t upedge = xpph->GetBinUpEdge(xpph->GetNbins());
7511 TF1 *fPH = new TF1("fPH",AliTRDCalibra::PH,-0.05,3.2,6);
7512 fPH->SetParameter(0,0.469); // Scaling
7513 fPH->SetParameter(1,0.18); // Start
7514 fPH->SetParameter(2,0.0857325); // AR
7515 fPH->SetParameter(3,1.89); // DR
7516 fPH->SetParameter(4,0.08); // QA/QD
7517 fPH->SetParameter(5,0.0); // Baseline
7519 TLine *line = new TLine();
7521 fVdriftCoef[0] = 0.0;
7524 if (idect%fFitPHPeriode == 0) {
7526 AliInfo(Form("<AliTRDCalibra::FitPH> The detector %d will be fitted",idect));
7527 fPH->SetParameter(0,(projPH->Integral()-(projPH->GetBinContent(1)*projPH->GetNbinsX())) * 0.00028); // Scaling
7528 fPH->SetParameter(1,fPhd[0] - 0.1); // Start
7529 fPH->SetParameter(2,fPhd[1] - fPhd[0]); // AR
7530 fPH->SetParameter(3,fPhd[2] - fPhd[1]); // DR
7531 fPH->SetParameter(4,0.225); // QA/QD
7532 fPH->SetParameter(5,(Float_t) projPH->GetBinContent(1));
7535 projPH->Fit(fPH,"0M","",0.0,upedge);
7539 TCanvas *cpente = new TCanvas("cpente","cpente",50,50,600,800);
7541 projPH->Fit(fPH,"M+","",0.0,upedge);
7543 line->SetLineColor(4);
7544 line->DrawLine(fPH->GetParameter(1)
7546 ,fPH->GetParameter(1)
7547 ,projPH->GetMaximum());
7548 line->DrawLine(fPH->GetParameter(1)+fPH->GetParameter(2)
7550 ,fPH->GetParameter(1)+fPH->GetParameter(2)
7551 ,projPH->GetMaximum());
7552 line->DrawLine(fPH->GetParameter(1)+fPH->GetParameter(2)+fPH->GetParameter(3)
7554 ,fPH->GetParameter(1)+fPH->GetParameter(2)+fPH->GetParameter(3)
7555 ,projPH->GetMaximum());
7558 if (fPH->GetParameter(3) != 0) {
7559 fVdriftCoef[0] = kDrWidth / (fPH->GetParameter(3));
7560 fT0Coef[0] = fPH->GetParameter(1);
7563 fVdriftCoef[0] = -TMath::Abs(fVdriftCoef[2]);
7564 fT0Coef[0] = -TMath::Abs(fT0Coef[2]);
7567 if ((fDebug == 1) ||
7569 fCoefVdrift[0]->SetBinContent(idect+1,fVdriftCoef[0]);
7570 fCoefT0[0]->SetBinContent(idect+1,fT0Coef[0]);
7571 if (fVdriftCoef[0] > 0.0){
7572 if (fVdriftCoef[2] != 0.0) {
7573 fDeltaVdrift[0]->SetBinContent(idect+1,(fVdriftCoef[0]-fVdriftCoef[2])/fVdriftCoef[2]);
7575 fDeltaT0[0]->SetBinContent(idect+1,(fT0Coef[0]-fT0Coef[2]));
7579 AliInfo(Form("fVdriftCoef[0]: %f",(Float_t) fVdriftCoef[0]));
7586 // Put the default value
7587 if ((fDebug <= 1) ||
7589 fCoefVdrift[0]->SetBinContent(idect+1,fVdriftCoef[2]);
7590 fCoefT0[0]->SetBinContent(idect+1,fT0Coef[2]);
7601 //_____________________________________________________________________________
7602 void AliTRDCalibra::FitPRF(TH1 *projPRF, Int_t idect)
7605 // Fit methode for the sigma of the pad response function
7611 projPRF->Fit("gaus","0M","",-fRangeFitPRF,fRangeFitPRF);
7615 TCanvas *cfit = new TCanvas("cfit","cfit",50,50,600,800);
7617 projPRF->Fit("gaus","M+","",-fRangeFitPRF,fRangeFitPRF);
7621 fPRFCoef[0] = projPRF->GetFunction("gaus")->GetParameter(2);
7622 if(fPRFCoef[0] <= 0.0) fPRFCoef[0] = -fPRFCoef[1];
7624 if ((fDebug == 1) ||
7626 fCoefPRF[0]->SetBinContent(idect+1,fPRFCoef[0]);
7627 if (fPRFCoef[1] != 0.0) {
7628 fDeltaPRF->SetBinContent(idect+1,(fPRFCoef[0]-fPRFCoef[1])/fPRFCoef[1]);
7632 AliInfo(Form("fPRFCoef[0]: %f",(Float_t) fPRFCoef[0]));
7637 //_____________________________________________________________________________
7638 void AliTRDCalibra::FitCH(TH1 *projch, Int_t idect)
7641 // Fit methode for the gain factor
7644 fChargeCoef[0] = 0.0;
7645 fChargeCoef[1] = 0.0;
7646 TF1 *fLandauGaus = new TF1("fLandauGaus",FuncLandauGaus,0,300,5);
7648 fChargeCoef[1] = projch->GetMean();
7649 projch->Fit("landau","0",""
7650 ,(Float_t) fChargeCoef[1]/fBeginFitCharge
7651 ,projch->GetBinCenter(projch->GetNbinsX()));
7652 fL3P0 = projch->GetFunction("landau")->GetParameter(0);
7653 Double_t l3P1 = projch->GetFunction("landau")->GetParameter(1);
7654 fL3P2 = projch->GetFunction("landau")->GetParameter(2);
7656 projch->Fit("gaus","0",""
7657 ,(Float_t) fChargeCoef[1]/fBeginFitCharge
7658 ,projch->GetBinCenter(projch->GetNbinsX()));
7659 Double_t g3P0 = projch->GetFunction("gaus")->GetParameter(0);
7660 fG3P2 = projch->GetFunction("gaus")->GetParameter(2);
7662 fLandauGaus->SetParameters(fL3P0,l3P1,fL3P2,g3P0,fG3P2);
7663 if ((fDebug <= 1) ||
7665 projch->Fit("fLandauGaus","0",""
7666 ,(Float_t) fChargeCoef[1]/fBeginFitCharge
7667 ,projch->GetBinCenter(projch->GetNbinsX()));
7671 TCanvas *cp = new TCanvas("cp","cp",50,50,600,800);
7673 projch->Fit("fLandauGaus","+",""
7674 ,(Float_t) fChargeCoef[1]/fBeginFitCharge
7675 ,projch->GetBinCenter(projch->GetNbinsX()));
7677 fLandauGaus->Draw("same");
7680 if (projch->GetFunction("fLandauGaus")->GetParameter(1) > 0) {
7681 // Calcul of "real" coef
7682 CalculChargeCoefMean(fCountDet[0],(Int_t) idect,kTRUE);
7683 fChargeCoef[0] = projch->GetFunction("fLandauGaus")->GetParameter(1);
7686 // Calcul of "real" coef
7687 CalculChargeCoefMean(fCountDet[0],(Int_t) idect,kFALSE);
7688 fChargeCoef[0] = -TMath::Abs(fChargeCoef[3]);
7692 AliInfo(Form("fChargeCoef[0]: %f",(Float_t) fChargeCoef[0]));
7693 AliInfo(Form("fChargeCoef[1]: %f",(Float_t) fChargeCoef[1]));
7696 if ((fDebug == 1) ||
7698 if (fChargeCoef[0] > 0.0) {
7699 fCoefCharge[0]->SetBinContent(idect+1,fChargeCoef[0]);
7700 fCoefCharge[1]->SetBinContent(idect+1,fChargeCoef[1]);
7701 fDeltaCharge[0]->SetBinContent(idect+1,fChargeCoef[0]);
7702 fDeltaCharge[1]->SetBinContent(idect+1,fChargeCoef[1]);
7705 fL3P0 = fLandauGaus->Integral(0.3*projch->GetMean(),3*projch->GetMean());
7706 fG3P2 = fLandauGaus->GetParameter(2);
7707 fL3P2 = fLandauGaus->GetParameter(4);
7715 //_____________________________________________________________________________
7716 void AliTRDCalibra::FitBisCH(TH1* projch, Int_t idect)
7719 // Fit methode for the gain factor more time consuming
7722 // Setting fit range and start values
7724 //Double_t sv[4] = { l3P2, fChargeCoef[1], projch->Integral("width"), fG3P2 };
7725 Double_t sv[4] = { fL3P2, fChargeCoef[1], fL3P0, fG3P2 };
7726 Double_t pllo[4] = { 0.001, 0.001, 0.001, 0.001 };
7727 Double_t plhi[4] = { 300.0, 300.0, 100000000.0, 300.0 };
7728 Double_t fp[4] = { 1.0, 1.0, 1.0, 1.0 };
7729 Double_t fpe[4] = { 1.0, 1.0, 1.0, 1.0 };
7730 fr[0] = 0.3 * projch->GetMean();
7731 fr[1] = 3.0 * projch->GetMean();
7732 fChargeCoef[2] = 0.0;
7736 TF1 *fitsnr = LanGauFit(projch,&fr[0],&sv[0]
7741 Double_t projchPeak;
7742 Double_t projchFWHM;
7743 LanGauPro(fp,projchPeak,projchFWHM);
7746 fChargeCoef[2] = fp[1];
7749 fChargeCoef[2] = -TMath::Abs(fChargeCoef[3]);
7753 AliInfo(Form("fChargeCoef[2]: %f",(Float_t) fChargeCoef[2]));
7754 TCanvas *cpy = new TCanvas("cpy","cpy",50,50,600,800);
7757 fitsnr->Draw("same");
7760 if ((fDebug == 1) ||
7762 if (fChargeCoef[2] > 0.0) {
7763 fCoefCharge[2]->SetBinContent(idect+1,fChargeCoef[2]);
7764 fDeltaCharge[2]->SetBinContent(idect+1,fChargeCoef[2]);
7774 //_____________________________________________________________________________
7775 void AliTRDCalibra::NormierungCharge()
7778 // Normalisation of the gain factor resulting for the fits
7781 // Calcul of the mean of the fit
7783 Float_t scalefactor = 1.0;
7784 for (Int_t k = 0; k < (Int_t) fVectorFitCH->GetEntriesFast(); k++) {
7786 Int_t detector = ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetDetector();
7787 Float_t *coef = ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetCoef();
7788 if (GetChamber(detector) == 2) {
7791 if (GetChamber(detector) != 2) {
7794 for (Int_t j = 0; j < total; j++) {
7802 fScaleFitFactor = fScaleFitFactor / sum;
7805 fScaleFitFactor = 1.0;
7809 if ((fDebug == 1) ||
7811 if ((fCoefCharge[0]->GetEntries() > 0.0) &&
7812 (fCoefCharge[0]->GetSumOfWeights() > 0.0)) {
7813 scalefactor = fCoefCharge[0]->GetEntries() / fCoefCharge[0]->GetSumOfWeights();
7814 fCoefCharge[0]->Scale(scalefactor);
7815 fDeltaCharge[0]->Scale(scalefactor);
7817 if ((fMeanChargeOn) &&
7818 (fCoefCharge[1]->GetEntries() > 0.0) &&
7819 (fCoefCharge[1]->GetSumOfWeights() > 0.0)) {
7820 fCoefCharge[1]->Scale(fCoefCharge[1]->GetEntries() / fCoefCharge[1]->GetSumOfWeights());
7822 if ((fFitChargeBisOn) &&
7823 (fCoefCharge[2]->GetEntries() > 0.0) &&
7824 (fCoefCharge[2]->GetSumOfWeights() > 0.0)) {
7825 fCoefCharge[2]->Scale(fCoefCharge[2]->GetEntries() / fCoefCharge[2]->GetSumOfWeights());
7827 if ((fMeanChargeOn) &&
7828 (fDeltaCharge[1]->GetEntries() > 0.0) &&
7829 (fDeltaCharge[1]->GetSumOfWeights() > 0.0)) {
7830 fDeltaCharge[1]->Scale(fDeltaCharge[1]->GetEntries() / fDeltaCharge[1]->GetSumOfWeights());
7832 if ((fFitChargeBisOn) &&
7833 (fDeltaCharge[2]->GetEntries() > 0.0) &&
7834 (fDeltaCharge[2]->GetSumOfWeights() > 0.0)) {
7835 fDeltaCharge[2]->Scale(fDeltaCharge[2]->GetEntries() / fDeltaCharge[2]->GetSumOfWeights());
7839 if ((fDebug == 3) ||
7841 fCoefChargeDB[0]->Scale(scalefactor);
7842 if ((fMeanChargeOn) &&
7843 (fCoefChargeDB[1]->GetEntries() > 0.0) &&
7844 (fCoefChargeDB[1]->GetSumOfWeights() > 0.0)) {
7845 fCoefChargeDB[1]->Scale(fCoefChargeDB[1]->GetEntries() / fCoefChargeDB[1]->GetSumOfWeights());
7847 if ((fFitChargeBisOn) &&
7848 (fCoefChargeDB[2]->GetEntries() > 0.0) &&
7849 (fCoefChargeDB[2]->GetSumOfWeights() > 0.0)) {
7850 fCoefChargeDB[2]->Scale(fCoefChargeDB[2]->GetEntries() / fCoefChargeDB[2]->GetSumOfWeights());
7854 if ((fDebug == 1) ||
7857 fDeltaCharge[0]->Add(fCoefCharge[3],-1);
7858 fDeltaCharge[0]->Divide(fCoefCharge[3]);
7860 if (fMeanChargeOn) {
7861 fDeltaCharge[1]->Add(fCoefCharge[3],-1);
7863 if (fMeanChargeOn) {
7864 fDeltaCharge[1]->Divide(fCoefCharge[3]);
7867 if (fFitChargeBisOn) {
7868 fDeltaCharge[2]->Add(fCoefCharge[3],-1);
7870 if (fFitChargeBisOn) {
7871 fDeltaCharge[2]->Divide(fCoefCharge[3]);
7878 //_____________________________________________________________________________
7879 TH1I *AliTRDCalibra::ReBin(TH1I *hist) const
7882 // Rebin of the 1D histo for the gain calibration if needed.
7883 // you have to choose fRebin, divider of fNumberBinCharge
7886 TAxis *xhist = hist->GetXaxis();
7887 TH1I *rehist = new TH1I("projrebin","",(Int_t) xhist->GetNbins()/fRebin
7888 ,xhist->GetBinLowEdge(1)
7889 ,xhist->GetBinUpEdge(xhist->GetNbins()));
7891 AliInfo(Form("fRebin: %d",fRebin));
7893 for (Int_t k = 1; k <= (Int_t) xhist->GetNbins()/fRebin; k++) {
7895 for (Int_t ji = i; ji < i+fRebin; ji++) {
7896 sum += hist->GetBinContent(ji);
7899 rehist->SetBinContent(k,sum);
7904 TCanvas *crebin = new TCanvas("crebin","",50,50,600,800);
7913 //_____________________________________________________________________________
7914 TH1F *AliTRDCalibra::ReBin(TH1F *hist) const
7917 // Rebin of the 1D histo for the gain calibration if needed
7918 // you have to choose fRebin divider of fNumberBinCharge
7921 TAxis *xhist = hist->GetXaxis();
7922 TH1F *rehist = new TH1F("projrebin","",(Int_t) xhist->GetNbins()/fRebin
7923 ,xhist->GetBinLowEdge(1)
7924 ,xhist->GetBinUpEdge(xhist->GetNbins()));
7926 AliInfo(Form("fRebin: %d",fRebin));
7928 for (Int_t k = 1; k <= (Int_t) xhist->GetNbins()/fRebin; k++) {
7930 for (Int_t ji = i; ji < i+fRebin; ji++) {
7931 sum += hist->GetBinContent(ji);
7934 rehist->SetBinContent(k,sum);
7939 TCanvas *crebin = new TCanvas("crebin","",50,50,600,800);
7948 //_____________________________________________________________________________
7949 TH1F *AliTRDCalibra::CorrectTheError(TGraphErrors *hist)
7952 // In the case of the vectors method the trees contains TGraphErrors for PH and PRF
7953 // to be able to add them after
7954 // We convert it to a TH1F to be able to applied the same fit function method
7955 // After having called this function you can not add the statistics anymore
7960 Int_t nbins = hist->GetN();
7961 Double_t *x = hist->GetX();
7962 Double_t *entries = hist->GetEX();
7963 Double_t *mean = hist->GetY();
7964 Double_t *square = hist->GetEY();
7965 fEntriesCurrent = 0;
7971 Double_t step = x[1] - x[0];
7972 Double_t minvalue = x[0] - step/2;
7973 Double_t maxvalue = x[(nbins-1)] + step/2;
7975 rehist = new TH1F("projcorrecterror","",nbins,minvalue,maxvalue);
7977 for (Int_t k = 0; k < nbins; k++) {
7978 rehist->SetBinContent(k+1,mean[k]);
7979 if (entries[k] > 0.0) {
7980 fEntriesCurrent += (Int_t) entries[k];
7981 Double_t d = TMath::Abs(square[k] - (mean[k]*mean[k]));
7982 rehist->SetBinError(k+1,TMath::Sqrt(d/entries[k]));
7985 rehist->SetBinError(k+1,0.0);
7993 //_____________________________________________________________________________
7994 TGraphErrors *AliTRDCalibra::AddProfiles(TGraphErrors *hist1
7995 , TGraphErrors *hist2) const
7998 // In the case of the vectors method we use TGraphErrors for PH and PRF
7999 // to be able to add the them after
8000 // Here we add the TGraphErrors
8003 // First TGraphErrors
8004 Int_t nbins1 = hist1->GetN();
8005 Double_t *x1 = hist1->GetX();
8006 Double_t *ex1 = hist1->GetEX();
8007 Double_t *y1 = hist1->GetY();
8008 Double_t *ey1 = hist1->GetEY();
8010 TGraphErrors *rehist = new TGraphErrors(nbins1);
8012 // Second TGraphErrors
8013 Double_t *ex2 = hist2->GetEX();
8014 Double_t *y2 = hist2->GetY();
8015 Double_t *ey2 = hist2->GetEY();
8017 // Define the Variables for the new TGraphErrors
8023 for (Int_t k = 0; k < nbins1; k++) {
8024 Double_t nentries = 0.0;
8029 if ((ex2[k] == 0.0) &&
8033 if ((ex2[k] == 0.0) &&
8040 if ((ex2[k] > 0.0) &&
8047 if ((ex2[k] > 0.0) &&
8049 nentries = ex1[k] + ex2[k];
8050 y = ( y1[k]*ex1[k]+ y2[k]*ex2[k]) / nentries;
8051 ey = (ey1[k]*ex1[k]+ey2[k]*ex2[k]) / nentries;
8054 rehist->SetPoint(k,x,y);
8055 rehist->SetPointError(k,ex,ey);
8063 //____________Some basic geometry function_____________________________________
8066 //_____________________________________________________________________________
8067 Int_t AliTRDCalibra::GetPlane(Int_t d) const
8070 // Reconstruct the plane number from the detector number
8073 return ((Int_t) (d % 6));
8077 //_____________________________________________________________________________
8078 Int_t AliTRDCalibra::GetChamber(Int_t d) const
8081 // Reconstruct the chamber number from the detector number
8085 return ((Int_t) (d % 30) / fgkNplan);
8089 //_____________________________________________________________________________
8090 Int_t AliTRDCalibra::GetSector(Int_t d) const
8093 // Reconstruct the sector number from the detector number
8097 return ((Int_t) (d / fg));
8102 //____________Fill and Init tree Gain, PRF, Vdrift and T0______________________
8105 //_____________________________________________________________________________
8106 void AliTRDCalibra::InitTreePRF()
8109 // Init the tree where the coefficients from the fit methods can be stored
8113 fPRFPad = new Float_t[2304];
8114 fPRF = new TTree("PRF","PRF");
8115 fPRF->Branch("detector",&fPRFDetector,"detector/I");
8116 fPRF->Branch("width" ,fPRFPad ,"width[2304]/F");
8118 // Set to default value for the plane 0 supposed to be the first one
8119 for (Int_t k = 0; k < 2304; k++) {
8126 //_____________________________________________________________________________
8127 void AliTRDCalibra::FillTreePRF(Int_t countdet)
8130 // Fill the tree with the sigma of the pad response function for the detector countdet
8133 Int_t numberofgroup = 0;
8134 fPRFDetector = countdet;
8137 if (GetChamber((Int_t)(countdet+1)) == 2) {
8138 numberofgroup = 1728;
8141 numberofgroup = 2304;
8144 // Reset to default value for the next
8145 for (Int_t k = 0; k < numberofgroup; k++) {
8146 if (GetPlane((Int_t) (countdet+1)) == 0) {
8149 if (GetPlane((Int_t) (countdet+1)) == 1) {
8152 if (GetPlane((Int_t) (countdet+1)) == 2) {
8155 if (GetPlane((Int_t) (countdet+1)) == 3) {
8158 if (GetPlane((Int_t) (countdet+1)) == 4) {
8161 if (GetPlane((Int_t) (countdet+1)) == 5) {
8170 //_____________________________________________________________________________
8171 void AliTRDCalibra::ConvertVectorFitCHTree()
8174 // Convert the vector stuff to a tree of 1D histos if the user
8175 // want to write it after the fill functions
8178 Int_t detector = -1;
8179 Int_t numberofgroup = 1;
8180 Float_t gainPad[2304];
8182 fGain = new TTree("Gain","Gain");
8183 fGain->Branch("detector",&detector,"detector/I");
8184 fGain->Branch("gainPad" ,gainPad ,"gainPad[2304]/F");
8186 Int_t loop = (Int_t) fVectorFitCH->GetEntriesFast();
8187 for (Int_t k = 0; k < loop; k++) {
8188 detector = ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetDetector();
8189 if (GetChamber((Int_t) ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetDetector()) == 2) {
8190 numberofgroup = 1728;
8193 numberofgroup = 2304;
8195 for (Int_t i = 0; i < numberofgroup; i++) {
8196 if (((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetCoef()[i] >= 0) {
8197 gainPad[i] = ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetCoef()[i] * fScaleFitFactor;
8200 gainPad[i] = (Float_t) ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetCoef()[i];
8208 //_____________________________________________________________________________
8209 void AliTRDCalibra::FillTreeVdrift(Int_t countdet)
8212 // Fill the tree with the drift velocities for the detector countdet
8215 Int_t numberofgroup = 0;
8216 fVdriftDetector = countdet;
8218 for (Int_t k = 0; k < 2304; k++) {
8222 if (GetChamber((Int_t)(countdet+1)) == 2) {
8223 numberofgroup = 1728;
8226 numberofgroup = 2304;
8228 // Reset to default value the gain coef
8229 for (Int_t k = 0; k < numberofgroup; k++) {
8230 fVdriftPad[k] = -1.5;
8232 fVdriftDetector = -1;
8236 //_____________________________________________________________________________
8237 void AliTRDCalibra::InitTreePH()
8240 // Init the tree where the coefficients from the fit methods can be stored
8244 fVdriftPad = new Float_t[2304];
8245 fVdrift = new TTree("Vdrift","Vdrift");
8246 fVdrift->Branch("detector",&fVdriftDetector,"detector/I");
8247 fVdrift->Branch("vdrift" ,fVdriftPad ,"vdrift[2304]/F");
8248 // Set to default value for the plane 0 supposed to be the first one
8249 for (Int_t k = 0; k < 2304; k++) {
8250 fVdriftPad[k] = -1.5;
8252 fVdriftDetector = -1;
8256 //_____________________________________________________________________________
8257 void AliTRDCalibra::FillTreeT0(Int_t countdet)
8260 // Fill the tree with the t0 value for the detector countdet
8263 Int_t numberofgroup = 0;
8265 fT0Detector = countdet;
8266 for (Int_t k = 0; k < 2304; k++) {
8270 if (GetChamber((Int_t) (countdet+1)) == 2) {
8271 numberofgroup = 1728;
8274 numberofgroup = 2304;
8276 // Reset to default value the gain coef
8277 for (Int_t k = 0; k < numberofgroup; k++) {
8284 //_____________________________________________________________________________
8285 void AliTRDCalibra::InitTreeT0()
8288 // Init the tree where the coefficients from the fit methods can be stored
8292 fT0Pad = new Float_t[2304];
8293 fT0 = new TTree("T0","T0");
8294 fT0->Branch("detector",&fT0Detector,"detector/I");
8295 fT0->Branch("t0",fT0Pad,"t0[2304]/F");
8296 //Set to default value for the plane 0 supposed to be the first one
8297 for(Int_t k = 0; k < 2304; k++){
8305 //____________Private Functions________________________________________________
8308 //_____________________________________________________________________________
8309 Double_t AliTRDCalibra::PH(Double_t *x, Double_t *par)
8312 // Function for the fit
8315 //TF1 *fAsymmGauss = new TF1("fAsymmGauss",AsymmGauss,0,4,6);
8317 //PARAMETERS FOR FIT PH
8319 //fAsymmGauss->SetParameter(0,0.113755);
8320 //fAsymmGauss->SetParameter(1,0.350706);
8321 //fAsymmGauss->SetParameter(2,0.0604244);
8322 //fAsymmGauss->SetParameter(3,7.65596);
8323 //fAsymmGauss->SetParameter(4,1.00124);
8324 //fAsymmGauss->SetParameter(5,0.870597); // No tail cancelation
8332 Double_t dx = 0.005;
8333 Double_t xs = par[1];
8335 Double_t paras[2] = { 0.0, 0.0 };
8338 if ((xs >= par[1]) &&
8339 (xs < (par[1]+par[2]))) {
8340 //fAsymmGauss->SetParameter(0,par[0]);
8341 //fAsymmGauss->SetParameter(1,xs);
8342 //ss += fAsymmGauss->Eval(xx);
8345 ss += AsymmGauss(&xx,paras);
8347 if ((xs >= (par[1]+par[2])) &&
8348 (xs < (par[1]+par[2]+par[3]))) {
8349 //fAsymmGauss->SetParameter(0,par[0]*par[4]);
8350 //fAsymmGauss->SetParameter(1,xs);
8351 //ss += fAsymmGauss->Eval(xx);
8352 paras[0] = par[0]*par[4];
8354 ss += AsymmGauss(&xx,paras);
8363 //_____________________________________________________________________________
8364 Double_t AliTRDCalibra::AsymmGauss(Double_t *x, Double_t *par)
8367 // Function for the fit
8370 //par[0] = normalization
8378 Double_t par1save = par[1];
8379 //Double_t par2save = par[2];
8380 Double_t par2save = 0.0604244;
8381 //Double_t par3save = par[3];
8382 Double_t par3save = 7.65596;
8383 //Double_t par5save = par[5];
8384 Double_t par5save = 0.870597;
8385 Double_t dx = x[0] - par1save;
8387 Double_t sigma2 = par2save*par2save;
8388 Double_t sqrt2 = TMath::Sqrt(2.0);
8389 Double_t exp1 = par3save * TMath::Exp(-par3save * (dx - 0.5 * par3save * sigma2))
8390 * (1.0 - TMath::Erf((par3save * sigma2 - dx) / (sqrt2 * par2save)));
8391 Double_t exp2 = par5save * TMath::Exp(-par5save * (dx - 0.5 * par5save * sigma2))
8392 * (1.0 - TMath::Erf((par5save * sigma2 - dx) / (sqrt2 * par2save)));
8394 //return par[0]*(exp1+par[4]*exp2);
8395 return par[0] * (exp1 + 1.00124 * exp2);
8399 //_____________________________________________________________________________
8400 Double_t AliTRDCalibra::FuncLandauGaus(Double_t *x, Double_t *par)
8403 // Sum Landau + Gaus with identical mean
8406 Double_t valLandau = par[0] * TMath::Landau(x[0],par[1],par[2]);
8407 //Double_t valGaus = par[3] * TMath::Gaus(x[0],par[4],par[5]);
8408 Double_t valGaus = par[3] * TMath::Gaus(x[0],par[1],par[4]);
8409 Double_t val = valLandau + valGaus;
8415 //_____________________________________________________________________________
8416 Double_t AliTRDCalibra::LanGauFun(Double_t *x, Double_t *par)
8419 // Function for the fit
8422 // par[0]=Width (scale) parameter of Landau density
8423 // par[1]=Most Probable (MP, location) parameter of Landau density
8424 // par[2]=Total area (integral -inf to inf, normalization constant)
8425 // par[3]=Width (sigma) of convoluted Gaussian function
8427 // In the Landau distribution (represented by the CERNLIB approximation),
8428 // the maximum is located at x=-0.22278298 with the location parameter=0.
8429 // This shift is corrected within this function, so that the actual
8430 // maximum is identical to the MP parameter.
8433 // Numeric constants
8434 Double_t invsq2pi = 0.3989422804014; // (2 pi)^(-1/2)
8435 Double_t mpshift = -0.22278298; // Landau maximum location
8437 // Control constants
8438 Double_t np = 100.0; // Number of convolution steps
8439 Double_t sc = 5.0; // Convolution extends to +-sc Gaussian sigmas
8451 // MP shift correction
8452 mpc = par[1] - mpshift * par[0];
8454 // Range of convolution integral
8455 xlow = x[0] - sc * par[3];
8456 xupp = x[0] + sc * par[3];
8458 step = (xupp - xlow) / np;
8460 // Convolution integral of Landau and Gaussian by sum
8461 for (i = 1.0; i <= np/2; i++) {
8463 xx = xlow + (i-.5) * step;
8464 fland = TMath::Landau(xx,mpc,par[0]) / par[0];
8465 sum += fland * TMath::Gaus(x[0],xx,par[3]);
8467 xx = xupp - (i-.5) * step;
8468 fland = TMath::Landau(xx,mpc,par[0]) / par[0];
8469 sum += fland * TMath::Gaus(x[0],xx,par[3]);
8473 return (par[2] * step * sum * invsq2pi / par[3]);
8477 //_____________________________________________________________________________
8478 TF1 *AliTRDCalibra::LanGauFit(TH1 *his, Double_t *fitrange, Double_t *startvalues
8479 , Double_t *parlimitslo, Double_t *parlimitshi
8480 , Double_t *fitparams, Double_t *fiterrors
8481 , Double_t *chiSqr, Int_t *ndf)
8484 // Function for the fit
8488 Char_t funname[100];
8490 AliInfo(Form(funname,"Fitfcn_%s",his->GetName()));
8492 TF1 *ffitold = (TF1 *) gROOT->GetListOfFunctions()->FindObject(funname);
8497 TF1 *ffit = new TF1(funname,LanGauFun,fitrange[0],fitrange[1],4);
8498 ffit->SetParameters(startvalues);
8499 ffit->SetParNames("Width","MP","Area","GSigma");
8501 for (i = 0; i < 4; i++) {
8502 ffit->SetParLimits(i,parlimitslo[i],parlimitshi[i]);
8505 his->Fit(funname,"RB0"); // Fit within specified range, use ParLimits, do not plot
8507 ffit->GetParameters(fitparams); // Obtain fit parameters
8508 for (i = 0; i < 4; i++) {
8509 fiterrors[i] = ffit->GetParError(i); // Obtain fit parameter errors
8511 chiSqr[0] = ffit->GetChisquare(); // Obtain chi^2
8512 ndf[0] = ffit->GetNDF(); // Obtain ndf
8514 return (ffit); // Return fit function
8518 //_____________________________________________________________________________
8519 Int_t AliTRDCalibra::LanGauPro(Double_t *params, Double_t &maxx, Double_t &fwhm)
8522 // Function for the fit
8535 Int_t maxcalls = 10000;
8537 // Search for maximum
8538 p = params[1] - 0.1 * params[0];
8539 step = 0.05 * params[0];
8543 while ((l != lold) && (i < maxcalls)) {
8547 l = LanGauFun(&x,params);
8549 step = -step / 10.0;
8554 if (i == maxcalls) {
8560 // Search for right x location of fy
8561 p = maxx + params[0];
8567 while ( (l != lold) && (i < maxcalls) ) {
8572 l = TMath::Abs(LanGauFun(&x,params) - fy);
8586 // Search for left x location of fy
8588 p = maxx - 0.5 * params[0];
8594 while ((l != lold) && (i < maxcalls)) {
8598 l = TMath::Abs(LanGauFun(&x,params) - fy);
8600 step = -step / 10.0;
8605 if (i == maxcalls) {