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->IsPadMasked(detector,col,row)) {
2537 //____________Functions for plotting the 2D____________________________________
2539 //_____________________________________________________________________________
2540 void AliTRDCalibra::Plot2d()
2543 // Plot the 2D histos
2558 //____________Writing the 2D___________________________________________________
2560 //_____________________________________________________________________________
2561 Bool_t AliTRDCalibra::Write2d()
2564 // Write the 2D histograms or the vectors converted in trees in the file
2565 // "TRD.calibration.root"
2568 TFile *fout = TFile::Open(fWriteName,"RECREATE");
2569 // Check if the file could be opened
2570 if (!fout || !fout->IsOpen()) {
2571 AliInfo("No File found!");
2574 AliInfo(Form("Numbertrack: %d Numberusedch[0]: %d, Numberusedch[1]: %d Numberusedph[0]: %d, Numberusedph[1]: %d"
2579 ,fNumberUsedPh[1]));
2581 TStopwatch stopwatch;
2585 if ((fCH2dOn ) && (fWrite[0])) {
2587 fout->WriteTObject(fCH2d);
2594 TTree *treeCH2d = ConvertVectorCTTreeHisto(fVectorCH,fPlaCH,"treeCH2d",(const char *) name);
2595 fout->WriteTObject(treeCH2d);
2598 if ((fPH2dOn ) && (fWrite[1])) {
2600 fout->WriteTObject(fPH2d);
2607 TTree *treePH2d = ConvertVectorPTreeHisto(fVectorPH,fPlaPH,"treePH2d",(const char *) name);
2608 fout->WriteTObject(treePH2d);
2611 if ((fPRF2dOn ) && (fWrite[2])) {
2613 fout->WriteTObject(fPRF2d);
2620 TTree *treePRF2d = ConvertVectorPTreeHisto(fVectorPRF,fPlaPRF,"treePRF2d",(const char *) name);
2621 fout->WriteTObject(treePRF2d);
2627 AliInfo(Form("Execution time Write2d: R:%.2fs C:%.2fs"
2628 ,stopwatch.RealTime(),stopwatch.CpuTime()));
2634 //_____________________________________________________________________________
2635 AliTRDCalDet *AliTRDCalibra::CreateDetObjectTree(TTree *tree, Int_t i)
2638 // It creates the AliTRDCalDet object from the tree of the coefficient
2639 // for the calibration i (i != 2)
2640 // It takes the mean value of the coefficients per detector
2641 // This object has to be written in the database
2644 // Create the DetObject
2645 AliTRDCalDet *object = 0x0;
2647 object = new AliTRDCalDet("ChamberGainFactor","GainFactor (detector value)");
2650 object = new AliTRDCalDet("ChamberVdrift","TRD drift velocities (detector value)");
2653 object = new AliTRDCalDet("ChamberT0","T0 (detector value)");
2657 Int_t detector = -1;
2658 Float_t values[2304];
2659 tree->SetBranchAddress("detector",&detector);
2661 tree->SetBranchAddress("gainPad",values);
2664 tree->SetBranchAddress("vdrift" ,values);
2667 tree->SetBranchAddress("t0" ,values);
2670 // For calculating the mean
2673 Int_t numberofentries = tree->GetEntries();
2675 if (numberofentries != 540) {
2676 AliInfo("The tree is not complete");
2679 for (Int_t det = 0; det < numberofentries; ++det) {
2680 tree->GetEntry(det);
2681 if (GetChamber(detector) == 2) {
2688 for (Int_t k = 0; k < nto; k++) {
2689 mean += TMath::Abs(values[k]) / nto;
2691 object->SetValue(detector,mean);
2698 //_____________________________________________________________________________
2699 TObject *AliTRDCalibra::CreatePadObjectTree(TTree *tree, Int_t i
2700 , AliTRDCalDet *detobject)
2703 // It Creates the AliTRDCalPad object from the tree of the
2704 // coefficient for the calibration i (i != 2)
2705 // You need first to create the object for the detectors,
2706 // where the mean value is put.
2707 // This object has to be written in the database
2710 // Create the DetObject
2711 AliTRDCalPad *object = 0x0;
2713 object = new AliTRDCalPad("GainFactor","GainFactor (local variations)");
2716 object = new AliTRDCalPad("LocalVdrift","TRD drift velocities (local variations)");
2719 object = new AliTRDCalPad("LocalT0","T0 (local variations)");
2723 Int_t detector = -1;
2724 Float_t values[2304];
2725 tree->SetBranchAddress("detector",&detector);
2727 tree->SetBranchAddress("gainPad",values);
2730 tree->SetBranchAddress("vdrift" ,values);
2733 tree->SetBranchAddress("t0" ,values);
2738 Int_t numberofentries = tree->GetEntries();
2740 if (numberofentries != 540) {
2741 AliInfo("The tree is not complete");
2744 for (Int_t det = 0; det < numberofentries; ++det) {
2745 tree->GetEntry(det);
2746 AliTRDCalROC *calROC = object->GetCalROC(detector);
2747 mean = detobject->GetValue(detector);
2751 Int_t rowMax = calROC->GetNrows();
2752 Int_t colMax = calROC->GetNcols();
2753 for (Int_t row = 0; row < rowMax; ++row) {
2754 for (Int_t col = 0; col < colMax; ++col) {
2756 calROC->SetValue(col,row,TMath::Abs(values[(Int_t) (col*rowMax+row)])/mean);
2766 //_____________________________________________________________________________
2767 TObject *AliTRDCalibra::CreatePadObjectTree(TTree *tree)
2770 // It Creates the AliTRDCalPad object from the tree of the
2771 // coefficient for the calibration PRF (i = 2)
2772 // This object has to be written in the database
2775 // Create the DetObject
2776 AliTRDCalPad *object = new AliTRDCalPad("PRFWidth","PRFWidth");
2779 Int_t detector = -1;
2780 Float_t values[2304];
2781 tree->SetBranchAddress("detector",&detector);
2782 tree->SetBranchAddress("width" ,values);
2785 Int_t numberofentries = tree->GetEntries();
2787 if (numberofentries != 540) {
2788 AliInfo("The tree is not complete");
2791 for (Int_t det = 0; det < numberofentries; ++det) {
2792 tree->GetEntry(det);
2793 AliTRDCalROC *calROC = object->GetCalROC(detector);
2794 Int_t rowMax = calROC->GetNrows();
2795 Int_t colMax = calROC->GetNcols();
2796 for (Int_t row = 0; row < rowMax; ++row) {
2797 for (Int_t col = 0; col < colMax; ++col) {
2798 calROC->SetValue(col,row,TMath::Abs(values[(Int_t) (col*rowMax+row)]));
2807 //_____________________________________________________________________________
2808 void AliTRDCalibra::SetRelativeScale(Float_t RelativeScale)
2811 // Set the factor that will divide the deposited charge
2812 // to fit in the histo range [0,300]
2815 if (RelativeScale > 0.0) {
2816 fRelativeScale = RelativeScale;
2819 AliInfo("RelativeScale must be strict positif!");
2824 //_____________________________________________________________________________
2825 void AliTRDCalibra::SetNz(Int_t i, Short_t Nz)
2828 // Set the mode of calibration group in the z direction for the parameter i
2836 AliInfo("You have to choose between 0 and 4");
2841 //_____________________________________________________________________________
2842 void AliTRDCalibra::SetNrphi(Int_t i, Short_t Nrphi)
2845 // Set the mode of calibration group in the rphi direction for the parameter i
2853 AliInfo("You have to choose between 0 and 6");
2858 //_____________________________________________________________________________
2859 void AliTRDCalibra::SetPeriodeFitPH(Int_t periodeFitPH)
2862 // Set FitPH if 1 then each detector will be fitted
2865 if (periodeFitPH > 0) {
2866 fFitPHPeriode = periodeFitPH;
2869 AliInfo("periodeFitPH must be higher than 0!");
2874 //_____________________________________________________________________________
2875 void AliTRDCalibra::SetBeginFitCharge(Float_t beginFitCharge)
2878 // The fit of the deposited charge distribution begins at
2879 // histo->Mean()/beginFitCharge
2880 // You can here set beginFitCharge
2883 if (beginFitCharge > 0) {
2884 fBeginFitCharge = beginFitCharge;
2887 AliInfo("beginFitCharge must be strict positif!");
2892 //_____________________________________________________________________________
2893 void AliTRDCalibra::SetT0Shift(Float_t t0Shift)
2896 // The t0 calculated with the maximum positif slope is shift from t0Shift
2897 // You can here set t0Shift
2904 AliInfo("t0Shift must be strict positif!");
2909 //_____________________________________________________________________________
2910 void AliTRDCalibra::SetRangeFitPRF(Float_t rangeFitPRF)
2913 // The fit of the PRF is from -rangeFitPRF to rangeFitPRF
2914 // You can here set rangeFitPRF
2917 if ((rangeFitPRF > 0) &&
2918 (rangeFitPRF <= 1.0)) {
2919 fRangeFitPRF = rangeFitPRF;
2922 AliInfo("rangeFitPRF must be between 0 and 1.0");
2927 //_____________________________________________________________________________
2928 void AliTRDCalibra::SetRebin(Short_t rebin)
2931 // Rebin with rebin time less bins the Ch histo
2932 // You can set here rebin that should divide the number of bins of CH histo
2937 AliInfo("You have to be sure that fRebin divides fNumberBinCharge used!");
2940 AliInfo("You have to choose a positiv value!");
2945 //_____________________________________________________________________________
2946 TTree *AliTRDCalibra::Sum2Trees(const Char_t *filename1
2947 , const Char_t *filename2
2948 , const Char_t *variablecali)
2951 // It returns the sum of two trees with the name variablecali
2952 // in the files filenam1 and filename2 equivalent of merging two 2D histos
2953 // The name of the resulting tree is the same as the two input trees
2954 // variablecali can be treeCH2d, treePH2d or treePRF2d
2958 TChain *treeChain = new TChain(variablecali);
2959 TObjArray *vectorplace = new TObjArray();
2960 TObjArray *where = new TObjArray();
2964 TFile *file1 = new TFile(filename1,"READ");
2965 TTree *tree1 = (TTree *) file1->Get(variablecali);
2970 vectorplace = ConvertTreeVector(tree1);
2972 // Say where it is in tree 1
2973 for (Int_t jui = 0; jui < (Int_t) vectorplace->GetEntriesFast(); jui++) {
2974 AliTRDPlace *placejui = new AliTRDPlace();
2975 placejui->SetPlace(jui);
2976 TObjArray *chainplace = new TObjArray();
2977 chainplace->Add((TObject *) placejui);
2978 where->Add((TObject *) chainplace);
2982 treeChain->Add(filename1);
2987 TFile *file2 = new TFile(filename2,"READ");
2988 TTree *tree2 = (TTree *) file2->Get(variablecali);
2993 TObjArray *vector2 = ConvertTreeVector(tree2);
2994 Int_t j = treeChain->GetEntries();
2996 for (Int_t jui = 0; jui < (Int_t) vector2->GetEntriesFast(); jui++) {
2997 // Search if already found
2998 Int_t place = SearchInTreeVector(vectorplace,((AliTRDPlace *) vector2->At(jui))->GetPlace());
2999 // Create a new element in the two std vectors
3001 AliTRDPlace *placejjui = new AliTRDPlace();
3002 placejjui->SetPlace((j+jui));
3003 TObjArray *chainplace = new TObjArray();
3004 chainplace->Add((TObject *) placejjui);
3005 vectorplace->Add((TObject *) (vector2->At(jui)));
3006 where->Add((TObject *) chainplace);
3008 // Update the element at the place "place" in the std vector whereinthechain
3010 AliTRDPlace *placejjui = new AliTRDPlace();
3011 placejjui->SetPlace((j+jui));
3012 TObjArray *chainplace = ((TObjArray *) where->At(place));
3013 chainplace->Add((TObject *) placejjui);
3014 where->AddAt((TObject *) chainplace,place);
3019 treeChain->Add(filename2);
3022 // Take care of the profile
3023 const Char_t *pattern = "P";
3026 if (!strstr(variablecali,pattern)) {
3028 // Ready to read the chain
3030 treeChain->SetBranchAddress("histo",&his);
3032 // Initialise the final tree
3034 TH1F *histsum = 0x0;
3036 tree = new TTree(variablecali,variablecali);
3037 tree->Branch("groupnumber",&group,"groupnumber/I");
3038 tree->Branch("histo","TH1F",&histsum,32000,0);
3041 if (treeChain->GetEntries() < 1) {
3045 for (Int_t h = 0; h < (Int_t) vectorplace->GetEntriesFast(); h++) {
3046 group = ((AliTRDPlace *) vectorplace->At(h))->GetPlace();
3047 TObjArray *chainplace = ((TObjArray *) where->At(h));
3048 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(0))->GetPlace());
3049 //Init for the first time
3051 histsum = new TH1F("","",his->GetXaxis()->GetNbins()
3052 ,his->GetXaxis()->GetBinLowEdge(1)
3053 ,his->GetXaxis()->GetBinUpEdge(his->GetXaxis()->GetNbins()));
3056 // Reset for each new group
3057 histsum->SetEntries(0.0);
3058 for (Int_t l = 0; l <= histsum->GetXaxis()->GetNbins(); l++) {
3059 histsum->SetBinContent(l,0.0);
3060 histsum->SetBinError(l,0.0);
3062 histsum->Add(his,1);
3063 if ((Int_t) chainplace->GetEntriesFast() > 1) {
3064 for (Int_t s = 1; s < (Int_t) chainplace->GetEntriesFast(); s++) {
3065 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(s))->GetPlace());
3066 histsum->Add(his,1);
3075 // Ready to read the chain
3076 TGraphErrors *his = 0x0;
3077 treeChain->SetBranchAddress("histo",&his);
3079 // Initialise the final tree
3081 TGraphErrors *histsum = 0x0;
3082 Double_t *xref = 0x0;
3084 tree = new TTree(variablecali,variablecali);
3085 tree->Branch("groupnumber",&group,"groupnumber/I");
3086 tree->Branch("histo","TGraphErrors",&histsum,32000,0);
3089 if (treeChain->GetEntries() < 1) {
3093 for (Int_t h = 0; h < (Int_t) vectorplace->GetEntriesFast(); h++) {
3095 group = ((AliTRDPlace *) vectorplace->At(h))->GetPlace();
3096 TObjArray *chainplace = ((TObjArray *) where->At(h));
3097 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(0))->GetPlace());
3098 //Init or reset for a new group
3099 Int_t nbins = his->GetN();
3101 x = new Double_t[nbins];
3104 ex = new Double_t[nbins];
3106 y = new Double_t[nbins];
3108 ey = new Double_t[nbins];
3110 for (Int_t lo = 0; lo < nbins; lo++) {
3117 histsum = new TGraphErrors(nbins,x,y,ex,ey);
3120 histsum = AddProfiles(his,histsum);
3121 if ((Int_t) chainplace->GetEntriesFast() > 1) {
3122 for (Int_t s = 1; s < (Int_t) chainplace->GetEntriesFast(); s++) {
3123 treeChain->GetEntry(((AliTRDPlace *) chainplace->At(s))->GetPlace());
3124 histsum = AddProfiles(his,histsum);
3138 //____________Function fill 2D for the moment out of the code__________________
3140 //____________Function fill 2D all objects from digits_________________________
3141 Bool_t AliTRDCalibra::Create2DDiSimOnline(Int_t iev1, Int_t iev2)
3144 // Only for simulations, after the simulation, create the 2D histos
3145 // from the digits stored in the file "TRD.Digits.root"
3146 // Only for CH and PH
3149 const Int_t kNplan = 6;
3150 const Int_t kNcham = 5;
3152 // RunLoader and so on
3154 delete gAlice->GetRunLoader();
3159 AliRunLoader *rl = AliRunLoader::Open("galice.root");
3165 gAlice = rl->GetAliRun();
3170 // Import the Trees for the event nEvent in the file
3171 rl->LoadKinematics();
3175 AliLoader *loader = rl->GetLoader("TRDLoader");
3177 AliInfo("No TRDLLoader found!");
3181 // Get the pointer to the TRD detector
3182 AliTRD *trd = (AliTRD *) gAlice->GetDetector("TRD");
3184 AliInfo("No TRD detector found");
3188 // Get the pointer to the geometry object
3189 AliTRDgeometry *geo;
3191 geo = trd->GetGeometry();
3194 AliInfo("No TRD geometry found");
3199 AliCDBManager *man = AliCDBManager::Instance();
3201 AliInfo("Could not get CDB Manager");
3205 // Get the parameter object
3206 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
3208 AliInfo("Could not get CommonParam");
3211 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
3213 AliInfo("Could not get calibDB");
3218 fTimeMax = cal->GetNumberOfTimeBins();
3219 fSf = (Float_t) cal->GetSamplingFrequency();
3220 if (fRelativeScaleAuto) {
3224 if (fRelativeScale <= 0.0) {
3225 AliInfo("You have to set the relativescale factor per hand!");
3230 // Create the 2D histos corresponding to the pad group calibration mode
3233 AliInfo(Form("We will fill the CH2d histo with the pad calibration mode: Nz %d, and Nrphi %d"
3237 // Calcul the number of Xbins
3239 ModePadCalibration(2,0);
3240 ModePadFragmentation(0,2,0,0);
3241 fDetChamb2[0] = fNfragZ[0] * fNfragRphi[0];
3242 fNtotal[0] += 6 * 18 * fDetChamb2[0];
3243 ModePadCalibration(0,0);
3244 ModePadFragmentation(0,0,0,0);
3245 fDetChamb0[0] = fNfragZ[0] * fNfragRphi[0];
3246 fNtotal[0] += 6 * 4 * 18 * fDetChamb0[0];
3247 AliInfo(Form("Total number of Xbins: %d",fNtotal[0]));
3249 // Create the 2D histo
3251 CreateCH2d(fNtotal[0]);
3254 fVectorCH = new TObjArray();
3255 fPlaCH = new TObjArray();
3262 AliInfo(Form("We will fill the PH2d histo with the pad calibration mode: Nz %d, and Nrphi %d"
3266 // Calcul the number of Xbins
3268 ModePadCalibration(2,1);
3269 ModePadFragmentation(0,2,0,1);
3270 fDetChamb2[1] = fNfragZ[1] * fNfragRphi[1];
3271 fNtotal[1] += 6 * 18 * fDetChamb2[1];
3272 ModePadCalibration(0,1);
3273 ModePadFragmentation(0,0,0,1);
3274 fDetChamb0[1] = fNfragZ[1] * fNfragRphi[1];
3275 fNtotal[1] += 6 * 4 * 18 * fDetChamb0[1];
3276 AliInfo(Form("Total number of Xbins: %d",fNtotal[1]));
3278 // Create the 2D histo
3280 CreatePH2d(fNtotal[1]);
3283 fVectorPH = new TObjArray();
3284 fPlaPH = new TObjArray();
3289 loader->LoadDigits();
3290 AliInfo("LoadDigits ");
3291 AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager();
3293 //iev2 egal to the max if 0
3295 iev2 = rl->GetNumberOfEvents();
3296 AliInfo(Form("Total number of events: %d",iev2));
3300 for (Int_t ievent = iev1; ievent < iev2; ievent++) {
3301 AliInfo(Form("Process event %d",ievent));
3302 rl->GetEvent(ievent);
3303 if (!loader->TreeD()) {
3304 AliInfo("loader Loading Digits ... ");
3305 loader->LoadDigits();
3307 digitsManager->ReadDigits(loader->TreeD());
3308 AliInfo("digitsManager Read Digits Done");
3309 // Read the digits from the file
3310 if (!(digitsManager->ReadDigits(loader->TreeD()))) {
3315 for (Int_t iSect = 0; iSect < 18; iSect++) {
3316 for (Int_t iChamb = 0; iChamb < kNcham; iChamb++) {
3317 for (Int_t iPlane = 0; iPlane < kNplan; iPlane++) {
3319 // A little geometry:
3320 Int_t iDet = geo->GetDetector(iPlane,iChamb,iSect);
3321 Int_t rowMax = parCom->GetRowMax(iPlane,iChamb,iSect);
3322 Int_t colMax = parCom->GetColMax(iPlane);
3324 // Variables for the group
3325 LocalisationDetectorXbins(iDet);
3327 // In the cas of charge
3329 amptotal = new Float_t[fNfragRphi[0]*fNfragZ[0]];
3331 for (Int_t k = 0; k < fNfragRphi[0]*fNfragZ[0]; k++) {
3336 // Loop through the detector pixel
3337 for (Int_t time = 0; time < fTimeMax; time++) {
3338 for (Int_t col = 0; col < colMax; col++) {
3339 for (Int_t row = 0; row < rowMax; row++) {
3341 // Amplitude and position in pad group
3342 AliTRDdigit *digit = digitsManager->GetDigit(row,col,time,iDet);
3343 Int_t amp = digit->GetAmp();
3344 Int_t posr[2] = {0,0};
3345 Int_t posc[2] = {0,0};
3348 posr[0] = (Int_t) row / fNnZ[0];
3351 (fNnRphi[0] != 0)) {
3352 posc[0] = (Int_t) col / fNnRphi[0];
3356 posr[1] = (Int_t) row / fNnZ[1];
3359 (fNnRphi[1] != 0)) {
3360 posc[1] = (Int_t) col / fNnRphi[1];
3365 if (amp < fThresholdDigit) {
3368 amptotal[(Int_t) (posc[0]*fNfragZ[0]+posr[0])] += amp;
3372 fPH2d->Fill((fXbins[1]+posc[1]*fNfragZ[1]+posr[1])+0.5,(Float_t) time/fSf,(Double_t) amp);
3375 UpdateVectorPH((fXbins[1]+posc[1]*fNfragZ[1]+posr[1]),time,(Double_t) amp);
3388 // If automatic scale
3389 if ((fCountRelativeScale < 100) && (fRelativeScaleAuto)) {
3390 // Take only the one zone track
3391 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3392 if ((fCountRelativeScale < 100) && (amptotal[k] > 2.0)) {
3393 fRelativeScale += amptotal[k]*0.014*0.01;
3394 fCountRelativeScale++;
3399 // We fill the CH2d after having scale with the first 100
3400 if ((fCountRelativeScale >= 100) && (fRelativeScaleAuto)) {
3402 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3404 (amptotal[k] > 0.0)) {
3405 fCH2d->Fill(fXbins[0]+k+0.5,amptotal[k]/fRelativeScale);
3408 (amptotal[k] > 0.0)) {
3409 UpdateVectorCH(fXbins[0]+k ,amptotal[k]/fRelativeScale);
3414 // No relative salce
3415 if (!fRelativeScaleAuto) {
3416 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3418 (amptotal[k] > 0.0)) {
3419 fCH2d->Fill(fXbins[0]+k+0.5, amptotal[k]/fRelativeScale);
3422 (amptotal[k] > 0.0)) {
3423 UpdateVectorCH(fXbins[0]+k, amptotal[k]/fRelativeScale);
3436 loader->UnloadDigits();
3441 if (fPH2dOn && fHisto2d) {
3444 if (fCH2dOn && fHisto2d) {
3449 if (fWrite[0] || fWrite[1]) {
3451 TFile *fout = TFile::Open(fWriteName,"RECREATE");
3452 // Check if the file could be opened
3453 if (!fout || !fout->IsOpen()) {
3454 AliInfo("<No File found!");
3458 if (fCH2dOn && fHisto2d && fWrite[0]) {
3459 fout->WriteTObject(fCH2d);
3461 if (fPH2dOn && fHisto2d && fWrite[1]) {
3462 fout->WriteTObject(fPH2d);
3465 if (fVector2d && fCH2dOn && fWrite[0]) {
3470 TTree *treeCH2d = ConvertVectorCTTreeHisto(fVectorCH,fPlaCH,"treeCH2d",(const char *) name);
3471 fout->WriteTObject(treeCH2d);
3474 if (fVector2d && fPH2dOn && fWrite[1]) {
3479 TTree *treePH2d = ConvertVectorPTreeHisto(fVectorPH,fPlaPH,"treePH2d",(const char *) name);
3480 fout->WriteTObject(treePH2d);
3491 //____________Function fill 2D all objects from Raw Data_______________________
3492 Bool_t AliTRDCalibra::Create2DRaDaOnline(Int_t iev1, Int_t iev2)
3495 // After having written the RAW DATA in the current directory, create the
3496 // 2D histos from these RAW DATA
3497 // Only for CH and PH
3500 const Int_t kNplan = 6;
3501 const Int_t kNcham = 5;
3502 TString dirname(".");
3505 AliCDBManager *man = AliCDBManager::Instance();
3507 AliInfo("Could not get CDB Manager");
3511 // Get the parameter object
3512 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
3514 AliInfo("Could not get CommonParam");
3518 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
3520 AliInfo("Could not get calibDB");
3525 fTimeMax = cal->GetNumberOfTimeBins();
3526 fSf = (Float_t) cal->GetSamplingFrequency();
3527 if (fRelativeScaleAuto) {
3531 if (fRelativeScale <= 0.0) {
3532 AliInfo("You have to set the relativescale factor per hand!");
3537 // Create the 2D histo corresponding to the pad group calibration mode
3540 AliInfo(Form("We will fill the CH2d histo with the pad calibration mode: Nz %d, and Nrphi %d"
3544 // Calcul the number of Xbins
3546 ModePadCalibration(2,0);
3547 ModePadFragmentation(0,2,0,0);
3548 fDetChamb2[0] = fNfragZ[0] * fNfragRphi[0];
3549 fNtotal[0] += 6 * 18 * fDetChamb2[0];
3550 ModePadCalibration(0,0);
3551 ModePadFragmentation(0,0,0,0);
3552 fDetChamb0[0] = fNfragZ[0] * fNfragRphi[0];
3553 fNtotal[0] += 6 * 4 * 18 * fDetChamb0[0];
3554 AliInfo(Form("Total number of Xbins: %d",fNtotal[0]));
3556 // Create the 2D histo
3558 CreateCH2d(fNtotal[0]);
3561 fVectorCH = new TObjArray();
3562 fPlaCH = new TObjArray();
3569 AliInfo(Form("We will fill the PH2d histo with the pad calibration mode: Nz %d, and Nrphi %d"
3573 // Calcul the number of Xbins
3575 ModePadCalibration(2,1);
3576 ModePadFragmentation(0,2,0,1);
3577 fDetChamb2[1] = fNfragZ[1] * fNfragRphi[1];
3578 fNtotal[1] += 6 * 18 * fDetChamb2[1];
3579 ModePadCalibration(0,1);
3580 ModePadFragmentation(0,0,0,1);
3581 fDetChamb0[1] = fNfragZ[1] * fNfragRphi[1];
3582 fNtotal[1] += 6 * 4 * 18 * fDetChamb0[1];
3583 AliInfo(Form("Total number of Xbins: %d",fNtotal[1]));
3585 // Create the 2D histo
3587 CreatePH2d(fNtotal[1]);
3590 fVectorPH = new TObjArray();
3591 fPlaPH = new TObjArray();
3596 AliTRDrawData *rawdata = new AliTRDrawData();
3597 AliInfo("AliTRDrawData object created ");
3600 for (Int_t ievent = iev1; ievent < iev2; ievent++) {
3603 AliRawReaderFile *readerfile = new AliRawReaderFile(dirname,ievent);
3605 AliInfo("No readerfile found!");
3609 AliTRDdigitsManager *digitsManager = rawdata->Raw2Digits((AliRawReader *) readerfile);
3610 if (!digitsManager) {
3611 AliInfo("No DigitsManager done!");
3615 // Loop on detectors
3616 for (Int_t iSect = 0; iSect < 18; iSect++) {
3617 for (Int_t iPlane = 0; iPlane < kNplan; iPlane++) {
3618 for (Int_t iChamb = 0; iChamb < kNcham; iChamb++) {
3620 // A little geometry:
3621 Int_t iDet = AliTRDgeometry::GetDetector(iPlane,iChamb,iSect);
3622 Int_t rowMax = parCom->GetRowMax(iPlane,iChamb,iSect);
3623 Int_t colMax = parCom->GetColMax(iPlane);
3625 // Variables for the group
3626 LocalisationDetectorXbins(iDet);
3628 // In the cas of charge
3630 amptotal = new Float_t[fNfragRphi[0]*fNfragZ[0]];
3632 for (Int_t k = 0; k < fNfragRphi[0]*fNfragZ[0]; k++) {
3637 // Loop through the detector pixel
3638 for (Int_t time = 0; time < fTimeMax; time++) {
3639 for (Int_t col = 0; col < colMax; col++) {
3640 for (Int_t row = 0; row < rowMax; row++) {
3642 // Amplitude and position of the digit
3643 AliTRDdigit *digit = digitsManager->GetDigit(row,col,time,iDet);
3644 Int_t amp = digit->GetAmp();
3645 Int_t posr[2] = { 0, 0 };
3646 Int_t posc[2] = { 0, 0 };
3649 posr[0] = (Int_t) row / fNnZ[0];
3652 (fNnRphi[0] != 0)) {
3653 posc[0] = (Int_t) col / fNnRphi[0];
3657 posr[1] = (Int_t) row / fNnZ[1];
3660 (fNnRphi[1] != 0)) {
3661 posc[1] = (Int_t) col / fNnRphi[1];
3666 if (amp < fThresholdDigit) {
3669 amptotal[(Int_t) (posc[0]*fNfragZ[0]+posr[0])] += amp;
3674 fPH2d->Fill((fXbins[1]+posc[1]*fNfragZ[1]+posr[1])+0.5,(Float_t)time/fSf,amp);
3677 UpdateVectorPH(fXbins[1]+posc[1]*fNfragZ[1]+posr[1],time,amp);
3689 // If automatic scale
3690 if ((fCountRelativeScale < 100) && (fRelativeScaleAuto)) {
3691 // Take only the one zone track
3692 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3693 if ((fCountRelativeScale < 100) && (amptotal[k] > 2.0)) {
3694 fRelativeScale += amptotal[k] * 0.014 * 0.01;
3695 fCountRelativeScale++;
3700 // We fill the CH2d after having scale with the first 100
3701 if ((fCountRelativeScale >= 100) && (fRelativeScaleAuto)) {
3703 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3704 if (fHisto2d && (amptotal[k] > 0.0)) {
3705 fCH2d->Fill(fXbins[0]+k+0.5,amptotal[k]/fRelativeScale);
3707 if (fVector2d && (amptotal[k] > 0.0)) {
3708 UpdateVectorCH(fXbins[0]+k, amptotal[k]/fRelativeScale);
3713 // No relative salce
3714 if (!fRelativeScaleAuto) {
3715 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3717 (amptotal[k] > 0.0)) {
3718 fCH2d->Fill(fXbins[0]+k+0.5,amptotal[k]/fRelativeScale);
3721 (amptotal[k] > 0.0)) {
3722 UpdateVectorCH(fXbins[0]+k, amptotal[k]/fRelativeScale);
3735 delete digitsManager;
3741 if (fPH2dOn && fHisto2d) {
3744 if (fCH2dOn && fHisto2d) {
3749 if (fWrite[0] || fWrite[1]) {
3751 TFile *fout = TFile::Open(fWriteName,"UPDATE");
3752 // Check if the file could be opened
3753 if (!fout || !fout->IsOpen()) {
3754 AliInfo("<No File found!");
3758 if (fCH2dOn && fHisto2d && fWrite[0]) {
3759 fout->WriteTObject(fCH2d);
3761 if (fPH2dOn && fHisto2d && fWrite[1]) {
3762 fout->WriteTObject(fPH2d);
3765 if (fVector2d && fCH2dOn && fWrite[0]) {
3770 TTree *treeCH2d = ConvertVectorCTTreeHisto(fVectorCH,fPlaCH,"treeCH2d",(const Char_t *) name);
3771 fout->WriteTObject(treeCH2d);
3774 if (fVector2d && fPH2dOn && fWrite[1]) {
3779 TTree *treePH2d = ConvertVectorPTreeHisto(fVectorPH,fPlaPH,"treePH2d",(const Char_t *) name);
3780 fout->WriteTObject(treePH2d);
3789 //____________Pad Calibration Public___________________________________________
3791 //____________Define the number of pads per group for one detector and one calibration
3792 void AliTRDCalibra::ModePadCalibration(Int_t iChamb, Int_t i)
3795 // Definition of the calibration mode
3796 // from Nz and Nrphi, the number of row and col pads per calibration groups are setted
3803 if ((fNz[i] == 0) && (iChamb == 2)) {
3806 if ((fNz[i] == 0) && (iChamb != 2)) {
3809 if ((fNz[i] == 1) && (iChamb == 2)) {
3812 if ((fNz[i] == 1) && (iChamb != 2)) {
3815 if ((fNz[i] == 2) && (iChamb == 2)) {
3818 if ((fNz[i] == 2) && (iChamb != 2)) {
3828 if (fNrphi[i] == 0) {
3831 if (fNrphi[i] == 1) {
3834 if (fNrphi[i] == 2) {
3837 if (fNrphi[i] == 3) {
3840 if (fNrphi[i] == 4) {
3843 if (fNrphi[i] == 5) {
3846 if (fNrphi[i] == 6) {
3852 //____________Define the number of pad groups in one detector for one calibration
3853 Bool_t AliTRDCalibra::ModePadFragmentation(Int_t iPlane,Int_t iChamb, Int_t iSect, Int_t i)
3856 // Definition of the calibration mode
3857 // From the number of row and col pads per calibration groups the
3858 // number of calibration groups are setted
3864 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
3866 AliInfo("Could not get CommonParam Manager");
3870 // A little geometry:
3871 Int_t rowMax = parCom->GetRowMax(iPlane,iChamb,iSect);
3872 Int_t colMax = parCom->GetColMax(iPlane);
3874 // The fragmentation
3876 fNfragZ[i] = (Int_t) rowMax / fNnZ[i];
3879 if (fNnRphi[i] != 0) {
3880 fNfragRphi[i] = (Int_t) colMax / fNnRphi[i];
3887 //____________Protected Functions______________________________________________
3888 //____________Create the 2D histo to be filled online__________________________
3891 //_____________________________________________________________________________
3892 void AliTRDCalibra::CreatePRF2d(Int_t nn)
3895 // Create the 2D histos
3903 fPRF2d = new TProfile2D("PRF2d",(const Char_t *) name
3904 ,nn,0,nn,fNumberBinPRF,-1.0,1.0);
3905 fPRF2d->SetXTitle("Det/pad groups");
3906 fPRF2d->SetYTitle("Position x/W [pad width units]");
3907 fPRF2d->SetZTitle("Q_{i}/Q_{total}");
3908 fPRF2d->SetStats(0);
3912 //_____________________________________________________________________________
3913 void AliTRDCalibra::CreatePH2d(Int_t nn)
3916 // Create the 2D histos
3924 fPH2d = new TProfile2D("PH2d",(const Char_t *) name
3926 ,-0.5/fSf,(Float_t) (fTimeMax-0.5)/fSf);
3927 fPH2d->SetXTitle("Det/pad groups");
3928 fPH2d->SetYTitle("time [#mus]");
3929 fPH2d->SetZTitle("<PH> [a.u.]");
3934 //_____________________________________________________________________________
3935 void AliTRDCalibra::CreateCH2d(Int_t nn)
3938 // Create the 2D histos
3946 fCH2d = new TH2I("CH2d",(const Char_t *) name
3947 ,nn,0,nn,fNumberBinCharge,0,300);
3948 fCH2d->SetXTitle("Det/pad groups");
3949 fCH2d->SetYTitle("charge deposit [a.u]");
3950 fCH2d->SetZTitle("counts");
3956 //____________Offine tracking in the AliTRDtracker_____________________________
3957 void AliTRDCalibra::FillTheInfoOfTheTrackCH()
3960 // For the offline tracking or mcm tracklets
3961 // This function will be called in the functions UpdateHistogram...
3962 // to fill the info of a track for the relativ gain calibration
3965 Int_t nb = 0; // Nombre de zones traversees
3966 Int_t fd = -1; // Premiere zone non nulle
3969 // See if the track goes through different zones
3970 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
3971 if (fAmpTotal[k] > 0.0) {
3979 // If automatic scale
3980 if ((fCountRelativeScale < 100) && (fRelativeScaleAuto)) {
3981 // Take only the one zone track
3983 fRelativeScale += fAmpTotal[fd] * 0.014 * 0.01;
3984 fCountRelativeScale++;
3988 // We fill the CH2d after having scale with the first 100
3989 if ((fCountRelativeScale >= 100) && (fRelativeScaleAuto)) {
3990 // Case of track with only one zone
3993 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
3996 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
3999 // Case of track with two zones
4001 // Two zones voisines sinon rien!
4002 if ((fAmpTotal[fd] > 0.0) &&
4003 (fAmpTotal[fd+1] > 0.0)) {
4004 // One of the two very big
4005 if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+1]) {
4007 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
4010 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4013 if (fAmpTotal[fd+1] > fProcent*fAmpTotal[fd]) {
4015 fCH2d->Fill(fXbins[0]+fd+1.5,fAmpTotal[fd+1]/fRelativeScale);
4018 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd+1]/fRelativeScale);
4025 // Fill with no automatic scale
4026 if (!fRelativeScaleAuto) {
4027 // Case of track with only one zone
4031 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
4034 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4037 // Case of track with two zones
4039 // Two zones voisines sinon rien!
4041 if ((fAmpTotal[fd] > 0.0) &&
4042 (fAmpTotal[fd+1] > 0.0)) {
4043 // One of the two very big
4044 if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+1]) {
4046 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
4049 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4053 if (fAmpTotal[fd+1] > fProcent*fAmpTotal[fd]) {
4055 fCH2d->Fill(fXbins[0]+fd+1.5,fAmpTotal[fd+1]/fRelativeScale);
4058 UpdateVectorCH(fXbins[0]+fd+1,fAmpTotal[fd+1]/fRelativeScale);
4064 if (fNfragZ[0] > 1) {
4065 if (fAmpTotal[fd] > 0.0) {
4066 if ((fd+fNfragZ[0]) < (fNfragZ[0]*fNfragRphi[0])) {
4067 if (fAmpTotal[fd+fNfragZ[0]] > 0.0) {
4068 // One of the two very big
4069 if (fAmpTotal[fd] > fProcent*fAmpTotal[fd+fNfragZ[0]]) {
4071 fCH2d->Fill(fXbins[0]+fd+0.5,fAmpTotal[fd]/fRelativeScale);
4074 UpdateVectorCH(fXbins[0]+fd,fAmpTotal[fd]/fRelativeScale);
4078 if (fAmpTotal[fd+fNfragZ[0]] > fProcent*fAmpTotal[fd]) {
4080 fCH2d->Fill(fXbins[0]+fd+fNfragZ[0]+0.5,fAmpTotal[fd+fNfragZ[0]]/fRelativeScale);
4084 UpdateVectorCH(fXbins[0]+fd+fNfragZ[0],fAmpTotal[fd+fNfragZ[0]]/fRelativeScale);
4097 //____________Offine tracking in the AliTRDtracker_____________________________
4098 void AliTRDCalibra::ResetfVariables()
4101 // Reset values of fAmpTotal, fPHValue and fPHPlace for
4102 // the updateHistogram... functions
4105 // Reset the good track
4108 // Reset the fAmpTotal where we put value
4110 for (Int_t k = 0; k < fNfragZ[0]*fNfragRphi[0]; k++) {
4115 // Reset the fPHValue
4117 for (Int_t k = 0; k < fTimeMax; k++) {
4125 //____________Offine tracking in the AliTRDtracker_____________________________
4126 void AliTRDCalibra::FillTheInfoOfTheTrackPH()
4129 // For the offline tracking or mcm tracklets
4130 // This function will be called in the functions UpdateHistogram...
4131 // to fill the info of a track for the drift velocity calibration
4134 Int_t nb = 1; // Nombre de zones traversees 1, 2 ou plus de 3
4135 Int_t fd1 = -1; // Premiere zone non nulle
4136 Int_t fd2 = -1; // Deuxieme zone non nulle
4137 Int_t k1 = -1; // Debut de la premiere zone
4138 Int_t k2 = -1; // Debut de la seconde zone
4140 // See if the track goes through different zones
4141 for (Int_t k = 0; k < fTimeMax; k++) {
4142 if (fPHValue[k] > 0.0) {
4147 if (fPHPlace[k] != fd1) {
4153 if (fPHPlace[k] != fd2) {
4161 // Case of track with only one zone
4164 for (Int_t i = 0; i < fTimeMax; i++) {
4165 if (fPHValue[i] > 0.0) {
4167 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4170 AliInfo(Form("WRITE nb %d ,place final: %d, fPHPlace[i]: %d, i: %d, fPHValue[i]: %f"
4171 ,nb,fXbins[1]+fPHPlace[i],fPHPlace[i],i,fPHValue[i]));
4174 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4179 // Case of track with two zones
4181 // Two zones voisines sinon rien!
4183 if ((fd1 == fd2+1) ||
4185 // One of the two fast all the think
4186 if (k2 > (k1+fDifference)) {
4188 for (Int_t i = k1; i < k2; i++) {
4189 if (fPHValue[i] > 0.0) {
4191 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4194 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4199 if ((k2+fDifference) < fTimeMax) {
4201 for (Int_t i = k2; i < fTimeMax; i++) {
4202 if (fPHValue[i] > 0.0) {
4204 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4207 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4213 // Two zones voisines sinon rien!
4214 if (fNfragZ[1] > 1) {
4216 if ((fd1+fNfragZ[1]) < (fNfragZ[1]*fNfragRphi[1])) {
4217 if (fd2 == (fd1+fNfragZ[1])) {
4218 // One of the two fast all the think
4219 if (k2 > (k1+fDifference)) {
4221 for (Int_t i = k1; i < k2; i++) {
4222 if (fPHValue[i] > 0.0) {
4224 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4227 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4232 if ((k2+fDifference) < fTimeMax) {
4234 for (Int_t i = k2; i < fTimeMax; i++) {
4235 if (fPHValue[i] > 0.0) {
4237 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4240 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4247 // Two zones voisines sinon rien!
4249 if ((fd1 - fNfragZ[1]) >= 0) {
4250 if (fd2 == (fd1 - fNfragZ[1])) {
4251 // One of the two fast all the think
4252 if (k2 > (k1 + fDifference)) {
4254 for (Int_t i = k1; i < k2; i++) {
4255 if (fPHValue[i] > 0.0) {
4257 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4260 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4265 if ((k2+fDifference) < fTimeMax) {
4267 for (Int_t i = k2; i < fTimeMax; i++) {
4268 if (fPHValue[i] > 0.0) {
4270 fPH2d->Fill((fXbins[1]+fPHPlace[i])+0.5,(Float_t) i/fSf,(Float_t) fPHValue[i]);
4273 UpdateVectorPH(fXbins[1]+fPHPlace[i],i,fPHValue[i]);
4286 //____________Set the pad calibration variables for the detector_______________
4287 Bool_t AliTRDCalibra::LocalisationDetectorXbins(Int_t detector)
4290 // For the detector calcul the first Xbins and set the number of row
4291 // and col pads per calibration groups, the number of calibration
4292 // groups in the detector.
4295 // first Xbins of the detector
4297 CalculXBins(detector,0);
4300 CalculXBins(detector,1);
4303 CalculXBins(detector,2);
4306 // fragmentation of idect
4307 for (Int_t i = 0; i < 3; i++) {
4308 ModePadCalibration((Int_t) GetChamber(detector),i);
4309 ModePadFragmentation((Int_t) GetPlane(detector)
4310 , (Int_t) GetChamber(detector)
4311 , (Int_t) GetSector(detector),i);
4318 //____________Plot the 2D histos filled Online_________________________________
4320 //_____________________________________________________________________________
4321 void AliTRDCalibra::PlotPH2d()
4324 // Plot the 2D histo
4327 TCanvas *cph2d = new TCanvas("cph2d","",50,50,600,800);
4329 fPH2d->Draw("LEGO");
4333 //_____________________________________________________________________________
4334 void AliTRDCalibra::PlotCH2d()
4337 // Plot the 2D histos
4340 TCanvas *cch2d = new TCanvas("cch2d","",50,50,600,800);
4342 fCH2d->Draw("LEGO");
4346 //_____________________________________________________________________________
4347 void AliTRDCalibra::PlotPRF2d()
4350 // Plot the 2D histos
4353 TCanvas *cPRF2d = new TCanvas("cPRF2d","",50,50,600,800);
4355 fPRF2d->Draw("LEGO");
4359 //____________Fit______________________________________________________________
4361 //____________Create histos if fDebug == 1 or fDebug >= 3______________________
4363 //_____________________________________________________________________________
4364 void AliTRDCalibra::CreateFitHistoPH(Int_t nbins, Double_t low, Double_t high)
4367 // Create the histos for fDebug = 1 and fDebug = 4 (Fit functions)
4370 // Histograms to store the coef
4371 fCoefVdrift[0] = new TH1F("coefvdrift0" ,"",nbins,low ,high);
4372 fCoefVdrift[1] = new TH1F("coefvdrift1" ,"",nbins,low ,high);
4373 fCoefVdrift[2] = new TH1F("coefvdrift2" ,"",nbins,low ,high);
4375 // Histograms for Debug
4376 fDeltaVdrift[0] = new TH1F("deltavdrift0","",nbins,low ,high);
4377 fDeltaVdrift[1] = new TH1F("deltavdrift1","",nbins,low ,high);
4378 fErrorVdrift[0] = new TH1I("errorvdrift0","",300 ,-0.5,0.5);
4379 fErrorVdrift[1] = new TH1I("errorvdrift1","",300 ,-0.5,0.5);
4381 fCoefVdrift[0]->SetXTitle("Det/pad groups");
4382 fCoefVdrift[0]->SetYTitle("Vdrift [cm/#mus]");
4383 fCoefVdrift[1]->SetXTitle("Det/pad groups");
4384 fCoefVdrift[1]->SetYTitle("Vdrift [cm/#mus]");
4385 fCoefVdrift[2]->SetXTitle("Det/pad groups");
4386 fCoefVdrift[2]->SetYTitle("Vdrift [cm/#mus]");
4388 fDeltaVdrift[0]->SetXTitle("Det/pad groups");
4389 fDeltaVdrift[0]->SetYTitle("#Deltav/v_{sim}");
4390 fDeltaVdrift[1]->SetXTitle("Det/pad groups");
4391 fDeltaVdrift[1]->SetYTitle("#Deltav/v_{sim}");
4393 fErrorVdrift[0]->SetXTitle("#Deltav/v_{sim}");
4394 fErrorVdrift[0]->SetYTitle("counts");
4395 fErrorVdrift[1]->SetXTitle("#Deltav/v_{sim}");
4396 fErrorVdrift[1]->SetYTitle("counts");
4398 fCoefVdrift[0]->SetStats(0);
4399 fCoefVdrift[1]->SetStats(0);
4400 fCoefVdrift[2]->SetStats(0);
4401 fDeltaVdrift[0]->SetStats(0);
4402 fDeltaVdrift[1]->SetStats(0);
4403 fErrorVdrift[0]->SetStats(0);
4404 fErrorVdrift[1]->SetStats(0);
4406 fCoefVdrift[0]->SetMarkerColor(6);
4407 fCoefVdrift[0]->SetMarkerStyle(26);
4408 fCoefVdrift[0]->SetLineColor(6);
4409 fCoefVdrift[1]->SetMarkerColor(2);
4410 fCoefVdrift[1]->SetMarkerStyle(24);
4411 fCoefVdrift[1]->SetLineColor(2);
4412 fCoefVdrift[2]->SetLineColor(4);
4414 fDeltaVdrift[1]->SetMarkerColor(2);
4415 fDeltaVdrift[1]->SetMarkerStyle(24);
4416 fDeltaVdrift[1]->SetLineColor(2);
4417 fDeltaVdrift[0]->SetMarkerColor(6);
4418 fDeltaVdrift[0]->SetMarkerStyle(26);
4419 fDeltaVdrift[0]->SetLineColor(6);
4421 fErrorVdrift[1]->SetLineColor(2);
4422 fErrorVdrift[1]->SetLineStyle(2);
4423 fErrorVdrift[0]->SetLineColor(6);
4424 fErrorVdrift[0]->SetLineStyle(1);
4428 //_____________________________________________________________________________
4429 void AliTRDCalibra::CreateFitHistoT0(Int_t nbins, Double_t low, Double_t high)
4432 // Create the histos for fDebug = 1 and fDebug = 4 (Fit functions)
4435 // Histograms to store the coef
4436 fCoefT0[0] = new TH1F("coefT00" ,"",nbins,low ,high);
4437 fCoefT0[1] = new TH1F("coefT01" ,"",nbins,low ,high);
4438 fCoefT0[2] = new TH1F("coefT02" ,"",nbins,low ,high);
4440 // Histograms for Debug
4441 fDeltaT0[0] = new TH1F("deltaT00","",nbins,low ,high);
4442 fDeltaT0[1] = new TH1F("deltaT01","",nbins,low ,high);
4443 fErrorT0[0] = new TH1I("errorT00","",100,-0.2,0.2);
4444 fErrorT0[1] = new TH1I("errorT01","",100,-0.2,0.2);
4446 fCoefT0[0]->SetXTitle("Det/pad groups");
4447 fCoefT0[0]->SetYTitle("t0 [timebin]");
4448 fCoefT0[1]->SetXTitle("Det/pad groups");
4449 fCoefT0[1]->SetYTitle("t0 [timebin]");
4450 fCoefT0[2]->SetXTitle("Det/pad groups");
4451 fCoefT0[2]->SetYTitle("t0 [timebin]");
4453 fDeltaT0[0]->SetXTitle("Det/pad groups");
4454 fDeltaT0[0]->SetYTitle("#Deltat0 [timebin]");
4455 fDeltaT0[1]->SetXTitle("Det/pad groups");
4456 fDeltaT0[1]->SetYTitle("#Deltat0 [timebin]");
4458 fErrorT0[0]->SetXTitle("#Deltat0 [timebin]");
4459 fErrorT0[0]->SetYTitle("counts");
4460 fErrorT0[1]->SetXTitle("#Deltat0 [timebin]");
4461 fErrorT0[1]->SetYTitle("counts");
4463 fCoefT0[0]->SetStats(0);
4464 fCoefT0[1]->SetStats(0);
4465 fCoefT0[2]->SetStats(0);
4466 fDeltaT0[0]->SetStats(0);
4467 fDeltaT0[1]->SetStats(0);
4468 fErrorT0[0]->SetStats(0);
4469 fErrorT0[1]->SetStats(0);
4471 fCoefT0[0]->SetMarkerColor(6);
4472 fCoefT0[0]->SetMarkerStyle(26);
4473 fCoefT0[0]->SetLineColor(6);
4474 fCoefT0[1]->SetMarkerColor(2);
4475 fCoefT0[1]->SetMarkerStyle(24);
4476 fCoefT0[1]->SetLineColor(2);
4477 fCoefT0[2]->SetLineColor(4);
4479 fDeltaT0[1]->SetMarkerColor(2);
4480 fDeltaT0[1]->SetMarkerStyle(24);
4481 fDeltaT0[1]->SetLineColor(2);
4482 fDeltaT0[0]->SetMarkerColor(6);
4483 fDeltaT0[0]->SetMarkerStyle(26);
4484 fDeltaT0[0]->SetLineColor(6);
4486 fErrorT0[1]->SetLineColor(2);
4487 fErrorT0[1]->SetLineStyle(2);
4488 fErrorT0[0]->SetLineColor(6);
4489 fErrorT0[0]->SetLineStyle(1);
4493 //_____________________________________________________________________________
4494 void AliTRDCalibra::CreateFitHistoCH(Int_t nbins, Double_t low, Double_t high)
4497 // Create the histos for fDebug = 1 and fDebug = 4 (Fit functions)
4500 // Histograms to store the coef
4501 fCoefCharge[0] = new TH1F("coefcharge0" ,"",nbins,low ,high);
4502 fCoefCharge[1] = new TH1F("coefcharge1" ,"",nbins,low ,high);
4503 fCoefCharge[2] = new TH1F("coefcharge2" ,"",nbins,low ,high);
4504 fCoefCharge[3] = new TH1F("coefcharge3" ,"",nbins,low ,high);
4506 // Histograms for Debug
4507 fDeltaCharge[0] = new TH1F("deltacharge0","",nbins,low ,high);
4508 fDeltaCharge[1] = new TH1F("deltacharge1","",nbins,low ,high);
4509 fDeltaCharge[2] = new TH1F("deltacharge2","",nbins,low ,high);
4511 fErrorCharge[0] = new TH1I("errorcharge0","",100 ,-0.5,0.5);
4512 fErrorCharge[1] = new TH1I("errorcharge1","",100 ,-0.5,0.5);
4513 fErrorCharge[2] = new TH1I("errorcharge2","",100 ,-0.5,0.5);
4515 fCoefCharge[0]->SetXTitle("Det/Pad groups");
4516 fCoefCharge[0]->SetYTitle("gain factor");
4517 fCoefCharge[1]->SetXTitle("Det/Pad groups");
4518 fCoefCharge[1]->SetYTitle("gain factor");
4519 fCoefCharge[2]->SetXTitle("Det/Pad groups");
4520 fCoefCharge[2]->SetYTitle("gain factor");
4521 fCoefCharge[3]->SetXTitle("Det/Pad groups");
4522 fCoefCharge[3]->SetYTitle("gain factor");
4524 fDeltaCharge[0]->SetXTitle("Det/Pad groups");
4525 fDeltaCharge[0]->SetYTitle("#Deltag/g_{sim}");
4526 fDeltaCharge[1]->SetXTitle("Det/Pad groups");
4527 fDeltaCharge[1]->SetYTitle("#Deltag/g_{sim}");
4528 fDeltaCharge[2]->SetXTitle("Det/Pad groups");
4529 fDeltaCharge[2]->SetYTitle("#Deltag/g_{sim}");
4530 fDeltaCharge[0]->SetAxisRange(-0.5,0.5,"Y");
4531 fDeltaCharge[1]->SetAxisRange(-0.5,0.5,"Y");
4532 fDeltaCharge[2]->SetAxisRange(-0.5,0.5,"Y");
4534 fErrorCharge[0]->SetXTitle("#Deltag/g_{sim}");
4535 fErrorCharge[0]->SetYTitle("counts");
4536 fErrorCharge[1]->SetXTitle("#Deltag/g_{sim}");
4537 fErrorCharge[1]->SetYTitle("counts");
4538 fErrorCharge[2]->SetXTitle("#Deltag/g_{sim}");
4539 fErrorCharge[2]->SetYTitle("counts");
4541 fDeltaCharge[1]->SetMarkerColor(2);
4542 fDeltaCharge[1]->SetMarkerStyle(24);
4543 fDeltaCharge[1]->SetLineColor(2);
4544 fErrorCharge[1]->SetLineColor(2);
4545 fErrorCharge[1]->SetLineStyle(2);
4546 fDeltaCharge[2]->SetMarkerColor(8);
4547 fDeltaCharge[2]->SetLineColor(8);
4548 fDeltaCharge[2]->SetMarkerStyle(9);
4549 fErrorCharge[2]->SetLineColor(8);
4550 fErrorCharge[2]->SetLineStyle(5);
4551 fDeltaCharge[0]->SetMarkerColor(6);
4552 fDeltaCharge[0]->SetLineColor(6);
4553 fDeltaCharge[0]->SetMarkerStyle(26);
4554 fErrorCharge[0]->SetLineColor(6);
4555 fErrorCharge[0]->SetLineStyle(1);
4557 fCoefCharge[3]->SetLineColor(4);
4558 fCoefCharge[1]->SetMarkerColor(2);
4559 fCoefCharge[1]->SetLineColor(2);
4560 fCoefCharge[1]->SetMarkerStyle(24);
4561 fCoefCharge[2]->SetMarkerColor(8);
4562 fCoefCharge[2]->SetLineColor(8);
4563 fCoefCharge[2]->SetMarkerStyle(9);
4564 fCoefCharge[0]->SetMarkerColor(6);
4565 fCoefCharge[0]->SetLineColor(6);
4566 fCoefCharge[0]->SetMarkerStyle(26);
4568 fErrorCharge[2]->SetLineWidth(3);
4570 fDeltaCharge[1]->SetStats(0);
4571 fDeltaCharge[2]->SetStats(0);
4572 fDeltaCharge[0]->SetStats(0);
4573 fErrorCharge[1]->SetStats(0);
4574 fErrorCharge[2]->SetStats(0);
4575 fErrorCharge[0]->SetStats(0);
4576 fCoefCharge[1]->SetStats(0);
4577 fCoefCharge[0]->SetStats(0);
4578 fCoefCharge[3]->SetStats(0);
4579 fCoefCharge[2]->SetStats(0);
4583 //_____________________________________________________________________________
4584 void AliTRDCalibra::CreateFitHistoPRF(Int_t nbins, Double_t low, Double_t high)
4587 // Create the histos for fDebug = 1 and fDebug = 4 (Fit functions)
4590 // Histograms to store the coef
4591 fCoefPRF[0] = new TH1F("coefPRF0","",nbins,low ,high);
4592 fCoefPRF[1] = new TH1F("coefPRF1","",nbins,low ,high);
4594 // Histograms for Debug
4595 fDeltaPRF = new TH1F("deltaPRF","",nbins,low ,high);
4596 fErrorPRF = new TH1I("errorPRF","",300,-0.5,0.5);
4598 fDeltaPRF->SetMarkerColor(6);
4599 fDeltaPRF->SetMarkerStyle(26);
4600 fDeltaPRF->SetLineColor(6);
4601 fErrorPRF->SetLineColor(6);
4602 fErrorPRF->SetLineStyle(2);
4604 fCoefPRF[1]->SetLineColor(4);
4605 fCoefPRF[0]->SetMarkerColor(6);
4606 fCoefPRF[0]->SetMarkerStyle(26);
4607 fCoefPRF[0]->SetLineColor(6);
4609 fCoefPRF[0]->SetXTitle("Det/Pad groups");
4610 fCoefPRF[0]->SetYTitle("#sigma_{PRF}");
4611 fCoefPRF[1]->SetXTitle("Det/Pad groups");
4612 fCoefPRF[1]->SetYTitle("#sigma_{PRF}");
4614 fDeltaPRF->SetXTitle("Det/Pad groups");
4615 fDeltaPRF->SetYTitle("#Delta#sigma/#sigma_{sim}");
4617 fErrorPRF->SetXTitle("#Delta#sigma/#sigma_{sim}");
4618 fErrorPRF->SetYTitle("counts");
4620 fDeltaPRF->SetStats(0);
4621 fErrorPRF->SetStats(0);
4622 fCoefPRF[1]->SetStats(0);
4623 fCoefPRF[0]->SetStats(0);
4627 //_____________________________________________________________________________
4628 void AliTRDCalibra::CreateFitHistoPRFDB(Int_t rowMax, Int_t colMax)
4631 // Create the histos for fDebug = 3 and fDebug = 4 (Fit functions)
4634 fCoefPRFDB = new TH2F("coefPRF","",rowMax,0,rowMax,colMax,0,colMax);
4636 fCoefPRFDB->SetStats(0);
4637 fCoefPRFDB->SetXTitle("row Number");
4638 fCoefPRFDB->SetYTitle("col Number");
4639 fCoefPRFDB->SetZTitle("PRF width [pad width units]");
4641 fCoefPRFDB->SetFillColor(6);
4642 fCoefPRFDB->SetLineColor(6);
4646 //_____________________________________________________________________________
4647 void AliTRDCalibra::CreateFitHistoCHDB(Int_t rowMax, Int_t colMax)
4650 // Create the histos for fDebug = 3 and fDebug = 4 (Fit functions)
4653 fCoefChargeDB[0] = new TH2F("coefchargedb0","",rowMax,0,rowMax,colMax,0,colMax);
4654 fCoefChargeDB[1] = new TH2F("coefchargedb1","",rowMax,0,rowMax,colMax,0,colMax);
4655 fCoefChargeDB[2] = new TH2F("coefchargedb2","",rowMax,0,rowMax,colMax,0,colMax);
4657 fCoefChargeDB[0]->SetStats(0);
4658 fCoefChargeDB[1]->SetStats(0);
4659 fCoefChargeDB[2]->SetStats(0);
4660 fCoefChargeDB[0]->SetXTitle("row Number");
4661 fCoefChargeDB[0]->SetYTitle("col Number");
4662 fCoefChargeDB[1]->SetXTitle("row Number");
4663 fCoefChargeDB[1]->SetYTitle("col Number");
4664 fCoefChargeDB[2]->SetXTitle("row Number");
4665 fCoefChargeDB[2]->SetYTitle("col Number");
4666 fCoefChargeDB[0]->SetZTitle("f_{g} Fit method");
4667 fCoefChargeDB[1]->SetZTitle("f_{g} Mean method");
4668 fCoefChargeDB[2]->SetZTitle("f_{g} Fitbis method");
4670 fCoefChargeDB[0]->SetFillColor(6);
4671 fCoefChargeDB[0]->SetLineColor(6);
4672 fCoefChargeDB[0]->SetLineColor(6);
4673 fCoefChargeDB[1]->SetFillColor(2);
4674 fCoefChargeDB[1]->SetLineColor(2);
4675 fCoefChargeDB[1]->SetLineColor(2);
4676 fCoefChargeDB[2]->SetFillColor(8);
4677 fCoefChargeDB[2]->SetLineColor(8);
4678 fCoefChargeDB[2]->SetLineColor(8);
4682 //_____________________________________________________________________________
4683 void AliTRDCalibra::CreateFitHistoPHDB(Int_t rowMax, Int_t colMax)
4686 // Create the histos for fDebug = 3 and fDebug = 4 (Fit functions)
4689 fCoefVdriftDB[0] = new TH2F("coefvdriftdb0","",rowMax,0,rowMax,colMax,0,colMax);
4690 fCoefVdriftDB[1] = new TH2F("coefvdriftdb1","",rowMax,0,rowMax,colMax,0,colMax);
4692 fCoefVdriftDB[0]->SetStats(0);
4693 fCoefVdriftDB[1]->SetStats(0);
4694 fCoefVdriftDB[0]->SetXTitle("row Number");
4695 fCoefVdriftDB[0]->SetYTitle("col Number");
4696 fCoefVdriftDB[1]->SetXTitle("row Number");
4697 fCoefVdriftDB[1]->SetYTitle("col Number");
4698 fCoefVdriftDB[0]->SetZTitle("v_{drift} Fit method");
4699 fCoefVdriftDB[1]->SetZTitle("v_{drift} slope method");
4701 fCoefVdriftDB[0]->SetFillColor(6);
4702 fCoefVdriftDB[0]->SetLineColor(6);
4703 fCoefVdriftDB[0]->SetLineColor(6);
4704 fCoefVdriftDB[1]->SetFillColor(2);
4705 fCoefVdriftDB[1]->SetLineColor(2);
4706 fCoefVdriftDB[1]->SetLineColor(2);
4710 //_____________________________________________________________________________
4711 void AliTRDCalibra::CreateFitHistoT0DB(Int_t rowMax, Int_t colMax)
4714 // Create the histos for fDebug = 3 and fDebug = 4 (Fit functions)
4717 fCoefT0DB[0] = new TH2F("coefT0db0","",rowMax,0,rowMax,colMax,0,colMax);
4718 fCoefT0DB[1] = new TH2F("coefT0db1","",rowMax,0,rowMax,colMax,0,colMax);
4720 fCoefT0DB[0]->SetStats(0);
4721 fCoefT0DB[1]->SetStats(0);
4722 fCoefT0DB[0]->SetXTitle("row Number");
4723 fCoefT0DB[0]->SetYTitle("col Number");
4724 fCoefT0DB[1]->SetXTitle("row Number");
4725 fCoefT0DB[1]->SetYTitle("col Number");
4726 fCoefT0DB[0]->SetZTitle("t0 Fit method");
4727 fCoefT0DB[1]->SetZTitle("t0 slope method");
4729 fCoefT0DB[0]->SetFillColor(6);
4730 fCoefT0DB[0]->SetLineColor(6);
4731 fCoefT0DB[0]->SetLineColor(6);
4732 fCoefT0DB[1]->SetFillColor(2);
4733 fCoefT0DB[1]->SetLineColor(2);
4734 fCoefT0DB[1]->SetLineColor(2);
4738 //_____________________________________________________________________________
4739 Bool_t AliTRDCalibra::FillVectorFitCH(Int_t countdet)
4742 // For the Fit functions fill the vector FitCH special for the gain calibration
4745 AliTRDFitCHInfo *fitCHInfo = new AliTRDFitCHInfo();
4748 if (GetChamber(countdet) == 2) {
4755 Float_t *coef = new Float_t[ntotal];
4756 for (Int_t i = 0; i < ntotal; i++) {
4757 coef[i] = fCoefCH[i];
4760 Int_t detector = countdet;
4762 fitCHInfo->SetCoef(coef);
4763 fitCHInfo->SetDetector(detector);
4764 fVectorFitCH->Add((TObject *) fitCHInfo);
4770 //____________Functions for initialising the AliTRDCalibra in the code_________
4771 Bool_t AliTRDCalibra::InitFit(Int_t nbins, Double_t lowedge
4772 , Double_t upedge, Int_t i)
4775 // Init the calibration mode (Nz, Nrphi), the histograms for
4776 // debugging the fit methods if fDebug > 0,
4779 gStyle->SetPalette(1);
4780 gStyle->SetOptStat(1111);
4781 gStyle->SetPadBorderMode(0);
4782 gStyle->SetCanvasColor(10);
4783 gStyle->SetPadLeftMargin(0.13);
4784 gStyle->SetPadRightMargin(0.01);
4786 // Get the parameter object
4787 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
4789 AliInfo("Could not get CommonParam");
4793 // Mode groups of pads: the total number of bins!
4794 Int_t numberofbinsexpected = 0;
4795 ModePadCalibration(2,i);
4796 ModePadFragmentation(0,2,0,i);
4797 fDetChamb2[i] = fNfragZ[i] * fNfragRphi[i];
4799 AliInfo(Form("For the chamber 2: %d",fDetChamb2[i]));
4801 numberofbinsexpected += 6 * 18 * fDetChamb2[i];
4802 ModePadCalibration(0,i);
4803 ModePadFragmentation(0,0,0,i);
4804 fDetChamb0[i] = fNfragZ[i] * fNfragRphi[i];
4806 AliInfo(Form("For the other chamber 0: %d",fDetChamb0[i]));
4808 numberofbinsexpected += 6 * 4 * 18 * fDetChamb0[i];
4810 // Quick verification that we have the good pad calibration mode if 2D histos!
4812 if (numberofbinsexpected != nbins) {
4813 AliInfo("It doesn't correspond to the mode of pad group calibration!");
4818 // Security for fDebug 3 and 4
4819 if ((fDebug >= 3) &&
4823 AliInfo("This detector doesn't exit!");
4827 // Determine fDet1 and fDet2
4831 fDect1[i] = fFitVoir;
4832 fDect2[i] = fDect1[i] +1;
4836 fDect2[i] = numberofbinsexpected;
4839 CalculXBins(AliTRDgeometry::GetDetector(fDet[0],fDet[1],fDet[2]),i);
4840 fDect1[i] = fXbins[i];
4841 CalculXBins((AliTRDgeometry::GetDetector(fDet[0],fDet[1],fDet[2])+1),i);
4842 fDect2[i] = fXbins[i];
4845 // Create the histos for debugging
4850 // Init the VectorFitCH
4851 fVectorFitCH = new TObjArray();
4852 fCoefCH = new Float_t[2304];
4853 for (Int_t k = 0; k < 2304; k++) {
4856 fScaleFitFactor = 0.0;
4858 // Number of Xbins(detectors or groups of pads) if Vector2d
4859 // Quick verification that we are not out of range!
4860 if (fVectorCH && fPlaCH) {
4862 (fVectorCH->GetEntriesFast() > 0) &&
4863 ((Int_t) fPlaCH->GetEntriesFast() > 0)) {
4864 if ((Int_t) fVectorCH->GetEntriesFast() > numberofbinsexpected) {
4865 AliInfo("ch doesn't correspond to the mode of pad group calibration!");
4868 if ((Int_t) fVectorCH->GetEntriesFast() != (Int_t) fPlaCH->GetEntriesFast()) {
4869 AliInfo("VectorCH doesn't correspond to PlaCH!");
4876 // Debugging: Create the histos
4879 // fDebug == 0 nothing
4884 // Create the histos replique de ch if histos2D
4885 CreateFitHistoCH(nbins,lowedge,upedge);
4888 // Ccreate the histos replique de ch vector2d
4889 CreateFitHistoCH(numberofbinsexpected,0,numberofbinsexpected);
4893 // fDebug == 2 and fFitVoir no histo
4895 if (fFitVoir < numberofbinsexpected) {
4896 AliInfo(Form("We will see the fit of the object %d",fFitVoir));
4899 AliInfo("fFitVoir is out of range of the histo!");
4904 // fDebug == 3 or 4 and fDet
4906 if ((fNz[0] == 0) && (fNrphi[0] == 0)) {
4907 AliInfo("Do you really want to see one detector without pad groups?");
4911 AliInfo(Form("You will see the detector: iPlane %d, iChamb %d, iSect %d"
4912 ,fDet[0],fDet[1],fDet[2]));
4913 // A little geometry:
4914 Int_t rowMax = parCom->GetRowMax(fDet[0],fDet[1],fDet[2]);
4915 Int_t colMax = parCom->GetColMax(fDet[0]);
4916 // Create the histos to visualise
4917 CreateFitHistoCHDB(rowMax,colMax);
4919 CreateFitHistoCH((Int_t) (fDect2[0]-fDect1[0]),fDect1[0],fDect2[0]);
4929 // Number of Xbins (detectors or groups of pads) if vector2d
4930 // Quick verification that we are not out of range!
4931 if (fVectorPH && fPlaPH) {
4933 (fVectorPH->GetEntriesFast() > 0) &&
4934 ((Int_t) fPlaPH->GetEntriesFast() > 0)) {
4935 if ((Int_t) fVectorPH->GetEntriesFast() > numberofbinsexpected) {
4936 AliInfo("ph doesn't correspond to the mode of pad group calibration!");
4939 if ((Int_t) fVectorPH->GetEntriesFast() != (Int_t) fPlaPH->GetEntriesFast()) {
4940 AliInfo("VectorPH doesn't correspond to PlaPH!");
4951 // Debugging: Create the histos
4954 // fDebug == 0 nothing
4959 // Create the histos replique de ch
4960 CreateFitHistoPH(nbins,lowedge,upedge);
4961 CreateFitHistoT0(nbins,lowedge,upedge);
4964 // Create the histos replique de ch if vector2d
4965 CreateFitHistoPH(numberofbinsexpected,0,numberofbinsexpected);
4966 CreateFitHistoT0(numberofbinsexpected,0,numberofbinsexpected);
4970 // fDebug == 2 and fFitVoir no histo
4972 if (fFitVoir < numberofbinsexpected) {
4973 AliInfo(Form("We will see the fit of the object %d",fFitVoir));
4976 AliInfo("fFitVoir is out of range of the histo!");
4981 // fDebug == 3 or 4 and fDet
4983 if ((fNz[1] == 0) &&
4985 AliInfo("Do you really want to see one detector without pad groups?");
4989 AliInfo(Form("You will see the detector: iPlane %d, iChamb %d, iSect %d"
4990 ,fDet[0],fDet[1],fDet[2]));
4991 // A little geometry:
4992 Int_t rowMax = parCom->GetRowMax(fDet[0],fDet[1],fDet[2]);
4993 Int_t colMax = parCom->GetColMax(fDet[0]);
4994 // Create the histos to visualise
4995 CreateFitHistoPHDB(rowMax,colMax);
4996 CreateFitHistoT0DB(rowMax,colMax);
4998 CreateFitHistoPH((Int_t) (fDect2[1]-fDect1[1]),fDect1[1],fDect2[1]);
4999 CreateFitHistoT0((Int_t) (fDect2[1]-fDect1[1]),fDect1[1],fDect2[1]);
5009 // Number of Xbins(detectors or groups of pads) if vector2d
5010 if (fVectorPRF && fPlaPRF){
5012 (fVectorPRF->GetEntriesFast() > 0) &&
5013 (fPlaPRF->GetEntriesFast() > 0)) {
5014 // Quick verification that we are not out of range!
5015 if ((Int_t) fVectorPRF->GetEntriesFast() > numberofbinsexpected) {
5016 AliInfo("ch doesn't correspond to the mode of pad group calibration!");
5019 if ((Int_t) fVectorPRF->GetEntriesFast() != (Int_t) fPlaPRF->GetEntriesFast()) {
5020 AliInfo("VectorPRF doesn't correspond to PlaCH!");
5030 // Debugging: Create the histos
5033 // fDebug == 0 nothing
5038 // Create the histos replique de ch
5039 CreateFitHistoPRF(nbins,lowedge,upedge);
5042 // Create the histos replique de ch
5043 CreateFitHistoPRF(numberofbinsexpected,0,numberofbinsexpected);
5047 // fDebug == 2 and fFitVoir no histo
5049 if (fFitVoir < numberofbinsexpected) {
5050 AliInfo(Form("We will see the fit of the object %d",fFitVoir));
5053 AliInfo("fFitVoir is out of range of the histo!");
5058 // fDebug == 3 or 4 and fDet
5060 if ((fNz[2] == 0) &&
5062 AliInfo("Do you really want to see one detector without pad groups?");
5066 AliInfo(Form("You will see the detector: iPlane %d, iChamb %d, iSect %d"
5067 ,fDet[0],fDet[1],fDet[2]));
5068 // A little geometry:
5069 Int_t rowMax = parCom->GetRowMax(fDet[0],fDet[1],fDet[2]);
5070 Int_t colMax = parCom->GetColMax(fDet[0]);
5071 // Create the histos to visualise
5072 CreateFitHistoPRFDB(rowMax,colMax);
5074 CreateFitHistoPRF((Int_t) (fDect2[2]-fDect1[2]),fDect1[2],fDect2[2]);
5085 //____________Functions for initialising the AliTRDCalibra in the code_________
5086 void AliTRDCalibra::InitfCountDetAndfCount(Int_t i)
5089 // Init the current detector where we are fCountDet and the
5090 // next fCount for the functions Fit...
5093 // Loop on the Xbins of ch!!
5094 fCountDet[i] = -1; // Current detector
5095 fCount[i] = 0; // To find the next detector
5100 // Set countdet to the detector
5101 fCountDet[i] = AliTRDgeometry::GetDetector(fDet[0],fDet[1],fDet[2]);
5103 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
5104 ModePadCalibration(fDet[1],i);
5105 ModePadFragmentation(fDet[0],fDet[1],fDet[2],i);
5107 // Set counter to write at the end of the detector
5108 fCount[i] = fDect1[i] + fNfragZ[i]*fNfragRphi[i];
5114 //____________Functions for initialising the AliTRDCalibra in the code_________
5115 void AliTRDCalibra::UpdatefCountDetAndfCount(Int_t idect, Int_t i)
5118 // See if we are in a new detector and update the
5119 // variables fNfragZ and fNfragRphi if yes
5122 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
5123 // If fDebug == 1 or 0
5124 if ((fDebug == 0) ||
5127 if (fCount[i] == idect) {
5129 // On en est au detector
5132 // Determination of fNnZ, fNnRphi, fNfragZ and fNfragRphi
5133 ModePadCalibration((Int_t) GetChamber(fCountDet[i]),i);
5134 ModePadFragmentation((Int_t) GetPlane(fCountDet[i])
5135 ,(Int_t) GetChamber(fCountDet[i])
5136 ,(Int_t) GetSector(fCountDet[i]),i);
5138 // Set for the next detector
5139 fCount[i] += fNfragZ[i]*fNfragRphi[i];
5147 //____________Functions for initialising the AliTRDCalibra in the code_________
5148 void AliTRDCalibra::ReconstructFitRowMinRowMax(Int_t idect, Int_t i)
5151 // Reconstruct the min pad row, max pad row, min pad col and
5152 // max pad col of the calibration group for the Fit functions
5156 ReconstructionRowPadGroup((Int_t) (idect-(fCount[i]-(fNfragZ[i]*fNfragRphi[i]))),i);
5159 ReconstructionRowPadGroup((Int_t) (idect-fDect1[i]),i);
5164 //____________Functions for initialising the AliTRDCalibra in the code_________
5165 Bool_t AliTRDCalibra::NotEnoughStatistic(Int_t idect, Int_t i)
5168 // For the case where there are not enough entries in the histograms
5169 // of the calibration group, the value present in the choosen database
5170 // will be put. A negativ sign enables to know that a fit was not possible.
5173 // Get the parameter object
5174 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
5176 AliInfo("Could not get CommonParam Manager");
5181 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
5183 AliInfo("Could not get calibDB");
5188 AliInfo(Form("The element %d in this detector %d has not enough statistic to be fitted"
5189 ,idect-(fCount[i]-(fNfragZ[i]*fNfragRphi[i])),fCountDet[i]));
5192 AliInfo("The element has not enough statistic to be fitted");
5195 if ((i == 0) && (fDebug != 2)) {
5197 // Calcul the coef from the database choosen
5198 CalculChargeCoefMean(fCountDet[0],(Int_t) (idect-fDect1[0]),kFALSE);
5200 // Fill the coefCH[2304] with negative value to say: not fitted
5201 for (Int_t k = fRowMin[0]; k < fRowMax[0]; k++) {
5202 for (Int_t j = fColMin[0]; j < fColMax[0]; j++) {
5203 if (GetChamber(fCountDet[0]) == 2) {
5204 fCoefCH[(Int_t)(j*12+k)] = -TMath::Abs(fChargeCoef[3]);
5206 if (GetChamber(fCountDet[0]) != 2) {
5207 fCoefCH[(Int_t)(j*16+k)] = -TMath::Abs(fChargeCoef[3]);
5212 // End of one detector
5213 if ((idect == (fCount[0]-1))) {
5214 FillVectorFitCH((Int_t) fCountDet[0]);
5216 for (Int_t k = 0; k < 2304; k++) {
5223 if ((i == 1) && (fDebug != 2)) {
5225 CalculVdriftCoefMean(fCountDet[1],(Int_t) (idect-fDect1[1]));
5226 CalculT0CoefMean(fCountDet[1],(Int_t) (idect-fDect1[1]));
5228 // Put the default value
5229 if ((fDebug == 1) ||
5233 fCoefVdrift[0]->SetBinContent(idect-fDect1[1]+1,-fVdriftCoef[2]);
5234 fCoefT0[0]->SetBinContent(idect-fDect1[1]+1,-fT0Coef[2]);
5237 fCoefVdrift[1]->SetBinContent(idect-fDect1[1]+1,-fVdriftCoef[2]);
5238 fCoefT0[1]->SetBinContent(idect-fDect1[1]+1,-fT0Coef[2]);
5242 // Put the default value
5244 fVdriftCoef[0] = fVdriftCoef[2];
5245 fVdriftCoef[1] = fVdriftCoef[2];
5247 fT0Coef[0] = fT0Coef[2];
5248 fT0Coef[1] = fT0Coef[2];
5252 // Fill the tree if end of a detector.
5253 // The pointer to the branch stays with the default value 1.5!!!
5255 // Pointer to the branch
5256 for (Int_t k = fRowMin[1]; k < fRowMax[1]; k++) {
5257 for (Int_t j = fColMin[1]; j < fColMax[1]; j++) {
5258 if (GetChamber(fCountDet[1]) == 2) {
5259 fVdriftPad[(Int_t)(j*12+k)] = -TMath::Abs(fVdriftCoef[2]);
5261 if (GetChamber(fCountDet[1]) != 2) {
5262 fVdriftPad[(Int_t)(j*16+k)] = -TMath::Abs(fVdriftCoef[2]);
5267 // End of one detector
5268 if ((idect == (fCount[1]-1)) && (fDebug != 2)) {
5269 FillTreeVdrift((Int_t) fCountDet[1]);
5273 // Fill the tree if end of a detector.
5274 // The pointer to the branch stays with the default value 1.5!!!
5275 // Pointer to the branch
5276 for (Int_t k = fRowMin[1]; k < fRowMax[1]; k++) {
5277 for (Int_t j = fColMin[1]; j < fColMax[1]; j++) {
5278 if (GetChamber(fCountDet[1]) == 2) {
5279 fT0Pad[(Int_t)(j*12+k)] = -TMath::Abs(fT0Coef[2]);
5281 if (GetChamber(fCountDet[1]) != 2) {
5282 fT0Pad[(Int_t)(j*16+k)] = -TMath::Abs(fT0Coef[2]);
5287 // End of one detector
5288 if ((idect == (fCount[1]-1)) && (fDebug != 2)) {
5289 FillTreeT0((Int_t) fCountDet[1]);
5294 if ((i == 2) && (fDebug != 2)) {
5296 CalculPRFCoefMean(fCountDet[2],(Int_t) (idect-fDect1[2]));
5298 if ((fDebug == 1) ||
5300 fCoefPRF[0]->SetBinContent(idect-fDect1[2]+1,fPRFCoef[1]);
5304 fPRFCoef[0] = fPRFCoef[1];
5308 // Fill the tree if end of a detector.
5309 // The pointer to the branch stays with the default value 1.5!!!
5310 // Pointer to the branch
5311 for (Int_t k = fRowMin[2]; k < fRowMax[2]; k++) {
5312 for (Int_t j = fColMin[2]; j < fColMax[2]; j++) {
5313 if((parCom->GetColMax(GetPlane(fCountDet[2])) != (j+1)) && (j != 0)){
5314 if (GetChamber(fCountDet[2]) == 2) {
5315 fPRFPad[(Int_t)(j*12+k)] = -fPRFCoef[1];
5317 if (GetChamber(fCountDet[2]) != 2) {
5318 fPRFPad[(Int_t)(j*16+k)] = -fPRFCoef[1];
5323 if (GetChamber(fCountDet[2]) == 2) {
5324 fPRFPad[(Int_t)(j*12+k)] = -((Float_t) cal->GetPRFWidth(fCountDet[2],j,k));
5326 if (GetChamber(fCountDet[2]) != 2) {
5327 fPRFPad[(Int_t)(j*16+k)] = -((Float_t) cal->GetPRFWidth(fCountDet[2],j,k));
5331 if (GetChamber(fCountDet[2]) == 2) {
5332 fPRFPad[(Int_t)(j*12+k)] = -((Float_t) GetPRFDefault(GetPlane(fCountDet[2])));
5334 if (GetChamber(fCountDet[2]) != 2) {
5335 fPRFPad[(Int_t)(j*16+k)] = -((Float_t) GetPRFDefault(GetPlane(fCountDet[2])));
5342 // End of one detector
5343 if ((idect == (fCount[2]-1)) && (fDebug != 2)) {
5344 FillTreePRF((Int_t) fCountDet[2]);
5353 //____________Functions for initialising the AliTRDCalibra in the code_________
5354 Bool_t AliTRDCalibra::FillInfosFit(Int_t idect, Int_t i)
5357 // Fill the coefficients found with the fits or other
5358 // methods from the Fit functions
5361 // Get the parameter object
5362 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
5364 AliInfo("Could not get CommonParam Manager");
5369 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
5371 AliInfo("Could not get calibDB");
5375 if ((i == 0) && (fDebug != 2)) {
5376 // Fill the coefCH[2304] with fChargeCoef[0]
5377 // that would be negativ only if the fit failed totally
5378 for (Int_t k = fRowMin[0]; k < fRowMax[0]; k++) {
5379 for (Int_t j = fColMin[0]; j < fColMax[0]; j++) {
5380 if (GetChamber(fCountDet[0]) == 2) {
5381 fCoefCH[(Int_t)(j*12+k)] = fChargeCoef[0];
5383 if (GetChamber(fCountDet[0]) != 2) {
5384 fCoefCH[(Int_t)(j*16+k)] = fChargeCoef[0];
5388 // End of one detector
5389 if ((idect == (fCount[0]-1))) {
5390 FillVectorFitCH((Int_t) fCountDet[0]);
5392 for (Int_t k = 0; k < 2304; k++) {
5398 if ((i == 1) && (fDebug != 2)) {
5401 // Pointer to the branch: fVdriftCoef[1] will ne negativ only if the fit failed totally
5402 for (Int_t k = fRowMin[1]; k < fRowMax[1]; k++) {
5403 for (Int_t j = fColMin[1]; j < fColMax[1]; j++) {
5404 if (GetChamber(fCountDet[1]) == 2) {
5405 fVdriftPad[(Int_t)(j*12+k)]=fVdriftCoef[1];
5407 if (GetChamber(fCountDet[1]) != 2) {
5408 fVdriftPad[(Int_t)(j*16+k)]=fVdriftCoef[1];
5412 // End of one detector
5413 if ((idect == (fCount[1]-1)) && (fDebug != 2)) {
5414 FillTreeVdrift((Int_t) fCountDet[1]);
5418 // Pointer to the branch: fT0Coef[1] will ne negativ only if the fit failed totally
5419 for (Int_t k = fRowMin[1]; k < fRowMax[1]; k++) {
5420 for (Int_t j = fColMin[1]; j < fColMax[1]; j++) {
5421 if (GetChamber(fCountDet[1]) == 2) {
5422 fT0Pad[(Int_t)(j*12+k)]=fT0Coef[1];
5424 if (GetChamber(fCountDet[1]) != 2) {
5425 fT0Pad[(Int_t)(j*16+k)]=fT0Coef[1];
5429 // End of one detector
5430 if ((idect == (fCount[1]-1)) && (fDebug != 2)) {
5431 FillTreeT0((Int_t) fCountDet[1]);
5436 if ((i == 2) && (fDebug != 2)) {
5437 // Pointer to the branch
5438 for (Int_t k = fRowMin[2]; k < fRowMax[2]; k++) {
5439 for (Int_t j = fColMin[2]; j < fColMax[2]; j++) {
5440 if ((parCom->GetColMax(GetPlane(fCountDet[2])) != (j+1)) && (j != 0)) {
5441 if (GetChamber(fCountDet[2]) == 2) {
5442 fPRFPad[(Int_t)(j*12+k)] = fPRFCoef[0];
5444 if (GetChamber(fCountDet[2]) != 2) {
5445 fPRFPad[(Int_t)(j*16+k)] = fPRFCoef[0];
5450 if (GetChamber(fCountDet[2]) == 2) {
5451 fPRFPad[(Int_t)(j*12+k)] = (Float_t) cal->GetPRFWidth(fCountDet[2],j,k);
5453 if (GetChamber(fCountDet[2]) != 2) {
5454 fPRFPad[(Int_t)(j*16+k)] = (Float_t) cal->GetPRFWidth(fCountDet[2],j,k);
5458 if (GetChamber(fCountDet[2]) == 2) {
5459 fPRFPad[(Int_t)(j*12+k)] = (Float_t) GetPRFDefault(GetPlane(fCountDet[2]));
5461 if (GetChamber(fCountDet[2]) != 2) {
5462 fPRFPad[(Int_t)(j*16+k)] = (Float_t) GetPRFDefault(GetPlane(fCountDet[2]));
5468 // End of one detector
5469 if ((idect == (fCount[2]-1)) && (fDebug != 2)) {
5470 FillTreePRF((Int_t) fCountDet[2]);
5478 //____________Functions for initialising the AliTRDCalibra in the code_________
5479 Bool_t AliTRDCalibra::WriteFitInfos(Int_t i)
5482 // In the case the user wants to write a file with a tree of the found
5483 // coefficients for the calibration before putting them in the database
5486 TFile *fout = TFile::Open(fWriteNameCoef,"UPDATE");
5487 // Check if the file could be opened
5488 if (!fout || !fout->IsOpen()) {
5489 AliInfo("No File found!");
5493 if ((i == 0) && (fDebug != 2)) {
5495 if ((fDebug == 1) ||
5500 if ((fDebug == 4) ||
5505 fout->WriteTObject(fGain,fGain->GetName(),(Option_t *) "writedelete");
5508 if ((i == 1) && (fDebug != 2)) {
5511 if ((fDebug == 1) ||
5516 if ((fDebug == 4) ||
5521 fout->WriteTObject(fVdrift,fVdrift->GetName(),(Option_t *) "writedelete");
5524 if ((fDebug == 1) ||
5529 if ((fDebug == 4) ||
5534 fout->WriteTObject(fT0,fT0->GetName(),(Option_t *) "writedelete");
5537 if ((i == 2) && (fDebug != 2)) {
5539 if ((fDebug == 1) ||
5544 if ((fDebug == 4) ||
5549 fout->WriteTObject(fPRF,fPRF->GetName(),(Option_t *) "writedelete");
5559 //____________Fill the Error histos in case of fDebug == 1_____________________
5562 //_____________________________________________________________________________
5563 void AliTRDCalibra::ErrorPRF()
5566 // Fill the error histos for fDebug = 1 and fDebug = 4 from the delta histos
5569 for (Int_t k= 0; k < fDeltaPRF->GetNbinsX(); k++) {
5570 if (fDeltaPRF->GetBinContent(k+1) != 0.0) {
5571 fErrorPRF->Fill(fDeltaPRF->GetBinContent(k+1));
5577 //_____________________________________________________________________________
5578 void AliTRDCalibra::ErrorCH()
5581 // Fill the error histos for fDebug = 1 and fDebug = 4 from the delta histos
5584 for (Int_t k= 0; k < fDeltaCharge[0]->GetNbinsX(); k++) {
5585 if (fDeltaCharge[0]->GetBinContent(k+1) != 0.0) {
5586 fErrorCharge[0]->Fill(fDeltaCharge[0]->GetBinContent(k+1));
5589 if (fMeanChargeOn) {
5590 for (Int_t k= 0; k < fDeltaCharge[1]->GetNbinsX(); k++) {
5591 if (fDeltaCharge[1]->GetBinContent(k+1) != 0.0) {
5592 fErrorCharge[1]->Fill(fDeltaCharge[1]->GetBinContent(k+1));
5596 if (fFitChargeBisOn ) {
5597 for (Int_t k= 0; k < fDeltaCharge[2]->GetNbinsX(); k++) {
5598 if (fDeltaCharge[2]->GetBinContent(k+1) != 0.0) {
5599 fErrorCharge[2]->Fill(fDeltaCharge[2]->GetBinContent(k+1));
5606 //_____________________________________________________________________________
5607 void AliTRDCalibra::ErrorPH()
5610 // Fill the error histos for fDebug = 1 and fDebug = 4 from the delta histos
5614 for (Int_t k = 0; k < fDeltaVdrift[0]->GetNbinsX(); k++) {
5615 if (fDeltaVdrift[0]->GetBinContent(k+1) != 0.0) {
5616 fErrorVdrift[0]->Fill(fDeltaVdrift[0]->GetBinContent(k+1));
5620 for (Int_t k = 0; k < fDeltaVdrift[1]->GetNbinsX(); k++) {
5621 if (fDeltaVdrift[1]->GetBinContent(k+1) != 0.0) {
5622 fErrorVdrift[1]->Fill(fDeltaVdrift[1]->GetBinContent(k+1));
5628 //_____________________________________________________________________________
5629 void AliTRDCalibra::ErrorT0()
5632 // Fill the error histos for fDebug = 1 and fDebug = 4 from the delta histos
5636 for (Int_t k = 0; k < fDeltaT0[0]->GetNbinsX(); k++) {
5637 if (fDeltaT0[0]->GetBinContent(k+1) != 0.0) {
5638 fErrorT0[0]->Fill(fDeltaT0[0]->GetBinContent(k+1));
5642 for (Int_t k = 0; k < fDeltaT0[1]->GetNbinsX(); k++) {
5643 if (fDeltaT0[1]->GetBinContent(k+1) != 0.0) {
5644 fErrorT0[1]->Fill(fDeltaT0[1]->GetBinContent(k+1));
5651 //____________Fill Coef DB in case of visualisation of one detector____________
5654 //_____________________________________________________________________________
5655 void AliTRDCalibra::FillCoefVdriftDB()
5658 // Fill the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5661 for (Int_t row = fRowMin[1]; row < fRowMax[1]; row++) {
5662 for (Int_t col = fColMin[1]; col < fColMax[1]; col++) {
5663 fCoefVdriftDB[1]->SetBinContent(row+1,col+1,TMath::Abs(fVdriftCoef[1]));
5665 fCoefVdriftDB[0]->SetBinContent(row+1,col+1,TMath::Abs(fVdriftCoef[0]));
5672 //_____________________________________________________________________________
5673 void AliTRDCalibra::FillCoefT0DB()
5676 // Fill the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5679 for (Int_t row = fRowMin[1]; row < fRowMax[1]; row++) {
5680 for (Int_t col = fColMin[1]; col < fColMax[1]; col++) {
5681 fCoefT0DB[1]->SetBinContent(row+1,col+1,TMath::Abs(fT0Coef[1]));
5683 fCoefT0DB[0]->SetBinContent(row+1,col+1,TMath::Abs(fT0Coef[0]));
5690 //_____________________________________________________________________________
5691 void AliTRDCalibra::FillCoefChargeDB()
5694 // Fill the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5697 for (Int_t row = fRowMin[0]; row < fRowMax[0]; row++) {
5698 for (Int_t col = fColMin[0]; col < fColMax[0]; col++) {
5699 if (fMeanChargeOn) {
5700 fCoefChargeDB[1]->SetBinContent(row+1,col+1,TMath::Abs(fChargeCoef[1]));
5702 if (fFitChargeBisOn) {
5703 fCoefChargeDB[2]->SetBinContent(row+1,col+1,TMath::Abs(fChargeCoef[2]));
5705 fCoefChargeDB[0]->SetBinContent(row+1,col+1,TMath::Abs(fChargeCoef[0]));
5711 //_____________________________________________________________________________
5712 void AliTRDCalibra::FillCoefPRFDB()
5715 // Fill the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5718 for (Int_t row = fRowMin[2]; row < fRowMax[2]; row++) {
5719 for (Int_t col = fColMin[2]; col < fColMax[2]; col++) {
5720 fCoefPRFDB->SetBinContent(row+1,col+1,fPRFCoef[0]);
5727 //____________Plot histos CoefPRF....__________________________________________
5730 //_____________________________________________________________________________
5731 void AliTRDCalibra::PlotCH()
5734 // Plot the histos for fDebug = 1 and fDebug = 4 for the errors
5737 TCanvas *cch1 = new TCanvas("cch1","",50,50,600,800);
5739 TLegend *legch1 = new TLegend(0.4,0.6,0.89,0.89);
5740 legch1->AddEntry(fCoefCharge[3],"f_{g} simulated","l");
5741 if (fMeanChargeOn) {
5742 legch1->AddEntry(fCoefCharge[1],"f_{g} mean","p");
5744 legch1->AddEntry(fCoefCharge[0],"f_{g} fit","p");
5745 if (fFitChargeBisOn ) {
5746 legch1->AddEntry(fCoefCharge[2],"f_{g} fitbis","p");
5749 fCoefCharge[0]->Draw("E2");
5750 if (fMeanChargeOn) {
5751 fCoefCharge[1]->Draw("E2 same");
5753 if (fFitChargeBisOn ) {
5754 fCoefCharge[2]->Draw("E2 same");
5756 fCoefCharge[3]->Draw("same");
5757 legch1->Draw("same");
5759 TCanvas *cch2 = new TCanvas("cch2","",50,50,600,800);
5762 TLegend *legch2 = new TLegend(0.4,0.6,0.89,0.89);
5763 if (fMeanChargeOn) {
5764 legch2->AddEntry(fErrorCharge[1],"f_{g} mean","l");
5766 legch2->AddEntry(fErrorCharge[0],"f_{g} fit","l");
5767 if (fFitChargeBisOn) {
5768 legch2->AddEntry(fErrorCharge[2],"f_{g} fitbis","l");
5770 fErrorCharge[0]->Draw();
5771 if (fMeanChargeOn) {
5772 fErrorCharge[1]->Draw("same");
5774 if (fFitChargeBisOn) {
5775 fErrorCharge[2]->Draw("same");
5777 legch2->Draw("same");
5779 TLegend *legch3 = new TLegend(0.4,0.6,0.89,0.89);
5780 if (fMeanChargeOn) {
5781 legch3->AddEntry(fDeltaCharge[1],"mean","p");
5783 legch3->AddEntry(fDeltaCharge[0],"fit","p");
5784 if (fFitChargeBisOn) {
5785 legch3->AddEntry(fDeltaCharge[2],"fit","p");
5787 fDeltaCharge[0]->Draw("E2");
5788 if (fMeanChargeOn) {
5789 fDeltaCharge[1]->Draw("E2 same");
5791 if (fFitChargeBisOn) {
5792 fDeltaCharge[2]->Draw("E2 same");
5794 legch3->Draw("same");
5798 //_____________________________________________________________________________
5799 void AliTRDCalibra::PlotPH()
5802 // Plot the histos for fDebug = 1 and fDebug = 4 for the errors
5805 TCanvas *cph1 = new TCanvas("cph1","",50,50,600,800);
5807 TLegend *legph1 = new TLegend(0.4,0.6,0.89,0.89);
5808 legph1->AddEntry(fCoefVdrift[2],"v_{real} simulated","l");
5809 legph1->AddEntry(fCoefVdrift[1],"v_{sm} slope 1 method","p");
5812 legph1->AddEntry(fCoefVdrift[0],"v_{fit} fit","p");
5814 fCoefVdrift[1]->Draw("E2");
5815 fCoefVdrift[2]->Draw("same");
5817 fCoefVdrift[0]->Draw("E2 same");
5819 legph1->Draw("same");
5821 TCanvas *cph2 = new TCanvas("cph2","",50,50,600,800);
5824 TLegend *legph2 = new TLegend(0.4,0.6,0.89,0.89);
5825 legph2->AddEntry(fErrorVdrift[1],"v_{sm} slope 1 method","l");
5827 legph2->AddEntry(fErrorVdrift[0],"v_{fit} fit","l");
5829 fErrorVdrift[1]->Draw();
5831 fErrorVdrift[0]->Draw("l,same");
5833 legph2->Draw("same");
5835 TLegend *legph3 = new TLegend(0.4,0.6,0.89,0.89);
5836 legph3->AddEntry(fDeltaVdrift[1],"v_{sm} slope 1 method","p");
5838 legph3->AddEntry(fDeltaVdrift[0],"v_{fit} fit","p");
5840 fDeltaVdrift[1]->Draw("E2");
5842 fDeltaVdrift[0]->Draw("E2 same");
5844 legph3->Draw("same");
5848 //_____________________________________________________________________________
5849 void AliTRDCalibra::PlotT0()
5852 // Plot the histos for fDebug = 1 and fDebug = 4 for the errors
5855 TCanvas *ct01 = new TCanvas("ct01","",50,50,600,800);
5857 TLegend *legt01 = new TLegend(0.4,0.6,0.89,0.89);
5858 legt01->AddEntry(fCoefT0[2],"t0 simulated","l");
5859 legt01->AddEntry(fCoefT0[1],"t0 slope 1 method","p");
5862 legt01->AddEntry(fCoefT0[0],"t0 fit","p");
5864 fCoefT0[1]->Draw("E2");
5865 fCoefT0[2]->Draw("same");
5867 fCoefT0[0]->Draw("E2 same");
5869 legt01->Draw("same");
5871 TCanvas *ct02 = new TCanvas("ct02","",50,50,600,800);
5874 TLegend *legt02 = new TLegend(0.4,0.6,0.89,0.89);
5875 legt02->AddEntry(fErrorT0[1],"t0 slope 1 method","l");
5877 legt02->AddEntry(fErrorT0[0],"t0 fit","l");
5879 fErrorT0[1]->Draw();
5881 fErrorT0[0]->Draw("l,same");
5883 legt02->Draw("same");
5885 TLegend *legt03 = new TLegend(0.4,0.6,0.89,0.89);
5886 legt03->AddEntry(fDeltaT0[1],"t0 slope 1 method","p");
5888 legt03->AddEntry(fDeltaT0[0],"t0 fit","p");
5890 fDeltaT0[1]->Draw("E2");
5892 fDeltaT0[0]->Draw("E2 same");
5894 legt03->Draw("same");
5898 //_____________________________________________________________________________
5899 void AliTRDCalibra::PlotPRF()
5902 // Plot the histos for fDebug = 1 and fDebug = 4 for the errors
5905 TCanvas *cprf1 = new TCanvas("cprf1","",50,50,600,800);
5907 TLegend *legprf1 = new TLegend(0.4,0.6,0.89,0.89);
5908 legprf1->AddEntry(fCoefPRF[1],"#sigma_{real} simulated","l");
5909 legprf1->AddEntry(fCoefPRF[0],"#sigma_{fit} reconstructed","p");
5911 fCoefPRF[0]->Draw("E2");
5912 fCoefPRF[1]->Draw("same");
5913 legprf1->Draw("same");
5915 TCanvas *cprf2 = new TCanvas("cprf2","",50,50,600,800);
5918 TLegend *legprf2 = new TLegend(0.4,0.6,0.89,0.89);
5919 legprf2->AddEntry(fErrorPRF,"#sigma_{fit} reconstructed","l");
5920 fErrorPRF->Draw("");
5921 legprf2->Draw("same");
5923 TLegend *legprf3 = new TLegend(0.4,0.6,0.89,0.89);
5924 legprf3->AddEntry(fDeltaPRF,"#sigma_{fit} reconstructed","p");
5925 fDeltaPRF->Draw("E2");
5926 legprf3->Draw("same");
5931 //____________Plot histos DB___________________________________________________
5934 //_____________________________________________________________________________
5935 void AliTRDCalibra::PlotCHDB()
5938 // Plot the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5941 TCanvas *cchdb = new TCanvas("cchdb","",50,50,600,800);
5942 if ((fFitChargeBisOn) && (fMeanChargeOn)) {
5945 fCoefChargeDB[0]->Draw("LEGO");
5947 fCoefChargeDB[1]->Draw("LEGO");
5949 fCoefChargeDB[2]->Draw("LEGO");
5951 if ((!fFitChargeBisOn) && (fMeanChargeOn)) {
5954 fCoefChargeDB[0]->Draw("LEGO");
5956 fCoefChargeDB[1]->Draw("LEGO");
5960 fCoefChargeDB[0]->Draw("LEGO");
5965 //_____________________________________________________________________________
5966 void AliTRDCalibra::PlotPHDB()
5969 // Plot the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5972 TCanvas *cphdb = new TCanvas("cphdb","",50,50,600,800);
5976 fCoefVdriftDB[0]->Draw("LEGO");
5978 fCoefVdriftDB[1]->Draw("LEGO");
5982 fCoefVdriftDB[1]->Draw("LEGO");
5987 //_____________________________________________________________________________
5988 void AliTRDCalibra::PlotT0DB()
5991 // Plot the histos for fDebug = 3 and fDebug = 4 to visualise the detector
5994 TCanvas *ct0db = new TCanvas("ct0db","",50,50,600,800);
5998 fCoefT0DB[0]->Draw("LEGO");
6000 fCoefT0DB[1]->Draw("LEGO");
6004 fCoefT0DB[1]->Draw("LEGO");
6009 //_____________________________________________________________________________
6010 void AliTRDCalibra::PlotPRFDB()
6013 // Plot the histos for fDebug = 3 and fDebug = 4 to visualise the detector
6016 TCanvas *cprfdb = new TCanvas("cprfdb","",50,50,600,800);
6018 fCoefPRFDB->Draw("LEGO");
6023 //____________Write histos Coef________________________________________________
6026 //_____________________________________________________________________________
6027 void AliTRDCalibra::WriteCH(TFile *fout)
6030 // If wanted, write the debug histos for fDebug = 1 and fDebug = 4
6033 fout->WriteTObject(fCoefCharge[0],fCoefCharge[0]->GetName(),(Option_t *) "OverWrite");
6034 if (fMeanChargeOn) {
6035 fout->WriteTObject(fCoefCharge[1],fCoefCharge[1]->GetName(),(Option_t *) "OverWrite");
6037 if (fFitChargeBisOn) {
6038 fout->WriteTObject(fCoefCharge[2],fCoefCharge[2]->GetName(),(Option_t *) "OverWrite");
6041 fout->WriteTObject(fCoefCharge[3],fCoefCharge[3]->GetName(),(Option_t *) "OverWrite");
6043 fout->WriteTObject(fDeltaCharge[0],fDeltaCharge[0]->GetName(),(Option_t *) "OverWrite");
6044 if (fMeanChargeOn) {
6045 fout->WriteTObject(fDeltaCharge[1],fDeltaCharge[1]->GetName(),(Option_t *) "OverWrite");
6047 if (fFitChargeBisOn) {
6048 fout->WriteTObject(fDeltaCharge[2],fDeltaCharge[2]->GetName(),(Option_t *) "OverWrite");
6051 fout->WriteTObject(fErrorCharge[0],fErrorCharge[0]->GetName(),(Option_t *) "OverWrite");
6052 if (fMeanChargeOn) {
6053 fout->WriteTObject(fErrorCharge[1],fErrorCharge[1]->GetName(),(Option_t *) "OverWrite");
6055 if (fFitChargeBisOn) {
6056 fout->WriteTObject(fErrorCharge[2],fErrorCharge[2]->GetName(),(Option_t *) "OverWrite");
6061 //_____________________________________________________________________________
6062 void AliTRDCalibra::WritePH(TFile *fout)
6065 // If wanted, write the debug histos for fDebug = 1 and fDebug = 4
6069 fout->WriteTObject(fCoefVdrift[0],fCoefVdrift[0]->GetName(),(Option_t *) "OverWrite");
6071 fout->WriteTObject(fCoefVdrift[1],fCoefVdrift[1]->GetName(),(Option_t *) "OverWrite");
6072 fout->WriteTObject(fCoefVdrift[2],fCoefVdrift[2]->GetName(),(Option_t *) "OverWrite");
6075 fout->WriteTObject(fDeltaVdrift[0],fDeltaVdrift[0]->GetName(),(Option_t *) "OverWrite");
6077 fout->WriteTObject(fDeltaVdrift[1],fDeltaVdrift[1]->GetName(),(Option_t *) "OverWrite");
6080 fout->WriteTObject(fErrorVdrift[0],fErrorVdrift[0]->GetName(),(Option_t *) "OverWrite");
6082 fout->WriteTObject(fErrorVdrift[1],fErrorVdrift[1]->GetName(),(Option_t *) "OverWrite");
6086 //_____________________________________________________________________________
6087 void AliTRDCalibra::WriteT0(TFile *fout)
6090 // If wanted, write the debug histos for fDebug = 1 and fDebug = 4
6094 fout->WriteTObject(fCoefT0[0],fCoefT0[0]->GetName(),(Option_t *) "OverWrite");
6096 fout->WriteTObject(fCoefT0[1],fCoefT0[1]->GetName(),(Option_t *) "OverWrite");
6097 fout->WriteTObject(fCoefT0[2],fCoefT0[2]->GetName(),(Option_t *) "OverWrite");
6100 fout->WriteTObject(fDeltaT0[0],fDeltaT0[0]->GetName(),(Option_t *) "OverWrite");
6102 fout->WriteTObject(fDeltaT0[1],fDeltaT0[1]->GetName(),(Option_t *) "OverWrite");
6105 fout->WriteTObject(fErrorT0[0],fErrorT0[0]->GetName(),(Option_t *) "OverWrite");
6107 fout->WriteTObject(fErrorT0[1],fErrorT0[1]->GetName(),(Option_t *) "OverWrite");
6111 //________________________________________________________________________________
6112 void AliTRDCalibra::WritePRF(TFile *fout)
6115 // If wanted, write the debug histos for fDebug = 1 and fDebug = 4
6118 fout->WriteTObject(fCoefPRF[0],fCoefPRF[0]->GetName(),(Option_t *) "OverWrite");
6119 fout->WriteTObject(fCoefPRF[1],fCoefPRF[1]->GetName(),(Option_t *) "OverWrite");
6121 fout->WriteTObject(fDeltaPRF,fDeltaPRF->GetName(), (Option_t *)"OverWrite");
6122 fout->WriteTObject(fErrorPRF,fErrorPRF->GetName(), (Option_t *)"OverWrite");
6127 //____________Write DB Histos__________________________________________________
6130 //_____________________________________________________________________________
6131 void AliTRDCalibra::WriteCHDB(TFile *fout)
6134 // If wanted, write the debug histos for fDebug = 3 and fDebug = 4
6137 fout->WriteTObject(fCoefChargeDB[0],fCoefChargeDB[0]->GetName(),(Option_t *) "OverWrite");
6138 if (fMeanChargeOn) {
6139 fout->WriteTObject(fCoefChargeDB[1],fCoefChargeDB[1]->GetName(),(Option_t *) "OverWrite");
6141 if (fFitChargeBisOn ) {
6142 fout->WriteTObject(fCoefChargeDB[2],fCoefChargeDB[2]->GetName(),(Option_t *) "OverWrite");
6147 //_____________________________________________________________________________
6148 void AliTRDCalibra::WritePHDB(TFile *fout)
6151 // If wanted, write the debug histos for fDebug = 3 and fDebug = 4
6155 fout->WriteTObject(fCoefVdriftDB[0],fCoefVdriftDB[0]->GetName(),(Option_t *) "OverWrite");
6157 fout->WriteTObject(fCoefVdriftDB[1],fCoefVdriftDB[1]->GetName(),(Option_t *) "OverWrite");
6161 //_____________________________________________________________________________
6162 void AliTRDCalibra::WriteT0DB(TFile *fout)
6165 // If wanted, write the debug histos for fDebug = 3 and fDebug = 4
6169 fout->WriteTObject(fCoefT0DB[0],fCoefT0DB[0]->GetName(),(Option_t *) "OverWrite");
6171 fout->WriteTObject(fCoefT0DB[1],fCoefT0DB[1]->GetName(),(Option_t *) "OverWrite");
6175 //_____________________________________________________________________________
6176 void AliTRDCalibra::WritePRFDB(TFile *fout)
6179 // If wanted, write the debug histos for fDebug = 3 and fDebug = 4
6182 fout->WriteTObject(fCoefPRFDB,fCoefPRFDB->GetName(),(Option_t *) "OverWrite");
6187 //____________Calcul Coef Mean_________________________________________________
6190 //_____________________________________________________________________________
6191 Bool_t AliTRDCalibra::CalculT0CoefMean(Int_t dect, Int_t idect)
6194 // For the detector Dect calcul the mean time 0
6195 // for the calibration group idect from the choosen database
6198 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
6200 AliInfo("Could not get calibDB Manager");
6206 if ((fDebug != 2) && fAccCDB) {
6208 for (Int_t row = fRowMin[1]; row < fRowMax[1]; row++) {
6209 for (Int_t col = fColMin[1]; col < fColMax[1]; col++) {
6213 fT0Coef[2] += (Float_t) cal->GetT0(dect,col,row);
6217 fT0Coef[2] += (Float_t) cal->GetT0Average(dect);
6222 fT0Coef[2] = fT0Coef[2] / ((fColMax[1]-fColMin[1])*(fRowMax[1]-fRowMin[1]));
6223 if ((fDebug == 1) ||
6225 fCoefT0[2]->SetBinContent(idect+1,fT0Coef[2]);
6234 //_____________________________________________________________________________
6235 Bool_t AliTRDCalibra::CalculChargeCoefMean(Int_t dect, Int_t idect, Bool_t vrai)
6238 // For the detector Dect calcul the mean gain factor
6239 // for the calibration group idect from the choosen database
6242 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
6244 AliInfo("Could not get calibDB Manager");
6247 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
6249 AliInfo("Could not get CommonParam Manager");
6253 fChargeCoef[3] = 0.0;
6257 for (Int_t row = fRowMin[0]; row < fRowMax[0]; row++) {
6258 for (Int_t col = fColMin[0]; col < fColMax[0]; col++) {
6263 fChargeCoef[3] += (Float_t) cal->GetGainFactor(dect,col,row);
6265 if (vrai && fAccCDB) {
6266 fScaleFitFactor += (Float_t) cal->GetGainFactor(dect,col,row);
6269 fChargeCoef[3] += 1.0;
6271 if (vrai && (!fAccCDB)) {
6272 fScaleFitFactor += 1.0;
6278 fChargeCoef[3] += (Float_t) cal->GetGainFactorAverage(dect);
6280 if (vrai && fAccCDB) {
6281 fScaleFitFactor += ((Float_t) cal->GetGainFactorAverage(dect));
6284 fChargeCoef[3] += 1.0;
6286 if (vrai && (!fAccCDB)) {
6287 fScaleFitFactor += 1.0;
6293 fChargeCoef[3] = fChargeCoef[3] / ((fColMax[0]-fColMin[0])*(fRowMax[0]-fRowMin[0]));
6294 if ((fDebug == 1) ||
6296 fCoefCharge[3]->SetBinContent(idect+1,fChargeCoef[3]);
6305 //_____________________________________________________________________________
6306 Bool_t AliTRDCalibra::CalculPRFCoefMean(Int_t dect, Int_t idect)
6309 // For the detector Dect calcul the mean sigma of pad response
6310 // function for the calibration group idect from the choosen database
6313 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
6315 AliInfo("Could not get calibDB Manager");
6319 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
6321 AliInfo("Could not get CommonParam Manager");
6330 for (Int_t row = fRowMin[2]; row < fRowMax[2]; row++) {
6331 for (Int_t col = fColMin[2]; col < fColMax[2]; col++) {
6332 if ((parCom->GetColMax(GetPlane(dect)) != (col+1)) && (col != 0)) {
6335 fPRFCoef[1] += (Float_t) cal->GetPRFWidth(dect,col,row);
6338 fPRFCoef[1] += GetPRFDefault(GetPlane(dect));
6345 fPRFCoef[1] = fPRFCoef[1]/cot;
6346 if ((fDebug == 1) ||
6348 fCoefPRF[1]->SetBinContent(idect+1,fPRFCoef[1]);
6352 if ((fDebug == 1) ||
6355 fCoefPRF[1]->SetBinContent(idect+1,cal->GetPRFWidth(dect,fColMin[2],fRowMin[2]));
6358 fCoefPRF[1]->SetBinContent(idect+1,GetPRFDefault(GetPlane(dect)));
6369 //_____________________________________________________________________________
6370 Bool_t AliTRDCalibra::CalculVdriftCoefMean(Int_t dect, Int_t idect)
6373 // For the detector dect calcul the mean drift velocity for the
6374 // calibration group idect from the choosen database
6377 AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
6379 AliInfo("Could not get calibDB Manager");
6383 fVdriftCoef[2] = 0.0;
6386 for (Int_t row = fRowMin[1]; row < fRowMax[1]; row++) {
6387 for (Int_t col = fColMin[1]; col < fColMax[1]; col++) {
6392 fVdriftCoef[2] += (Float_t) cal->GetVdrift(dect,col,row);
6395 fVdriftCoef[2] += 1.5;
6401 fVdriftCoef[2] += (Float_t) cal->GetVdriftAverage(dect);
6404 fVdriftCoef[2] += 1.5;
6409 fVdriftCoef[2] = fVdriftCoef[2] / ((fColMax[1]-fColMin[1])*(fRowMax[1]-fRowMin[1]));
6410 if ((fDebug == 1) ||
6412 fCoefVdrift[2]->SetBinContent(idect+1,fVdriftCoef[2]);
6420 //_____________________________________________________________________________
6421 Float_t AliTRDCalibra::GetPRFDefault(Int_t plane) const
6424 // Default width of the PRF if there is no database as reference
6452 //____________Pad group calibration mode_______________________________________
6455 //_____________________________________________________________________________
6456 void AliTRDCalibra::ReconstructionRowPadGroup(Int_t idect, Int_t i)
6459 // For the calibration group idect in a detector calculate the
6460 // first and last row pad and col pad.
6461 // The pads in the interval will have the same calibrated coefficients
6471 if (fNfragZ[i] != 0) {
6472 posc = (Int_t) idect / fNfragZ[i];
6474 if (fNfragRphi[i] != 0) {
6475 posr = (Int_t) idect % fNfragZ[i];
6477 fRowMin[i] = posr * fNnZ[i];
6478 fRowMax[i] = (posr+1) * fNnZ[i];
6479 fColMin[i] = posc * fNnRphi[i];
6480 fColMax[i] = (posc+1) * fNnRphi[i];
6484 //_____________________________________________________________________________
6485 void AliTRDCalibra::CalculXBins(Int_t idect, Int_t i)
6488 // For the detector idect calcul the first Xbins
6493 AliInfo(Form("detector: %d", idect));
6497 Int_t sector = GetSector(idect);
6498 fXbins[i] += sector*(6*fDetChamb2[i]+6*4*fDetChamb0[i]);
6500 // In which chamber?
6501 Int_t chamber = GetChamber(idect);
6503 while (kc < chamber) {
6505 fXbins[i] += 6 * fDetChamb2[i];
6508 fXbins[i] += 6 * fDetChamb0[i];
6514 Int_t plane = GetPlane(idect);
6516 fXbins[i] += plane*fDetChamb2[i];
6519 fXbins[i] += plane*fDetChamb0[i];
6524 //_____________________________________________________________________________
6525 Int_t AliTRDCalibra::SearchInVector(Int_t group, Int_t i) const
6528 // Search if the calibration group "group" has already been
6529 // initialised by a previous track in the vector
6533 for (Int_t k = 0; k < (Int_t) fPlaCH->GetEntriesFast(); k++) {
6534 if (((AliTRDPlace *) fPlaCH->At(k))->GetPlace() == group) {
6542 for (Int_t k = 0; k < (Int_t) fPlaPH->GetEntriesFast(); k++) {
6543 if (((AliTRDPlace *) fPlaPH->At(k))->GetPlace() == group) {
6551 for (Int_t k = 0; k < (Int_t) fPlaPRF->GetEntriesFast(); k++) {
6552 if (((AliTRDPlace *) fPlaPRF->At(k))->GetPlace() == group) {
6563 //_____________________________________________________________________________
6564 Int_t AliTRDCalibra::SearchInTreeVector(TObjArray *vectorplace, Int_t group) const
6567 // Search if the calibration group "group" is present in the tree
6570 for (Int_t k = 0; k < (Int_t) vectorplace->GetEntriesFast(); k++) {
6571 if (((AliTRDPlace *) vectorplace->At(k))->GetPlace() == group) {
6580 //_____________________________________________________________________________
6581 Int_t AliTRDCalibra::SearchBin(Float_t value, Int_t i) const
6589 Int_t fbinmax = (Int_t) value;
6590 Int_t fNumberOfBin = -1;
6596 fNumberOfBin = fNumberBinCharge;
6603 fNumberOfBin = fNumberBinPRF;
6607 if ((value >= fbinmax) ||
6608 (value < fbinmin)) {
6613 reponse = (Int_t) ((fNumberOfBin*(value-fbinmin)) / (fbinmax-fbinmin));
6620 //_____________________________________________________________________________
6621 Bool_t AliTRDCalibra::UpdateVectorCH(Int_t group, Float_t value)
6624 // Fill the vector if a new calibration group "group" or update the
6625 // values of the calibration group "group" if already here
6629 Int_t bin = SearchBin(value,0);
6631 if ((bin < 0) || (bin >= fNumberBinCharge)) {
6636 Int_t place = SearchInVector(group,0);
6640 AliTRDPlace *placegroup = new AliTRDPlace();
6641 placegroup->SetPlace(group);
6642 fPlaCH->Add((TObject *) placegroup);
6644 AliTRDCTInfo *fCHInfo = new AliTRDCTInfo();
6645 UShort_t *entries = new UShort_t[fNumberBinCharge];
6647 for(Int_t k = 0; k < fNumberBinCharge; k++) {
6653 fCHInfo->SetEntries(entries);
6654 // Set in the vector
6655 fVectorCH->Add((TObject *) fCHInfo);
6657 // Group already exits
6660 AliTRDCTInfo *fCHInfo = new AliTRDCTInfo();
6662 fCHInfo = ((AliTRDCTInfo *) fVectorCH->At(place));
6663 UShort_t *entries = fCHInfo->GetEntries();
6667 fCHInfo->SetEntries(entries);
6668 // Update the vector
6669 fVectorCH->AddAt((TObject *) fCHInfo,place);
6676 //_____________________________________________________________________________
6677 Bool_t AliTRDCalibra::UpdateVectorPRF(Int_t group, Float_t x, Float_t y)
6680 // Fill the vector if a new calibration group "group" or update the
6681 // values of the calibration group "group" if already here
6685 Int_t bin = SearchBin(x,2);
6687 if ((bin < 0) || (bin >= fNumberBinPRF)) {
6692 Int_t place = SearchInVector(group,2);
6697 AliTRDPlace *placegroup = new AliTRDPlace();
6698 placegroup->SetPlace(group);
6699 fPlaPRF->Add((TObject *) placegroup);
6700 AliTRDPInfo *fPRFInfo = new AliTRDPInfo();
6702 Float_t *sum = new Float_t[fNumberBinPRF];
6703 Float_t *sumsquare = new Float_t[fNumberBinPRF];
6704 UShort_t *entries = new UShort_t[fNumberBinPRF];
6707 for (Int_t k = 0; k < fNumberBinPRF; k++) {
6715 sumsquare[bin] += y*y;
6719 fPRFInfo->SetSum(sum);
6720 fPRFInfo->SetSumSquare(sumsquare);
6721 fPRFInfo->SetEntries(entries);
6723 // Set in the vector
6724 fVectorPRF->Add((TObject *) fPRFInfo);
6727 // Group already exits
6730 AliTRDPInfo *fPRFInfo = new AliTRDPInfo();
6732 fPRFInfo = (AliTRDPInfo *) fVectorPRF->At(place);
6734 Float_t *sum = fPRFInfo->GetSum();
6735 Float_t *sumsquare = fPRFInfo->GetSumSquare();
6736 UShort_t *entries = fPRFInfo->GetEntries();
6739 Double_t calcul = (((Double_t) fPRFInfo->GetEntries()[bin])
6740 * ((Double_t) fPRFInfo->GetSum()[bin]) + (Double_t) y)
6741 / (((Double_t) fPRFInfo->GetEntries()[bin]) + 1);
6742 sum[bin] = (Float_t) calcul;
6743 Double_t calculsquare = (((Double_t) fPRFInfo->GetSumSquare()[bin])
6744 * ((Double_t) fPRFInfo->GetEntries()[bin]) + ((Double_t) y)*((Double_t) y))
6745 / (((Double_t) fPRFInfo->GetEntries()[bin]) + 1);
6746 sumsquare[bin] = (Float_t) calculsquare;
6750 fPRFInfo->SetSum(sum);
6751 fPRFInfo->SetSumSquare(sumsquare);
6752 fPRFInfo->SetEntries(entries);
6754 // Update the vector
6755 fVectorPRF->AddAt((TObject *) fPRFInfo,place);
6763 //_____________________________________________________________________________
6764 Bool_t AliTRDCalibra::UpdateVectorPH(Int_t group, Int_t time, Float_t value)
6767 // Fill the vector if a new calibration group "group" or update
6768 // the values of the calibration group "group" if already here
6775 (bin >= fTimeMax)) {
6780 Int_t place = SearchInVector(group,1);
6785 AliTRDPlace *placegroup = new AliTRDPlace();
6786 placegroup->SetPlace(group);
6787 fPlaPH->Add((TObject *) placegroup);
6788 AliTRDPInfo *fPHInfo = new AliTRDPInfo();
6790 Float_t *sum = new Float_t[fTimeMax];
6791 Float_t *sumsquare = new Float_t[fTimeMax];
6792 UShort_t *entries = new UShort_t[fTimeMax];
6795 for (Int_t k = 0; k < fTimeMax; k++) {
6803 sumsquare[bin] += value*value;
6807 fPHInfo->SetSum(sum);
6808 fPHInfo->SetSumSquare(sumsquare);
6809 fPHInfo->SetEntries(entries);
6811 // Set in the vector
6812 fVectorPH->Add((TObject *) fPHInfo);
6815 // Group already exits
6818 AliTRDPInfo *fPHInfo = new AliTRDPInfo();
6820 fPHInfo = (AliTRDPInfo *) fVectorPH->At(place);
6822 Float_t *sum = fPHInfo->GetSum();
6823 Float_t *sumsquare = fPHInfo->GetSumSquare();
6824 UShort_t *entries = fPHInfo->GetEntries();
6827 Double_t calcul = (((Double_t) fPHInfo->GetEntries()[bin])
6828 * ((Double_t) fPHInfo->GetSum()[bin]) + (Double_t) value)
6829 / (((Double_t) fPHInfo->GetEntries()[bin]) + 1);
6830 sum[bin] = (Float_t) calcul;
6831 Double_t calculsquare = ((((Double_t) fPHInfo->GetSumSquare()[bin])
6832 * ((Double_t) fPHInfo->GetEntries()[bin]))
6833 + (((Double_t) value) * ((Double_t)value)))
6834 / (((Double_t) fPHInfo->GetEntries()[bin]) + 1);
6835 sumsquare[bin] = (Float_t) calculsquare;
6839 fPHInfo->SetSum(sum);
6840 fPHInfo->SetSumSquare(sumsquare);
6841 fPHInfo->SetEntries(entries);
6843 // Update the vector
6844 fVectorPH->AddAt((TObject *) fPHInfo,place);
6852 //_____________________________________________________________________________
6853 TGraphErrors *AliTRDCalibra::ConvertVectorPHisto(AliTRDPInfo *pInfo
6854 , const Char_t *name) const
6857 // Convert the PInfo in a 1D grapherror, name must contains "PRF"
6858 // if PRF calibration and not "PRF" for Vdrift calibration
6861 TGraphErrors *histo;
6862 const Char_t *pattern1 = "PRF";
6869 Double_t step = 0.0;
6874 if (strstr(name,pattern1)) {
6875 ntimes = fNumberBinPRF;
6880 x = new Double_t[ntimes]; // Xaxis
6881 y = new Double_t[ntimes]; // Mean
6882 ex = new Double_t[ntimes]; // Nentries
6883 ey = new Double_t[ntimes]; // Sum of square/nentries
6886 if (!strstr(name,pattern1)) {
6891 step = (1.0 - (-1.0)) / fNumberBinPRF;
6892 min = -1.0 + step / 2.0;
6896 for (Int_t k = 0; k < ntimes; k++) {
6897 x[k] = min + k*step;
6901 // Fill only if there is more than 0 something
6902 if (pInfo->GetEntries()[k] > 0) {
6903 ex[k] = pInfo->GetEntries()[k];
6904 y[k] = pInfo->GetSum()[k];
6905 ey[k] = pInfo->GetSumSquare()[k];
6909 // Define the TGraphErrors
6910 histo = new TGraphErrors(ntimes,x,y,ex,ey);
6911 histo->SetTitle(name);
6916 //_____________________________________________________________________________
6917 TH1F *AliTRDCalibra::ConvertVectorCTHisto(AliTRDCTInfo *cTInfo
6918 , const Char_t * name) const
6921 // Convert the CTInfo in a 1D histo
6926 Int_t ntimes = fNumberBinCharge;
6927 UShort_t *entries = cTInfo->GetEntries();
6930 histo = new TH1F(name,name,fNumberBinCharge,0,300);
6933 for (Int_t k = 0; k < ntimes; k++) {
6934 histo->SetBinContent(k+1,entries[k]);
6935 histo->SetBinError(k+1,TMath::Sqrt(TMath::Abs(entries[k])));
6942 //_____________________________________________________________________________
6943 TTree *AliTRDCalibra::ConvertVectorCTTreeHisto(TObjArray *vVectorCT
6945 , const Char_t *name
6946 , const Char_t *nametitle) const
6949 // Convert the vector in a tree with two branchs: the group number
6950 // and the TH1F histo reconstructed from the vector
6953 // Size of the things
6954 Int_t ntotal = (Int_t) pPlaCT->GetEntriesFast();
6956 AliInfo("nothing to write!");
6957 TTree *treeCT = new TTree(name,nametitle);
6961 // Variable of the tree
6962 Int_t groupnumber = -1; // Group calibration
6964 TObjArray vectorCT = *vVectorCT;
6965 TObjArray plaCT = *pPlaCT;
6968 TTree *treeCT = new TTree(name,nametitle);
6969 treeCT->Branch("groupnumber",&groupnumber,"groupnumber/I");
6970 treeCT->Branch("histo","TH1F",&histo,32000,0);
6974 while (k < ntotal) {
6976 groupnumber = ((AliTRDPlace *) plaCT.At(0))->GetPlace();
6977 nome += groupnumber;
6978 histo = ConvertVectorCTHisto(((AliTRDCTInfo *) vectorCT.At(0)),nome);
6980 vectorCT.RemoveAt(0);
6981 vectorCT.Compress();
6991 //_____________________________________________________________________________
6992 TTree *AliTRDCalibra::ConvertVectorPTreeHisto(TObjArray *vVectorP
6994 , const Char_t *name
6995 , const Char_t *nametitle) const
6998 // Convert the vector in a tree with two branchs: the group number
6999 // and the TGraphErrors histo reconstructed from the vector.
7000 // The name must contain "PRF" for PRF calibration and not "PRF"
7001 // for Vdrift calibration
7004 // Size of the things
7005 Int_t ntotal = (Int_t) pPlaP->GetEntriesFast();
7007 AliInfo("nothing to write!");
7008 TTree *treeP = new TTree(name,nametitle);
7012 // Variable of the tree
7013 Int_t groupnumber = -1; // Group calibration
7014 TGraphErrors *histo = 0x0;
7015 TObjArray vectorP = *vVectorP;
7016 TObjArray plaP = *pPlaP;
7019 TTree *treeP = new TTree(name,nametitle);
7020 treeP->Branch("groupnumber",&groupnumber,"groupnumber/I");
7021 treeP->Branch("histo","TGraphErrors",&histo,32000,0);
7025 while (k < ntotal) {
7027 groupnumber = ((AliTRDPlace *) plaP.At(0))->GetPlace();
7028 nome += groupnumber;
7029 histo = ConvertVectorPHisto((AliTRDPInfo *) vectorP.At(0),nome);
7031 vectorP.RemoveAt(0);
7042 //_____________________________________________________________________________
7043 TObjArray *AliTRDCalibra::ConvertTreeVector(TTree *tree) const
7046 // Convert the branch groupnumber of the tree taken from
7047 // TRD.calibration.root in case of vector method in a std::vector
7052 TObjArray *vectorplace = new TObjArray();
7054 // Variable of the tree
7055 Int_t groupnumber = -1; // Group calibration
7058 tree->SetBranchAddress("groupnumber",&groupnumber);
7061 Int_t ntotal = tree->GetEntries();
7062 for (Int_t k = 0; k < ntotal; k++) {
7064 AliTRDPlace *placegroupnumber = new AliTRDPlace();
7065 placegroupnumber->SetPlace(groupnumber);
7066 vectorplace->Add((TObject *) placegroupnumber);
7073 //_____________________________________________________________________________
7074 Bool_t AliTRDCalibra::MergeVectorCT(TObjArray *vVectorCT2, TObjArray *pPlaCT2)
7077 // Add the two vectors and place the result in the first
7080 if (((Int_t) pPlaCT2->GetEntriesFast()) != ((Int_t) vVectorCT2->GetEntriesFast())){
7081 AliInfo("VectorCT2 doesn't correspond to PlaCT2!");
7086 for (Int_t k = 0; k < (Int_t) fPlaCH->GetEntriesFast(); k++) {
7088 // Look if PlaCT1[k] it is also in the second vector
7090 for (Int_t j = 0; j < (Int_t) pPlaCT2->GetEntriesFast(); j++) {
7091 if (((AliTRDPlace *) pPlaCT2->At(j))->GetPlace() ==
7092 ((AliTRDPlace *) fPlaCH->At(k))->GetPlace()) {
7098 // If not in the second vector nothing to do
7100 // If in the second vector
7103 AliTRDCTInfo *fCTInfo = new AliTRDCTInfo();
7104 UShort_t *entries = new UShort_t[fNumberBinCharge];
7106 for (Int_t nu = 0; nu < fNumberBinCharge; nu++) {
7107 entries[nu] = ((AliTRDCTInfo *) fVectorCH->At(((AliTRDPlace *) fPlaCH->At(k))->GetPlace()))->GetEntries()[nu]
7108 + ((AliTRDCTInfo *) vVectorCT2->At(((AliTRDPlace *) fPlaCH->At(k))->GetPlace()))->GetEntries()[nu];
7112 fCTInfo->SetEntries(entries);
7114 // Nothing to do on PlaCT1
7116 // Update the vector
7117 fVectorCH->AddAt((TObject *) fCTInfo,((AliTRDPlace *) fPlaCH->At(k))->GetPlace());
7123 // And at the end the vector in CT2 but not in CH1
7124 for (Int_t k = 0; k < (Int_t) pPlaCT2->GetEntriesFast(); k++) {
7126 // Look if pPlaCT2[k] it is also in the second vector
7128 for (Int_t j = 0; j < (Int_t) fPlaCH->GetEntriesFast(); j++) {
7129 if (((AliTRDPlace *) fPlaCH->At(j))->GetPlace() == ((AliTRDPlace *) pPlaCT2->At(k))->GetPlace()) {
7135 // If not in the first vector
7138 AliTRDCTInfo *fCTInfo = new AliTRDCTInfo();
7139 fCTInfo = ((AliTRDCTInfo *) vVectorCT2->At(((AliTRDPlace *) pPlaCT2->At(k))->GetPlace()));
7142 fPlaCH->Add((TObject *) (pPlaCT2->At(k)));
7143 fVectorCH->Add((TObject *) fCTInfo);
7153 //_____________________________________________________________________________
7154 Bool_t AliTRDCalibra::MergeVectorP(TObjArray *vVectorP2
7159 // Add the two vectors and place the result in the first
7162 if (((Int_t) pPlaP2->GetEntriesFast()) != ((Int_t) vVectorP2->GetEntriesFast())) {
7163 AliInfo("VectorP2 doesn't correspond to PlaP2!");
7170 for (Int_t k = 0; k < (Int_t) fPlaPH->GetEntriesFast(); k++) {
7172 // Look if fPlaPH[k] it is also in the second vector
7174 for (Int_t j = 0; j < (Int_t) pPlaP2->GetEntriesFast(); j++) {
7175 if (((AliTRDPlace *) pPlaP2->At(j))->GetPlace() == ((AliTRDPlace *) fPlaPH->At(k))->GetPlace()) {
7181 // If not in the second vector nothing to do
7183 // If in the second vector
7186 AliTRDPInfo *fPInfo = new AliTRDPInfo();
7187 UShort_t *entries = new UShort_t[fTimeMax];
7188 Float_t *sum = new Float_t[fTimeMax];
7189 Float_t *sumsquare = new Float_t[fTimeMax];
7191 for (Int_t nu = 0; nu < fTimeMax; nu++) {
7193 entries[nu] = ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]
7194 + ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu];
7196 Double_t calcul = ((((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSum()[nu])
7197 * ((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]))
7198 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSum()[nu])
7199 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu])))
7200 / ((Double_t) fPInfo->GetEntries()[nu]);
7202 sum[nu] = (Float_t) calcul;
7204 Double_t calculsquare = ((((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSumSquare()[nu])
7205 * ((Double_t) ((AliTRDPInfo *) fVectorPH->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu]))
7206 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetSumSquare()[nu])
7207 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPH->At(k))->GetPlace()))->GetEntries()[nu])))
7208 / ((Double_t) fPInfo->GetEntries()[nu]);
7211 sumsquare[nu] = calculsquare;
7216 fPInfo->SetSum(sum);
7217 fPInfo->SetSumSquare(sumsquare);
7218 fPInfo->SetEntries(entries);
7220 // Nothing to do on PlaCT1
7222 // Update the vector VectorCT1
7223 fVectorPH->AddAt((TObject *) fPInfo,((AliTRDPlace *) fPlaPH->At(k))->GetPlace());
7229 // And at the end the vector in P2 but not in CH1
7230 for (Int_t k = 0; k < (Int_t) pPlaP2->GetEntriesFast(); k++) {
7232 // Look if PlaCT2[k] it is also in the second vector
7234 for (Int_t j = 0; j < (Int_t) fPlaPH->GetEntriesFast(); j++) {
7235 if (((AliTRDPlace *) fPlaPH->At(j))->GetPlace() == ((AliTRDPlace *) pPlaP2->At(k))->GetPlace()) {
7241 // If not in the first vector
7244 AliTRDPInfo *fPInfo = new AliTRDPInfo();
7245 fPInfo = (AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) pPlaP2->At(k))->GetPlace());
7247 // Add at the end of CH1
7248 fPlaPH->Add(((TObject *) pPlaP2->At(k)));
7249 fVectorPH->Add((TObject *) fPInfo);
7261 for (Int_t k = 0; k < (Int_t) fPlaPRF->GetEntriesFast(); k++) {
7263 // Look if fPlaPRF[k] it is also in the second vector
7265 for (Int_t j = 0; j < (Int_t) pPlaP2->GetEntriesFast(); j++) {
7266 if (((AliTRDPlace *) pPlaP2->At(j))->GetPlace() == ((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()) {
7272 // If not in the second vector nothing to do
7274 // If in the second vector
7277 AliTRDPInfo *fPInfo = new AliTRDPInfo();
7278 UShort_t *entries = new UShort_t[fNumberBinPRF];
7279 Float_t *sum = new Float_t[fNumberBinPRF];
7280 Float_t *sumsquare = new Float_t[fNumberBinPRF];
7282 for (Int_t nu = 0; nu < fNumberBinPRF; nu++) {
7284 entries[nu] = ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]
7285 + ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu];
7287 Double_t calcul = ((((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSum()[nu])
7288 * ((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]))
7289 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSum()[nu])
7290 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu])))
7291 / ((Double_t) fPInfo->GetEntries()[nu]);
7293 sum[nu] = (Float_t) calcul;
7295 Double_t calculsquare = ((((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSumSquare()[nu])
7296 * ((Double_t) ((AliTRDPInfo *) fVectorPRF->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu]))
7297 + (((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetSumSquare()[nu])
7298 * ((Double_t) ((AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) fPlaPRF->At(k))->GetPlace()))->GetEntries()[nu])))
7299 / ((Double_t) fPInfo->GetEntries()[nu]);
7301 sumsquare[nu] = calculsquare;
7306 fPInfo->SetSum(sum);
7307 fPInfo->SetSumSquare(sumsquare);
7308 fPInfo->SetEntries(entries);
7310 // Nothing to do on PlaCT1
7312 // Update the vector VectorCT1
7313 fVectorPRF->AddAt((TObject *) fPInfo,((AliTRDPlace *) fPlaPRF->At(k))->GetPlace());
7319 // And at the end the vector in P2 but not in CH1
7320 for (Int_t k = 0; k < (Int_t) pPlaP2->GetEntriesFast(); k++) {
7322 // Look if PlaCT2[k] it is also in the second vector
7324 for (Int_t j = 0; j < (Int_t) fPlaPRF->GetEntriesFast(); j++) {
7325 if (((AliTRDPlace *) fPlaPRF->At(j))->GetPlace() == ((AliTRDPlace *) pPlaP2->At(k))->GetPlace()) {
7331 // If not in the first vector
7334 AliTRDPInfo *fPInfo = new AliTRDPInfo();
7335 fPInfo = (AliTRDPInfo *) vVectorP2->At(((AliTRDPlace *) pPlaP2->At(k))->GetPlace());
7337 // Add at the end of CH1
7338 fPlaPRF->Add(((TObject *) pPlaP2->At(k)));
7339 fVectorPRF->Add((TObject *) fPInfo);
7351 //____________Fit Methods______________________________________________________
7353 //_____________________________________________________________________________
7354 void AliTRDCalibra::FitPente(TH1* projPH, Int_t idect)
7357 // Slope methode for the drift velocity
7361 const Float_t kDrWidth = AliTRDgeometry::DrThick();
7368 fVdriftCoef[1] = 0.0;
7370 TLine *line = new TLine();
7373 TAxis *xpph = projPH->GetXaxis();
7374 Int_t nbins = xpph->GetNbins();
7375 Double_t lowedge = xpph->GetBinLowEdge(1);
7376 Double_t upedge = xpph->GetBinUpEdge(xpph->GetNbins());
7377 Double_t widbins = (upedge - lowedge) / nbins;
7378 Double_t limit = upedge + 0.5 * widbins;
7380 // Beginning of the signal
7381 TH1D *pentea = new TH1D("pentea","pentea",projPH->GetNbinsX(),0,(Float_t) limit);
7382 for (Int_t k = 1; k < projPH->GetNbinsX(); k++) {
7383 pentea->SetBinContent(k,(Double_t) (projPH->GetBinContent(k+1) - projPH->GetBinContent(k)));
7386 binmax = (Int_t) pentea->GetMaximumBin();
7389 AliInfo("Put the binmax from 1 to 2 to enable the fit");
7391 pentea->Fit("pol2","0MR","",TMath::Max(pentea->GetBinCenter(binmax-1),0.0),pentea->GetBinCenter(binmax+1));
7392 Float_t l3P1am = pentea->GetFunction("pol2")->GetParameter(1);
7393 Float_t l3P2am = pentea->GetFunction("pol2")->GetParameter(2);
7395 fPhd[0] = -(l3P1am / (2 * l3P2am));
7398 // Amplification region
7401 for (Int_t kbin = 1; kbin < projPH->GetNbinsX(); kbin ++) {
7402 if (((projPH->GetBinContent(kbin+1) - projPH->GetBinContent(kbin)) <= 0.0) && (ju == 0)) {
7409 AliInfo("Put the binmax from 1 to 2 to enable the fit");
7411 projPH->Fit("pol2","0MR","",TMath::Max(projPH->GetBinCenter(binmax-1),0.0),projPH->GetBinCenter(binmax+1));
7412 Float_t l3P1amf = projPH->GetFunction("pol2")->GetParameter(1);
7413 Float_t l3P2amf = projPH->GetFunction("pol2")->GetParameter(2);
7416 fPhd[1] = -(l3P1amf / (2 * l3P2amf));
7420 TH1D *pente = new TH1D("pente","pente",projPH->GetNbinsX(),0,(Float_t) limit);
7421 for (Int_t k = binmax+4; k < projPH->GetNbinsX(); k++) {
7422 pente->SetBinContent(k,(Double_t) (projPH->GetBinContent(k+1) - projPH->GetBinContent(k)));
7424 binmin = (Int_t) pente->GetMinimumBin();
7427 AliInfo("Put the binmax from 1 to 2 to enable the fit");
7432 ,TMath::Max(pente->GetBinCenter(binmin-1), 0.0)
7433 ,TMath::Min(pente->GetBinCenter(binmin+2),(Double_t) limit));
7434 Float_t l3P1dr = pente->GetFunction("pol2")->GetParameter(1);
7435 Float_t l3P2dr = pente->GetFunction("pol2")->GetParameter(2);
7437 fPhd[2] = -(l3P1dr / (2 * l3P2dr));
7440 if ((fPhd[2] > fPhd[0]) &&
7441 (fPhd[2] > fPhd[1]) &&
7442 (fPhd[1] > fPhd[0])) {
7443 fVdriftCoef[1] = (kDrWidth) / (fPhd[2]-fPhd[1]);
7444 if (fPhd[0] >= 0.0) {
7445 fT0Coef[1] = (fPhd[0] - fT0Shift) / widbins;
7446 if (fT0Coef[1] < 0.0) {
7447 fT0Coef[1] = - TMath::Abs(fT0Coef[2]);
7451 fT0Coef[1] = -TMath::Abs(fT0Coef[2]);
7455 fVdriftCoef[1] = -TMath::Abs(fVdriftCoef[2]);
7456 fT0Coef[1] = -TMath::Abs(fT0Coef[2]);
7459 if ((fDebug == 1) ||
7461 fCoefVdrift[1]->SetBinContent(idect+1,fVdriftCoef[1]);
7462 fCoefT0[1]->SetBinContent(idect+1,fT0Coef[1]);
7463 if (fVdriftCoef[1] > 0.0) {
7464 if (fVdriftCoef[2] != 0.0) {
7465 fDeltaVdrift[1]->SetBinContent(idect+1,(fVdriftCoef[1] - fVdriftCoef[2]) / fVdriftCoef[2]);
7468 if(fT0Coef[1] >= 0.0){
7469 fDeltaT0[1]->SetBinContent(idect+1,(fT0Coef[1] - fT0Coef[2]));
7474 TCanvas *cpentei = new TCanvas("cpentei","cpentei",50,50,600,800);
7477 line->SetLineColor(2);
7478 line->DrawLine(fPhd[0],0,fPhd[0],projPH->GetMaximum());
7479 line->DrawLine(fPhd[1],0,fPhd[1],projPH->GetMaximum());
7480 line->DrawLine(fPhd[2],0,fPhd[2],projPH->GetMaximum());
7481 AliInfo(Form("fPhd[0] (beginning of the signal): %f" ,(Float_t) fPhd[0]));
7482 AliInfo(Form("fPhd[1] (end of the amplification region): %f" ,(Float_t) fPhd[1]));
7483 AliInfo(Form("fPhd[2] (end of the drift region): %f" ,(Float_t) fPhd[2]));
7484 AliInfo(Form("fVriftCoef[1] (with only the drift region(default)): %f",(Float_t) fVdriftCoef[1]));
7496 //_____________________________________________________________________________
7497 void AliTRDCalibra::FitPH(TH1* projPH, Int_t idect)
7500 // Fit methode for the drift velocity
7504 const Float_t kDrWidth = AliTRDgeometry::DrThick();
7507 TAxis *xpph = projPH->GetXaxis();
7508 Double_t upedge = xpph->GetBinUpEdge(xpph->GetNbins());
7510 TF1 *fPH = new TF1("fPH",AliTRDCalibra::PH,-0.05,3.2,6);
7511 fPH->SetParameter(0,0.469); // Scaling
7512 fPH->SetParameter(1,0.18); // Start
7513 fPH->SetParameter(2,0.0857325); // AR
7514 fPH->SetParameter(3,1.89); // DR
7515 fPH->SetParameter(4,0.08); // QA/QD
7516 fPH->SetParameter(5,0.0); // Baseline
7518 TLine *line = new TLine();
7520 fVdriftCoef[0] = 0.0;
7523 if (idect%fFitPHPeriode == 0) {
7525 AliInfo(Form("<AliTRDCalibra::FitPH> The detector %d will be fitted",idect));
7526 fPH->SetParameter(0,(projPH->Integral()-(projPH->GetBinContent(1)*projPH->GetNbinsX())) * 0.00028); // Scaling
7527 fPH->SetParameter(1,fPhd[0] - 0.1); // Start
7528 fPH->SetParameter(2,fPhd[1] - fPhd[0]); // AR
7529 fPH->SetParameter(3,fPhd[2] - fPhd[1]); // DR
7530 fPH->SetParameter(4,0.225); // QA/QD
7531 fPH->SetParameter(5,(Float_t) projPH->GetBinContent(1));
7534 projPH->Fit(fPH,"0M","",0.0,upedge);
7538 TCanvas *cpente = new TCanvas("cpente","cpente",50,50,600,800);
7540 projPH->Fit(fPH,"M+","",0.0,upedge);
7542 line->SetLineColor(4);
7543 line->DrawLine(fPH->GetParameter(1)
7545 ,fPH->GetParameter(1)
7546 ,projPH->GetMaximum());
7547 line->DrawLine(fPH->GetParameter(1)+fPH->GetParameter(2)
7549 ,fPH->GetParameter(1)+fPH->GetParameter(2)
7550 ,projPH->GetMaximum());
7551 line->DrawLine(fPH->GetParameter(1)+fPH->GetParameter(2)+fPH->GetParameter(3)
7553 ,fPH->GetParameter(1)+fPH->GetParameter(2)+fPH->GetParameter(3)
7554 ,projPH->GetMaximum());
7557 if (fPH->GetParameter(3) != 0) {
7558 fVdriftCoef[0] = kDrWidth / (fPH->GetParameter(3));
7559 fT0Coef[0] = fPH->GetParameter(1);
7562 fVdriftCoef[0] = -TMath::Abs(fVdriftCoef[2]);
7563 fT0Coef[0] = -TMath::Abs(fT0Coef[2]);
7566 if ((fDebug == 1) ||
7568 fCoefVdrift[0]->SetBinContent(idect+1,fVdriftCoef[0]);
7569 fCoefT0[0]->SetBinContent(idect+1,fT0Coef[0]);
7570 if (fVdriftCoef[0] > 0.0){
7571 if (fVdriftCoef[2] != 0.0) {
7572 fDeltaVdrift[0]->SetBinContent(idect+1,(fVdriftCoef[0]-fVdriftCoef[2])/fVdriftCoef[2]);
7574 fDeltaT0[0]->SetBinContent(idect+1,(fT0Coef[0]-fT0Coef[2]));
7578 AliInfo(Form("fVdriftCoef[0]: %f",(Float_t) fVdriftCoef[0]));
7585 // Put the default value
7586 if ((fDebug <= 1) ||
7588 fCoefVdrift[0]->SetBinContent(idect+1,fVdriftCoef[2]);
7589 fCoefT0[0]->SetBinContent(idect+1,fT0Coef[2]);
7600 //_____________________________________________________________________________
7601 void AliTRDCalibra::FitPRF(TH1 *projPRF, Int_t idect)
7604 // Fit methode for the sigma of the pad response function
7610 projPRF->Fit("gaus","0M","",-fRangeFitPRF,fRangeFitPRF);
7614 TCanvas *cfit = new TCanvas("cfit","cfit",50,50,600,800);
7616 projPRF->Fit("gaus","M+","",-fRangeFitPRF,fRangeFitPRF);
7620 fPRFCoef[0] = projPRF->GetFunction("gaus")->GetParameter(2);
7621 if(fPRFCoef[0] <= 0.0) fPRFCoef[0] = -fPRFCoef[1];
7623 if ((fDebug == 1) ||
7625 fCoefPRF[0]->SetBinContent(idect+1,fPRFCoef[0]);
7626 if (fPRFCoef[1] != 0.0) {
7627 fDeltaPRF->SetBinContent(idect+1,(fPRFCoef[0]-fPRFCoef[1])/fPRFCoef[1]);
7631 AliInfo(Form("fPRFCoef[0]: %f",(Float_t) fPRFCoef[0]));
7636 //_____________________________________________________________________________
7637 void AliTRDCalibra::FitCH(TH1 *projch, Int_t idect)
7640 // Fit methode for the gain factor
7643 fChargeCoef[0] = 0.0;
7644 fChargeCoef[1] = 0.0;
7645 TF1 *fLandauGaus = new TF1("fLandauGaus",FuncLandauGaus,0,300,5);
7647 fChargeCoef[1] = projch->GetMean();
7648 projch->Fit("landau","0",""
7649 ,(Float_t) fChargeCoef[1]/fBeginFitCharge
7650 ,projch->GetBinCenter(projch->GetNbinsX()));
7651 fL3P0 = projch->GetFunction("landau")->GetParameter(0);
7652 Double_t l3P1 = projch->GetFunction("landau")->GetParameter(1);
7653 fL3P2 = projch->GetFunction("landau")->GetParameter(2);
7655 projch->Fit("gaus","0",""
7656 ,(Float_t) fChargeCoef[1]/fBeginFitCharge
7657 ,projch->GetBinCenter(projch->GetNbinsX()));
7658 Double_t g3P0 = projch->GetFunction("gaus")->GetParameter(0);
7659 fG3P2 = projch->GetFunction("gaus")->GetParameter(2);
7661 fLandauGaus->SetParameters(fL3P0,l3P1,fL3P2,g3P0,fG3P2);
7662 if ((fDebug <= 1) ||
7664 projch->Fit("fLandauGaus","0",""
7665 ,(Float_t) fChargeCoef[1]/fBeginFitCharge
7666 ,projch->GetBinCenter(projch->GetNbinsX()));
7670 TCanvas *cp = new TCanvas("cp","cp",50,50,600,800);
7672 projch->Fit("fLandauGaus","+",""
7673 ,(Float_t) fChargeCoef[1]/fBeginFitCharge
7674 ,projch->GetBinCenter(projch->GetNbinsX()));
7676 fLandauGaus->Draw("same");
7679 if (projch->GetFunction("fLandauGaus")->GetParameter(1) > 0) {
7680 // Calcul of "real" coef
7681 CalculChargeCoefMean(fCountDet[0],(Int_t) idect,kTRUE);
7682 fChargeCoef[0] = projch->GetFunction("fLandauGaus")->GetParameter(1);
7685 // Calcul of "real" coef
7686 CalculChargeCoefMean(fCountDet[0],(Int_t) idect,kFALSE);
7687 fChargeCoef[0] = -TMath::Abs(fChargeCoef[3]);
7691 AliInfo(Form("fChargeCoef[0]: %f",(Float_t) fChargeCoef[0]));
7692 AliInfo(Form("fChargeCoef[1]: %f",(Float_t) fChargeCoef[1]));
7695 if ((fDebug == 1) ||
7697 if (fChargeCoef[0] > 0.0) {
7698 fCoefCharge[0]->SetBinContent(idect+1,fChargeCoef[0]);
7699 fCoefCharge[1]->SetBinContent(idect+1,fChargeCoef[1]);
7700 fDeltaCharge[0]->SetBinContent(idect+1,fChargeCoef[0]);
7701 fDeltaCharge[1]->SetBinContent(idect+1,fChargeCoef[1]);
7704 fL3P0 = fLandauGaus->Integral(0.3*projch->GetMean(),3*projch->GetMean());
7705 fG3P2 = fLandauGaus->GetParameter(2);
7706 fL3P2 = fLandauGaus->GetParameter(4);
7714 //_____________________________________________________________________________
7715 void AliTRDCalibra::FitBisCH(TH1* projch, Int_t idect)
7718 // Fit methode for the gain factor more time consuming
7721 // Setting fit range and start values
7723 //Double_t sv[4] = { l3P2, fChargeCoef[1], projch->Integral("width"), fG3P2 };
7724 Double_t sv[4] = { fL3P2, fChargeCoef[1], fL3P0, fG3P2 };
7725 Double_t pllo[4] = { 0.001, 0.001, 0.001, 0.001 };
7726 Double_t plhi[4] = { 300.0, 300.0, 100000000.0, 300.0 };
7727 Double_t fp[4] = { 1.0, 1.0, 1.0, 1.0 };
7728 Double_t fpe[4] = { 1.0, 1.0, 1.0, 1.0 };
7729 fr[0] = 0.3 * projch->GetMean();
7730 fr[1] = 3.0 * projch->GetMean();
7731 fChargeCoef[2] = 0.0;
7735 TF1 *fitsnr = LanGauFit(projch,&fr[0],&sv[0]
7740 Double_t projchPeak;
7741 Double_t projchFWHM;
7742 LanGauPro(fp,projchPeak,projchFWHM);
7745 fChargeCoef[2] = fp[1];
7748 fChargeCoef[2] = -TMath::Abs(fChargeCoef[3]);
7752 AliInfo(Form("fChargeCoef[2]: %f",(Float_t) fChargeCoef[2]));
7753 TCanvas *cpy = new TCanvas("cpy","cpy",50,50,600,800);
7756 fitsnr->Draw("same");
7759 if ((fDebug == 1) ||
7761 if (fChargeCoef[2] > 0.0) {
7762 fCoefCharge[2]->SetBinContent(idect+1,fChargeCoef[2]);
7763 fDeltaCharge[2]->SetBinContent(idect+1,fChargeCoef[2]);
7773 //_____________________________________________________________________________
7774 void AliTRDCalibra::NormierungCharge()
7777 // Normalisation of the gain factor resulting for the fits
7780 // Calcul of the mean of the fit
7782 Float_t scalefactor = 1.0;
7783 for (Int_t k = 0; k < (Int_t) fVectorFitCH->GetEntriesFast(); k++) {
7785 Int_t detector = ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetDetector();
7786 Float_t *coef = ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetCoef();
7787 if (GetChamber(detector) == 2) {
7790 if (GetChamber(detector) != 2) {
7793 for (Int_t j = 0; j < total; j++) {
7801 fScaleFitFactor = fScaleFitFactor / sum;
7804 fScaleFitFactor = 1.0;
7808 if ((fDebug == 1) ||
7810 if ((fCoefCharge[0]->GetEntries() > 0.0) &&
7811 (fCoefCharge[0]->GetSumOfWeights() > 0.0)) {
7812 scalefactor = fCoefCharge[0]->GetEntries() / fCoefCharge[0]->GetSumOfWeights();
7813 fCoefCharge[0]->Scale(scalefactor);
7814 fDeltaCharge[0]->Scale(scalefactor);
7816 if ((fMeanChargeOn) &&
7817 (fCoefCharge[1]->GetEntries() > 0.0) &&
7818 (fCoefCharge[1]->GetSumOfWeights() > 0.0)) {
7819 fCoefCharge[1]->Scale(fCoefCharge[1]->GetEntries() / fCoefCharge[1]->GetSumOfWeights());
7821 if ((fFitChargeBisOn) &&
7822 (fCoefCharge[2]->GetEntries() > 0.0) &&
7823 (fCoefCharge[2]->GetSumOfWeights() > 0.0)) {
7824 fCoefCharge[2]->Scale(fCoefCharge[2]->GetEntries() / fCoefCharge[2]->GetSumOfWeights());
7826 if ((fMeanChargeOn) &&
7827 (fDeltaCharge[1]->GetEntries() > 0.0) &&
7828 (fDeltaCharge[1]->GetSumOfWeights() > 0.0)) {
7829 fDeltaCharge[1]->Scale(fDeltaCharge[1]->GetEntries() / fDeltaCharge[1]->GetSumOfWeights());
7831 if ((fFitChargeBisOn) &&
7832 (fDeltaCharge[2]->GetEntries() > 0.0) &&
7833 (fDeltaCharge[2]->GetSumOfWeights() > 0.0)) {
7834 fDeltaCharge[2]->Scale(fDeltaCharge[2]->GetEntries() / fDeltaCharge[2]->GetSumOfWeights());
7838 if ((fDebug == 3) ||
7840 fCoefChargeDB[0]->Scale(scalefactor);
7841 if ((fMeanChargeOn) &&
7842 (fCoefChargeDB[1]->GetEntries() > 0.0) &&
7843 (fCoefChargeDB[1]->GetSumOfWeights() > 0.0)) {
7844 fCoefChargeDB[1]->Scale(fCoefChargeDB[1]->GetEntries() / fCoefChargeDB[1]->GetSumOfWeights());
7846 if ((fFitChargeBisOn) &&
7847 (fCoefChargeDB[2]->GetEntries() > 0.0) &&
7848 (fCoefChargeDB[2]->GetSumOfWeights() > 0.0)) {
7849 fCoefChargeDB[2]->Scale(fCoefChargeDB[2]->GetEntries() / fCoefChargeDB[2]->GetSumOfWeights());
7853 if ((fDebug == 1) ||
7856 fDeltaCharge[0]->Add(fCoefCharge[3],-1);
7857 fDeltaCharge[0]->Divide(fCoefCharge[3]);
7859 if (fMeanChargeOn) {
7860 fDeltaCharge[1]->Add(fCoefCharge[3],-1);
7862 if (fMeanChargeOn) {
7863 fDeltaCharge[1]->Divide(fCoefCharge[3]);
7866 if (fFitChargeBisOn) {
7867 fDeltaCharge[2]->Add(fCoefCharge[3],-1);
7869 if (fFitChargeBisOn) {
7870 fDeltaCharge[2]->Divide(fCoefCharge[3]);
7877 //_____________________________________________________________________________
7878 TH1I *AliTRDCalibra::ReBin(TH1I *hist) const
7881 // Rebin of the 1D histo for the gain calibration if needed.
7882 // you have to choose fRebin, divider of fNumberBinCharge
7885 TAxis *xhist = hist->GetXaxis();
7886 TH1I *rehist = new TH1I("projrebin","",(Int_t) xhist->GetNbins()/fRebin
7887 ,xhist->GetBinLowEdge(1)
7888 ,xhist->GetBinUpEdge(xhist->GetNbins()));
7890 AliInfo(Form("fRebin: %d",fRebin));
7892 for (Int_t k = 1; k <= (Int_t) xhist->GetNbins()/fRebin; k++) {
7894 for (Int_t ji = i; ji < i+fRebin; ji++) {
7895 sum += hist->GetBinContent(ji);
7898 rehist->SetBinContent(k,sum);
7903 TCanvas *crebin = new TCanvas("crebin","",50,50,600,800);
7912 //_____________________________________________________________________________
7913 TH1F *AliTRDCalibra::ReBin(TH1F *hist) const
7916 // Rebin of the 1D histo for the gain calibration if needed
7917 // you have to choose fRebin divider of fNumberBinCharge
7920 TAxis *xhist = hist->GetXaxis();
7921 TH1F *rehist = new TH1F("projrebin","",(Int_t) xhist->GetNbins()/fRebin
7922 ,xhist->GetBinLowEdge(1)
7923 ,xhist->GetBinUpEdge(xhist->GetNbins()));
7925 AliInfo(Form("fRebin: %d",fRebin));
7927 for (Int_t k = 1; k <= (Int_t) xhist->GetNbins()/fRebin; k++) {
7929 for (Int_t ji = i; ji < i+fRebin; ji++) {
7930 sum += hist->GetBinContent(ji);
7933 rehist->SetBinContent(k,sum);
7938 TCanvas *crebin = new TCanvas("crebin","",50,50,600,800);
7947 //_____________________________________________________________________________
7948 TH1F *AliTRDCalibra::CorrectTheError(TGraphErrors *hist)
7951 // In the case of the vectors method the trees contains TGraphErrors for PH and PRF
7952 // to be able to add them after
7953 // We convert it to a TH1F to be able to applied the same fit function method
7954 // After having called this function you can not add the statistics anymore
7959 Int_t nbins = hist->GetN();
7960 Double_t *x = hist->GetX();
7961 Double_t *entries = hist->GetEX();
7962 Double_t *mean = hist->GetY();
7963 Double_t *square = hist->GetEY();
7964 fEntriesCurrent = 0;
7970 Double_t step = x[1] - x[0];
7971 Double_t minvalue = x[0] - step/2;
7972 Double_t maxvalue = x[(nbins-1)] + step/2;
7974 rehist = new TH1F("projcorrecterror","",nbins,minvalue,maxvalue);
7976 for (Int_t k = 0; k < nbins; k++) {
7977 rehist->SetBinContent(k+1,mean[k]);
7978 if (entries[k] > 0.0) {
7979 fEntriesCurrent += (Int_t) entries[k];
7980 Double_t d = TMath::Abs(square[k] - (mean[k]*mean[k]));
7981 rehist->SetBinError(k+1,TMath::Sqrt(d/entries[k]));
7984 rehist->SetBinError(k+1,0.0);
7992 //_____________________________________________________________________________
7993 TGraphErrors *AliTRDCalibra::AddProfiles(TGraphErrors *hist1
7994 , TGraphErrors *hist2) const
7997 // In the case of the vectors method we use TGraphErrors for PH and PRF
7998 // to be able to add the them after
7999 // Here we add the TGraphErrors
8002 // First TGraphErrors
8003 Int_t nbins1 = hist1->GetN();
8004 Double_t *x1 = hist1->GetX();
8005 Double_t *ex1 = hist1->GetEX();
8006 Double_t *y1 = hist1->GetY();
8007 Double_t *ey1 = hist1->GetEY();
8009 TGraphErrors *rehist = new TGraphErrors(nbins1);
8011 // Second TGraphErrors
8012 Double_t *ex2 = hist2->GetEX();
8013 Double_t *y2 = hist2->GetY();
8014 Double_t *ey2 = hist2->GetEY();
8016 // Define the Variables for the new TGraphErrors
8022 for (Int_t k = 0; k < nbins1; k++) {
8023 Double_t nentries = 0.0;
8028 if ((ex2[k] == 0.0) &&
8032 if ((ex2[k] == 0.0) &&
8039 if ((ex2[k] > 0.0) &&
8046 if ((ex2[k] > 0.0) &&
8048 nentries = ex1[k] + ex2[k];
8049 y = ( y1[k]*ex1[k]+ y2[k]*ex2[k]) / nentries;
8050 ey = (ey1[k]*ex1[k]+ey2[k]*ex2[k]) / nentries;
8053 rehist->SetPoint(k,x,y);
8054 rehist->SetPointError(k,ex,ey);
8062 //____________Some basic geometry function_____________________________________
8065 //_____________________________________________________________________________
8066 Int_t AliTRDCalibra::GetPlane(Int_t d) const
8069 // Reconstruct the plane number from the detector number
8072 return ((Int_t) (d % 6));
8076 //_____________________________________________________________________________
8077 Int_t AliTRDCalibra::GetChamber(Int_t d) const
8080 // Reconstruct the chamber number from the detector number
8084 return ((Int_t) (d % 30) / fgkNplan);
8088 //_____________________________________________________________________________
8089 Int_t AliTRDCalibra::GetSector(Int_t d) const
8092 // Reconstruct the sector number from the detector number
8096 return ((Int_t) (d / fg));
8101 //____________Fill and Init tree Gain, PRF, Vdrift and T0______________________
8104 //_____________________________________________________________________________
8105 void AliTRDCalibra::InitTreePRF()
8108 // Init the tree where the coefficients from the fit methods can be stored
8112 fPRFPad = new Float_t[2304];
8113 fPRF = new TTree("PRF","PRF");
8114 fPRF->Branch("detector",&fPRFDetector,"detector/I");
8115 fPRF->Branch("width" ,fPRFPad ,"width[2304]/F");
8117 // Set to default value for the plane 0 supposed to be the first one
8118 for (Int_t k = 0; k < 2304; k++) {
8125 //_____________________________________________________________________________
8126 void AliTRDCalibra::FillTreePRF(Int_t countdet)
8129 // Fill the tree with the sigma of the pad response function for the detector countdet
8132 Int_t numberofgroup = 0;
8133 fPRFDetector = countdet;
8136 if (GetChamber((Int_t)(countdet+1)) == 2) {
8137 numberofgroup = 1728;
8140 numberofgroup = 2304;
8143 // Reset to default value for the next
8144 for (Int_t k = 0; k < numberofgroup; k++) {
8145 if (GetPlane((Int_t) (countdet+1)) == 0) {
8148 if (GetPlane((Int_t) (countdet+1)) == 1) {
8151 if (GetPlane((Int_t) (countdet+1)) == 2) {
8154 if (GetPlane((Int_t) (countdet+1)) == 3) {
8157 if (GetPlane((Int_t) (countdet+1)) == 4) {
8160 if (GetPlane((Int_t) (countdet+1)) == 5) {
8169 //_____________________________________________________________________________
8170 void AliTRDCalibra::ConvertVectorFitCHTree()
8173 // Convert the vector stuff to a tree of 1D histos if the user
8174 // want to write it after the fill functions
8177 Int_t detector = -1;
8178 Int_t numberofgroup = 1;
8179 Float_t gainPad[2304];
8181 fGain = new TTree("Gain","Gain");
8182 fGain->Branch("detector",&detector,"detector/I");
8183 fGain->Branch("gainPad" ,gainPad ,"gainPad[2304]/F");
8185 Int_t loop = (Int_t) fVectorFitCH->GetEntriesFast();
8186 for (Int_t k = 0; k < loop; k++) {
8187 detector = ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetDetector();
8188 if (GetChamber((Int_t) ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetDetector()) == 2) {
8189 numberofgroup = 1728;
8192 numberofgroup = 2304;
8194 for (Int_t i = 0; i < numberofgroup; i++) {
8195 if (((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetCoef()[i] >= 0) {
8196 gainPad[i] = ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetCoef()[i] * fScaleFitFactor;
8199 gainPad[i] = (Float_t) ((AliTRDFitCHInfo *) fVectorFitCH->At(k))->GetCoef()[i];
8207 //_____________________________________________________________________________
8208 void AliTRDCalibra::FillTreeVdrift(Int_t countdet)
8211 // Fill the tree with the drift velocities for the detector countdet
8214 Int_t numberofgroup = 0;
8215 fVdriftDetector = countdet;
8217 for (Int_t k = 0; k < 2304; k++) {
8221 if (GetChamber((Int_t)(countdet+1)) == 2) {
8222 numberofgroup = 1728;
8225 numberofgroup = 2304;
8227 // Reset to default value the gain coef
8228 for (Int_t k = 0; k < numberofgroup; k++) {
8229 fVdriftPad[k] = -1.5;
8231 fVdriftDetector = -1;
8235 //_____________________________________________________________________________
8236 void AliTRDCalibra::InitTreePH()
8239 // Init the tree where the coefficients from the fit methods can be stored
8243 fVdriftPad = new Float_t[2304];
8244 fVdrift = new TTree("Vdrift","Vdrift");
8245 fVdrift->Branch("detector",&fVdriftDetector,"detector/I");
8246 fVdrift->Branch("vdrift" ,fVdriftPad ,"vdrift[2304]/F");
8247 // Set to default value for the plane 0 supposed to be the first one
8248 for (Int_t k = 0; k < 2304; k++) {
8249 fVdriftPad[k] = -1.5;
8251 fVdriftDetector = -1;
8255 //_____________________________________________________________________________
8256 void AliTRDCalibra::FillTreeT0(Int_t countdet)
8259 // Fill the tree with the t0 value for the detector countdet
8262 Int_t numberofgroup = 0;
8264 fT0Detector = countdet;
8265 for (Int_t k = 0; k < 2304; k++) {
8269 if (GetChamber((Int_t) (countdet+1)) == 2) {
8270 numberofgroup = 1728;
8273 numberofgroup = 2304;
8275 // Reset to default value the gain coef
8276 for (Int_t k = 0; k < numberofgroup; k++) {
8283 //_____________________________________________________________________________
8284 void AliTRDCalibra::InitTreeT0()
8287 // Init the tree where the coefficients from the fit methods can be stored
8291 fT0Pad = new Float_t[2304];
8292 fT0 = new TTree("T0","T0");
8293 fT0->Branch("detector",&fT0Detector,"detector/I");
8294 fT0->Branch("t0",fT0Pad,"t0[2304]/F");
8295 //Set to default value for the plane 0 supposed to be the first one
8296 for(Int_t k = 0; k < 2304; k++){
8304 //____________Private Functions________________________________________________
8307 //_____________________________________________________________________________
8308 Double_t AliTRDCalibra::PH(Double_t *x, Double_t *par)
8311 // Function for the fit
8314 //TF1 *fAsymmGauss = new TF1("fAsymmGauss",AsymmGauss,0,4,6);
8316 //PARAMETERS FOR FIT PH
8318 //fAsymmGauss->SetParameter(0,0.113755);
8319 //fAsymmGauss->SetParameter(1,0.350706);
8320 //fAsymmGauss->SetParameter(2,0.0604244);
8321 //fAsymmGauss->SetParameter(3,7.65596);
8322 //fAsymmGauss->SetParameter(4,1.00124);
8323 //fAsymmGauss->SetParameter(5,0.870597); // No tail cancelation
8331 Double_t dx = 0.005;
8332 Double_t xs = par[1];
8334 Double_t paras[2] = { 0.0, 0.0 };
8337 if ((xs >= par[1]) &&
8338 (xs < (par[1]+par[2]))) {
8339 //fAsymmGauss->SetParameter(0,par[0]);
8340 //fAsymmGauss->SetParameter(1,xs);
8341 //ss += fAsymmGauss->Eval(xx);
8344 ss += AsymmGauss(&xx,paras);
8346 if ((xs >= (par[1]+par[2])) &&
8347 (xs < (par[1]+par[2]+par[3]))) {
8348 //fAsymmGauss->SetParameter(0,par[0]*par[4]);
8349 //fAsymmGauss->SetParameter(1,xs);
8350 //ss += fAsymmGauss->Eval(xx);
8351 paras[0] = par[0]*par[4];
8353 ss += AsymmGauss(&xx,paras);
8362 //_____________________________________________________________________________
8363 Double_t AliTRDCalibra::AsymmGauss(Double_t *x, Double_t *par)
8366 // Function for the fit
8369 //par[0] = normalization
8377 Double_t par1save = par[1];
8378 //Double_t par2save = par[2];
8379 Double_t par2save = 0.0604244;
8380 //Double_t par3save = par[3];
8381 Double_t par3save = 7.65596;
8382 //Double_t par5save = par[5];
8383 Double_t par5save = 0.870597;
8384 Double_t dx = x[0] - par1save;
8386 Double_t sigma2 = par2save*par2save;
8387 Double_t sqrt2 = TMath::Sqrt(2.0);
8388 Double_t exp1 = par3save * TMath::Exp(-par3save * (dx - 0.5 * par3save * sigma2))
8389 * (1.0 - TMath::Erf((par3save * sigma2 - dx) / (sqrt2 * par2save)));
8390 Double_t exp2 = par5save * TMath::Exp(-par5save * (dx - 0.5 * par5save * sigma2))
8391 * (1.0 - TMath::Erf((par5save * sigma2 - dx) / (sqrt2 * par2save)));
8393 //return par[0]*(exp1+par[4]*exp2);
8394 return par[0] * (exp1 + 1.00124 * exp2);
8398 //_____________________________________________________________________________
8399 Double_t AliTRDCalibra::FuncLandauGaus(Double_t *x, Double_t *par)
8402 // Sum Landau + Gaus with identical mean
8405 Double_t valLandau = par[0] * TMath::Landau(x[0],par[1],par[2]);
8406 //Double_t valGaus = par[3] * TMath::Gaus(x[0],par[4],par[5]);
8407 Double_t valGaus = par[3] * TMath::Gaus(x[0],par[1],par[4]);
8408 Double_t val = valLandau + valGaus;
8414 //_____________________________________________________________________________
8415 Double_t AliTRDCalibra::LanGauFun(Double_t *x, Double_t *par)
8418 // Function for the fit
8421 // par[0]=Width (scale) parameter of Landau density
8422 // par[1]=Most Probable (MP, location) parameter of Landau density
8423 // par[2]=Total area (integral -inf to inf, normalization constant)
8424 // par[3]=Width (sigma) of convoluted Gaussian function
8426 // In the Landau distribution (represented by the CERNLIB approximation),
8427 // the maximum is located at x=-0.22278298 with the location parameter=0.
8428 // This shift is corrected within this function, so that the actual
8429 // maximum is identical to the MP parameter.
8432 // Numeric constants
8433 Double_t invsq2pi = 0.3989422804014; // (2 pi)^(-1/2)
8434 Double_t mpshift = -0.22278298; // Landau maximum location
8436 // Control constants
8437 Double_t np = 100.0; // Number of convolution steps
8438 Double_t sc = 5.0; // Convolution extends to +-sc Gaussian sigmas
8450 // MP shift correction
8451 mpc = par[1] - mpshift * par[0];
8453 // Range of convolution integral
8454 xlow = x[0] - sc * par[3];
8455 xupp = x[0] + sc * par[3];
8457 step = (xupp - xlow) / np;
8459 // Convolution integral of Landau and Gaussian by sum
8460 for (i = 1.0; i <= np/2; i++) {
8462 xx = xlow + (i-.5) * step;
8463 fland = TMath::Landau(xx,mpc,par[0]) / par[0];
8464 sum += fland * TMath::Gaus(x[0],xx,par[3]);
8466 xx = xupp - (i-.5) * step;
8467 fland = TMath::Landau(xx,mpc,par[0]) / par[0];
8468 sum += fland * TMath::Gaus(x[0],xx,par[3]);
8472 return (par[2] * step * sum * invsq2pi / par[3]);
8476 //_____________________________________________________________________________
8477 TF1 *AliTRDCalibra::LanGauFit(TH1 *his, Double_t *fitrange, Double_t *startvalues
8478 , Double_t *parlimitslo, Double_t *parlimitshi
8479 , Double_t *fitparams, Double_t *fiterrors
8480 , Double_t *chiSqr, Int_t *ndf)
8483 // Function for the fit
8487 Char_t funname[100];
8489 AliInfo(Form(funname,"Fitfcn_%s",his->GetName()));
8491 TF1 *ffitold = (TF1 *) gROOT->GetListOfFunctions()->FindObject(funname);
8496 TF1 *ffit = new TF1(funname,LanGauFun,fitrange[0],fitrange[1],4);
8497 ffit->SetParameters(startvalues);
8498 ffit->SetParNames("Width","MP","Area","GSigma");
8500 for (i = 0; i < 4; i++) {
8501 ffit->SetParLimits(i,parlimitslo[i],parlimitshi[i]);
8504 his->Fit(funname,"RB0"); // Fit within specified range, use ParLimits, do not plot
8506 ffit->GetParameters(fitparams); // Obtain fit parameters
8507 for (i = 0; i < 4; i++) {
8508 fiterrors[i] = ffit->GetParError(i); // Obtain fit parameter errors
8510 chiSqr[0] = ffit->GetChisquare(); // Obtain chi^2
8511 ndf[0] = ffit->GetNDF(); // Obtain ndf
8513 return (ffit); // Return fit function
8517 //_____________________________________________________________________________
8518 Int_t AliTRDCalibra::LanGauPro(Double_t *params, Double_t &maxx, Double_t &fwhm)
8521 // Function for the fit
8534 Int_t maxcalls = 10000;
8536 // Search for maximum
8537 p = params[1] - 0.1 * params[0];
8538 step = 0.05 * params[0];
8542 while ((l != lold) && (i < maxcalls)) {
8546 l = LanGauFun(&x,params);
8548 step = -step / 10.0;
8553 if (i == maxcalls) {
8559 // Search for right x location of fy
8560 p = maxx + params[0];
8566 while ( (l != lold) && (i < maxcalls) ) {
8571 l = TMath::Abs(LanGauFun(&x,params) - fy);
8585 // Search for left x location of fy
8587 p = maxx - 0.5 * params[0];
8593 while ((l != lold) && (i < maxcalls)) {
8597 l = TMath::Abs(LanGauFun(&x,params) - fy);
8599 step = -step / 10.0;
8604 if (i == maxcalls) {