Extacting the OCDB in a separate module. The detectors have write permission in the...
[u/mrichter/AliRoot.git] / PHOS / macros / BadModules / makeBadModESD.C
CommitLineData
d76dd43e 1#if !defined(__CINT__) || defined(__MAKECINT__)
2#include "TH3.h"
3#include "TF1.h"
4#include "TFile.h"
5#include "TLine.h"
6#include "TCanvas.h"
7#include "AliPHOSEmcBadChannelsMap.h"
8#include "AliCDBId.h"
9#include "AliCDBMetaData.h"
10#include "AliCDBEntry.h"
11#include "AliCDBManager.h"
12#endif
13void makeBadModESD(Int_t mod=3){
14 //Here we use 3D distributions of pedestals and spectra
15 //produced by the macro scanESD.C
16 //Having these distributions we calculate distributions
17 //over number of clusters in "Soft" and "Hard" regions and
18 //exlude cells, havig considerable deviations from avrage
19 //Result is stored in TH2S in root file and as CDB root file
20 //
21 //Currect version is valid for one PHOS module (e.g. cosmic run)
22 //Author D.Peressounko
23
24 TFile * f = new TFile("scan.root") ;
25 TH3D * hSp = (TH3D*)f->Get("hEsd") ;
26
27 //Boundary of "Soft" and "Hard" regions
28 //Soft region used to exclude dead or almost dead cells
29 //Hard region used to exclude "singers"
30 //All energies in GeV
31 Double_t aSoftMin=0.22 ;
32 Double_t aSoftMax=0.5 ;
33 Double_t aHardMin=1. ;
34 Double_t aHardMax=2. ;
35
36 //Final bad map
37 TH2D * hBadMap = new TH2D("hBadMap","Bad Modules map",64,0.,64.,56,0.,56.) ;
38
39 TH2D * hSoft = new TH2D("hSoft","Number of soft clusters",64,0.,64.,56,0.,56.) ;
40 TH2D * hHard = new TH2D("hHard","Number of hard clusters",64,0.,64.,56,0.,56.) ;
41
42 //The range of these histos depends on statistics, one probaly will have to increase it
43 TH1D * hAllSoft = new TH1D("hallSoft","Summary of Soft",200,0.,200.) ;
44 TH1D * hAllHard = new TH1D("hallHard","Summary of Hard",100,0.,100.) ;
45
46 TAxis * ax = hSp->GetZaxis() ;
47 Int_t iSoftMin=ax->FindBin(aSoftMin) ;
48 Int_t iSoftMax=ax->FindBin(aSoftMax) ;
49 Int_t iHardMin=ax->FindBin(aHardMin) ;
50 Int_t iHardMax=ax->FindBin(aHardMax) ;
51
52 //Fill histograms with number of counts per cell
53 for(Int_t i=1; i<=hSp->GetXaxis()->GetNbins(); i++){
54 for(Int_t j=1; j<=hSp->GetYaxis()->GetNbins(); j++){
55 Double_t soft=0 ;
56 for(Int_t k=iSoftMin; k<=iSoftMax;k++){
57 soft+=hSp->GetBinContent(i,j,k) ;
58 }
59 Double_t hard=0. ;
60 for(Int_t k=iHardMin; k<=iHardMax;k++){
61 hard+=hSp->GetBinContent(i,j,k) ;
62 }
63 hSoft->SetBinContent(i,j,soft) ;
64 hHard->SetBinContent(i,j,hard) ;
65 hAllSoft->Fill(soft) ;
66 hAllHard->Fill(hard) ;
67 }
68 }
69
70
71 //=====================================================================================
72 //
73 //Usually the part below needs some manual intervention to properly select good modules region
74 //
75 //Fit overall distribution with Gauss and calculate limits of "goodness"
76 //usually +-3-4 sigma for soft and ~4-5 sigma for hard (less statistics)
77 TF1 * gaus = new TF1("gaus","[0]*exp(-(x-[1])*(x-[1])/2./[2]/[2])",0.,1000.) ;
78 gaus->SetParameters(100.,10.,20.) ;
79 gaus->FixParameter(1,0.) ;
80 hAllSoft->Fit(gaus,"","",10,80) ;
81 Double_t m = gaus->GetParameter(1) ;
82 Double_t s = TMath::Abs(gaus->GetParameter(2)) ;
83 //one can use either sigma or hardwired parameterization
84 //of limits of goodness
85 Double_t cutSoftMin= 8. ; //TMath::Max(0.,m-4.*s) ;
86 Double_t cutSoftMax= 80. ; //m+4.*s ;
87 Double_t softH= gaus->GetParameter(0) ;
88
89 //repeat same procedure for hard
90 gaus->SetParameters(100.,2.,10.) ;
91 gaus->FixParameter(1,0.) ;
92 hAllHard->Fit(gaus,"","",1.,20.) ;
93 m = gaus->GetParameter(1) ;
94 s = TMath::Abs(gaus->GetParameter(2)) ;
95 Double_t cutHardMin= TMath::Max(0.,m-5.*s) ;
96 Double_t cutHardMax= m+5.*s ;
97 Double_t hardH= gaus->GetParameter(0) ;
98
99 //=====================================================================================
100
101 //Draw results
102 TCanvas * cm = new TCanvas("Soft","Soft",15,22,500,300) ;
103 hSoft->Draw("colz") ;
104
105 TCanvas * cr = new TCanvas("Hard","Hard",520,22,500,300) ;
106 hHard->Draw("colz") ;
107
108 TCanvas * cma = new TCanvas("AllSoft","Soft Multiplicity",15,330,500,300) ;
109 cma->SetLogy() ;
110 hAllSoft->Draw() ;
111 TLine * l1 = new TLine(cutSoftMin,0.1,cutSoftMin,2.*softH) ;
112 l1->SetLineColor(2) ;
113 l1->Draw() ;
114 TLine * l2 = new TLine(cutSoftMax,0.1,cutSoftMax,2.*softH) ;
115 l2->SetLineColor(2) ;
116 l2->Draw() ;
117
118 TCanvas * cra = new TCanvas("AllHard","Hard Multiplicity",520,330,500,300) ;
119 cra->SetLogy() ;
120 hAllHard->Draw() ;
121 printf("Min=%f, max=%f \n",cutHardMin,cutHardMax) ;
122 TLine * l3 = new TLine(cutHardMin,0.1,cutHardMin,2.*hardH) ;
123 l3->SetLineColor(2) ;
124 l3->Draw() ;
125 TLine * l4 = new TLine(cutHardMax,0.1,cutHardMax,2.*hardH) ;
126 l4->SetLineColor(2) ;
127 l4->Draw() ;
128
129
130 //now calculate bad map
131 //Exclude everething beyond the limits
132 AliPHOSEmcBadChannelsMap badMap;
133 for(Int_t i=1; i<=hSp->GetXaxis()->GetNbins(); i++){
134 for(Int_t j=1; j<=hSp->GetYaxis()->GetNbins(); j++){
135 Double_t soft = hSoft->GetBinContent(i,j) ;
136 Double_t hard= hHard->GetBinContent(i,j) ;
137 if(soft>=cutSoftMin && soft <=cutSoftMax &&
138 hard>=cutHardMin && hard <=cutHardMax){
139 hBadMap->SetBinContent(i,j,1.) ;
140 hBadMap->SetBinContent(i,j,1.) ;
141 }
142 else{
143 hBadMap->SetBinContent(i,j,0.) ;
144 badMap.SetBadChannel(mod,i,j); //module, col,row
145 }
146 }
147 }
148
149 TCanvas * cb = new TCanvas("BadMap","Bad map",15,630,500,300) ;
150 hBadMap->DrawClone("colz") ;
151
152 TFile * fout = new TFile("BadMap.root","recreate") ;
153 hBadMap->Write() ;
154 fout->Close() ;
155
156 //put now result to local CDB
157 AliCDBManager *CDB = AliCDBManager::Instance();
162637e4 158 // CDB->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
d76dd43e 159 CDB->SetDefaultStorage("local://./");
160 // CDB->SetSpecificStorage("local://./","PHOS");
161
162 AliCDBMetaData *md= new AliCDBMetaData();
163 md->SetResponsible("Dmitri Peressounko");
164 md->SetComment("Dead channels for run 1234");
165 AliCDBId id("PHOS/Calib/EmcBadChannels",0,999999);
166 CDB->Put(&badMap,id, md);
167
168
169}