AliPhysics  64f4410 (64f4410)
LInfo.cxx
Go to the documentation of this file.
1 #include "LInfo.h"
2 #include <TCanvas.h>
3 #include <TDatime.h>
4 #include <TFile.h>
5 #include <TGrid.h>
6 #include <TH2F.h>
7 #include <TMap.h>
8 #include <TMath.h>
9 #include <TNtuple.h>
10 #include <TProfile2D.h>
11 
13 {
14  if (fIsComputed)
15  return;
16 
17  char id[100];
18  char title[100];
19  const char *sideStr[] = {"A", "C"};
20  const Int_t kNCol = NCol();
21  const Int_t kNRow = NRow();
22  const Int_t kNStrip = NStrip();
23 
24  for (Int_t iSM=0; iSM<kNSM; ++iSM) {
25  Int_t isector = iSM/2;
26  Int_t iside = iSM%2;
27  for (Int_t igain=0; igain<2; ++igain) {
28  if (!fhAmpOverMon[iSM][igain]) {
29  sprintf(id, "hAmpOverMon%02d%d", iSM, igain);
30  sprintf(title, "LED amplitude over LEDMON: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
31  fhAmpOverMon[iSM][igain] = new TH2F(id, title, kNCol, -0.5, kNCol-0.5, kNRow, -0.5, kNRow - 0.5);
32  fhAmpOverMon[iSM][igain]->SetDirectory(0);
33  sprintf(id, "hStripRmsOverMean%02d%d", iSM, igain);
34  sprintf(title, "LedMon rms over mean: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
35  fhStripRmsOverMean[iSM][igain] = new TH1F(id, title, kNCol, -0.5, kNCol-0.5);
36  fhStripRmsOverMean[iSM][igain]->SetDirectory(0);
37  sprintf(id, "hLedRmsOverMean%02d%d", iSM, igain);
38  sprintf(title, "Led rms over mean: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
39  fhLedRmsOverMean[iSM][igain] = new TH2F(id, title, kNCol, -0.5, kNCol-0.5, kNRow, -0.5, kNRow-0.5);
40  fhLedRmsOverMean[iSM][igain]->SetDirectory(0);
41  } else {
42  fhAmpOverMon[iSM][igain]->Reset();
43  fhStripRmsOverMean[iSM][igain]->Reset();
44  fhLedRmsOverMean[iSM][igain]->Reset();
45  }
46  }
47  }
48 
49  for (Int_t iSM=0; iSM<kNSM; ++iSM) {
50  Int_t isector = iSM/2;
51  Int_t iside = iSM%2;
52  for (Int_t igain=0; igain<2; ++igain) {
53  for (Int_t col=0;col<kNCol;++col) {
54  Int_t strip = col / 2;
55  Int_t mbin = fhStrip[iSM][igain]->FindBin(strip);
56  Double_t ledMonAmp = fhStrip[iSM][igain]->GetBinContent(mbin);
57  for (Int_t row=0;row<kNRow;++row) {
58  Int_t lbin = fhLed[iSM][igain]->FindBin(col,row);
59  Double_t ledAmp = fhLed[iSM][igain]->GetBinContent(lbin);
60  if (ledMonAmp!=0) {
61  Double_t weightf = ledAmp/ledMonAmp;
62  fhAmpOverMon[iSM][igain]->SetBinContent(lbin,weightf);
63  }
64  Double_t ledRms = fhLedCount[iSM][igain]->GetBinContent(lbin);
65  if (ledAmp>0) {
66  fhLedRmsOverMean[iSM][igain]->SetBinContent(lbin,ledRms/ledAmp);
67  }
68  }
69  }
70  for (Int_t strip=0;strip<kNStrip;++strip) {
71  Int_t bin = fhStrip[iSM][igain]->FindBin(strip);
72  Double_t mean = fhStrip[iSM][igain]->GetBinContent(bin);
73  Double_t rms = fhStripCount[iSM][igain]->GetBinContent(bin);
74  if (mean>0)
75  fhStripRmsOverMean[iSM][igain]->SetBinContent(bin,rms/mean);
76  }
77  }
78  }
79  fIsComputed = 1;
80 }
81 
83 {
84  // book histograms
85  char id[100];
86  char title[100];
87  const char *sideStr[] = {"A", "C"};
88  const Int_t kNCol = NCol();
89  const Int_t kNRow = NRow();
90  const Int_t kNStrip = NStrip();
91  const char *opt = "S";
92 
93  for (Int_t iSM=0; iSM<kNSM; ++iSM) {
94  Int_t isector = iSM/2;
95  Int_t iside = iSM%2;
96  for (Int_t igain=0; igain<2; ++igain) {
97  sprintf(id, "hStrip%02d%d", iSM, igain);
98  sprintf(title, "LEDMon Amplitude: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
99  fhStrip[iSM][igain] = new TProfile(id, title, kNStrip, -0.5, kNStrip-0.5, opt);
100  fhStrip[iSM][igain]->SetDirectory(0);
101 
102  sprintf(id, "hStripCount%02d%d", iSM, igain);
103  sprintf(title, "LEDMon Entries: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
104  fhStripCount[iSM][igain] = new TProfile(id, title, kNStrip, -0.5, kNStrip-0.5, opt);
105  fhStripCount[iSM][igain]->SetDirectory(0);
106 
107  sprintf(id, "hStripWeighted%02d%d", iSM, igain);
108  sprintf(title, "LEDMon Weighted Amplitude: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
109  fhStripWeighted[iSM][igain] = new TProfile(id, title, kNStrip, -0.5, kNStrip-0.5, opt);
110  fhStripWeighted[iSM][igain]->SetDirectory(0);
111 
112  sprintf(id, "hLed%02d%d", iSM, igain);
113  sprintf(title, "Led Amplitude: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
114  fhLed[iSM][igain] = new TProfile2D(id, title, kNCol, -0.5, kNCol-0.5, kNRow, -0.5, kNRow - 0.5, opt);
115  fhLed[iSM][igain]->SetDirectory(0);
116 
117  sprintf(id, "hLedCount%02d%d", iSM, igain);
118  sprintf(title, "Led Entries: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
119  fhLedCount[iSM][igain] = new TProfile2D(id, title, kNCol, -0.5, kNCol-0.5, kNRow, -0.5, kNRow - 0.5, opt);
120  fhLedCount[iSM][igain]->SetDirectory(0);
121 
122  sprintf(id, "hLedWeighted%02d%d", iSM, igain);
123  sprintf(title, "Led Weighted Amplitude: SM %d (%1s%d) gain %d", iSM, sideStr[iside], isector, igain);
124  fhLedWeighted[iSM][igain] = new TProfile2D(id, title, kNCol, -0.5, kNCol-0.5, kNRow, -0.5, kNRow - 0.5, opt);
125  fhLedWeighted[iSM][igain]->SetDirectory(0);
126 
127  fhAmpOverMon[iSM][igain] = 0;
128  fhLedRmsOverMean[iSM][igain] = 0;
129  fhStripRmsOverMean[iSM][igain] = 0;
130  }
131  }
132 }
133 
134 void LInfo::FillLed(Int_t mod, Int_t gain, Int_t col, Int_t row, Double_t amp, Double_t rms)
135 {
136  if ((amp<0)||(amp>1500))
137  return;
138  fhLed[mod][gain]->Fill(col, row, amp);
139  if (rms>0) {
140  fhLedCount[mod][gain]->Fill(col, row, rms);
141  fhLedWeighted[mod][gain]->Fill(col, row, amp, 1./TMath::Power(rms,2));
142  }
143 }
144 
145 TCanvas *LInfo::DrawHist(Int_t which, Int_t gain, const char *opt) const
146 {
147  TString name;
148  TString hopt;
149  if ((which<1||which>4) && !fIsComputed) {
150  cerr << "Execute Linfo::Compute first" << endl;
151  return 0;
152  }
153  TCanvas *c = new TCanvas("c","c",1600,1600);
154  c->Divide(2,10);
155  c->SetLeftMargin(0.02);
156  c->SetRightMargin(0.02);
157  c->SetTopMargin(0.02);
158  c->SetBottomMargin(0.02);
159  for (Int_t sm=0;sm<kNSM;++sm) {
160  c->cd(sm+1);
161  gPad->SetLeftMargin(0.02);
162  gPad->SetRightMargin(0.02);
163  gPad->SetTopMargin(0.02);
164  gPad->SetBottomMargin(0.02);
165  TH1 *h = 0;
166  switch (which) {
167  case 1:
168  h = GetStripHist(sm,gain);
169  name = Form("StripMeanHist_%d_%d",gain,fRunNo);
170  break;
171  case 2:
172  h = GetStripRmsHist(sm,gain);
173  name = Form("StripRmsHist_%d_%d",gain,fRunNo);
174  break;
175  case 3:
176  h = GetLedHist(sm,gain);
177  name = Form("LedMeanHist_%d_%d",gain,fRunNo);
178  hopt = "colz";
179  break;
180  case 4:
181  h = GetLedRmsHist(sm,gain);
182  name = Form("LedRmsHist_%d_%d",gain,fRunNo);
183  hopt = "colz";
184  break;
185  case 5:
186  h = GetLedMonDispHist(sm,gain);
187  name = Form("LedMonDispHist_%d_%d",gain,fRunNo);
188  break;
189  case 6:
190  h = GetLedDispHist(sm,gain);
191  name = Form("LedDispHist_%d_%d",gain,fRunNo);
192  hopt = "colz";
193  break;
194  default :
195  h = GetLedOverMonHist(sm,gain);
196  name = Form("LedOverMonHist_%d_%d",gain,fRunNo);
197  hopt = "colz";
198  }
199  if (opt)
200  h->Draw(opt);
201  else
202  h->Draw(hopt);
203  }
204  c->SetName(name);
205  c->SetTitle(name);
206  return c;
207 }
208 
209 
210 void LInfo::FillStrip(Int_t mod, Int_t gain, Int_t strip, Double_t amp, Double_t rms)
211 {
212  if ((amp<0)||(amp>1500))
213  return;
214  Int_t spos = strip;
215  if (mod%2==1)
216  spos = 23-strip;
217  fhStrip[mod][gain]->Fill(spos, amp);
218  if (rms>0) {
219  fhStripCount[mod][gain]->Fill(spos, rms);
220  fhStripWeighted[mod][gain]->Fill(spos, amp, 1./TMath::Power(rms,2));
221  }
222 }
223 
225 {
226  const Int_t kGain=gain;
227  Double_t ret=0,all=0;
228  for (Int_t iSM=0; iSM<kNSM; ++iSM) {
229  if (sm>=0&&iSM!=sm)
230  continue;
231  Int_t ncols=NCol();
232  if (iSM>11 && iSM<18)
233  ncols=32;
234  Int_t nrows=NRow();
235  if (iSM==10||iSM==11||iSM==18||iSM==19)
236  nrows=8;
237  for (Int_t col=0;col<ncols;++col) {
238  for (Int_t row=0;row<nrows;++row) {
239  ++all;
240  Int_t bin=fhLed[iSM][kGain]->GetBin(col,row);
241  if (fhLed[iSM][kGain]->GetBinContent(bin)>0.)
242  ++ret;
243  }
244  }
245  }
246  if (all>0)
247  return ret/all;
248  return 0;
249 }
250 
252 {
253  const Int_t kstripGain=gain;
254  Double_t ret=0,all=0;
255  for (Int_t iSM=0; iSM<kNSM; ++iSM) {
256  if (sm>=0&&iSM!=sm)
257  continue;
258  Int_t nstrips=NStrip();
259  if (iSM>11 && iSM<18)
260  nstrips=32/2;
261  for (Int_t strip=1; strip<=nstrips; ++strip) {
262  ++all;
263  if (fhStrip[iSM][kstripGain]->GetBinContent(strip)>0.)
264  ++ret;
265  }
266  }
267  if (all>0)
268  return ret/all;
269  return 0;
270 }
271 
272 void LInfo::Print(Option_t *option) const
273 {
274  cout << "Runno: " << fRunNo << endl;
275 }
TProfile * fhStripWeighted[kNSM][2]
Definition: LInfo.h:48
double Double_t
Definition: External.C:58
Definition: External.C:236
const char * title
Definition: MakeQAPdf.C:27
Double_t FracLeds(Int_t sm, Int_t gain=1) const
Definition: LInfo.cxx:224
TH1 * GetLedMonDispHist(Int_t sm, Int_t gain=1) const
Definition: LInfo.h:27
TCanvas * DrawHist(Int_t which, Int_t gain=1, const char *opt=0) const
Definition: LInfo.cxx:145
TProfile * fhStripCount[kNSM][2]
Definition: LInfo.h:47
TProfile * fhStrip[kNSM][2]
Definition: LInfo.h:46
TCanvas * c
Definition: TestFitELoss.C:172
void Print(Option_t *option="") const
Definition: LInfo.cxx:272
TH2 * fhAmpOverMon[kNSM][2]
Definition: LInfo.h:52
static Int_t NCol()
Definition: LInfo.h:39
TProfile2D * fhLedCount[kNSM][2]
Definition: LInfo.h:50
TProfile2D * fhLedWeighted[kNSM][2]
Definition: LInfo.h:51
int Int_t
Definition: External.C:63
Bool_t fIsComputed
RMS over Mean for Led.
Definition: LInfo.h:55
static Int_t NStrip()
Definition: LInfo.h:41
void FillLed(Int_t mod, Int_t gain, Int_t col, Int_t row, Double_t amp, Double_t rms)
Definition: LInfo.cxx:134
TProfile2D * GetLedRmsHist(Int_t sm, Int_t gain=1) const
Definition: LInfo.h:24
void Compute()
Definition: LInfo.cxx:12
Double_t FracStrips(Int_t sm, Int_t gain=1) const
Definition: LInfo.cxx:251
TProfile * GetStripHist(Int_t sm, Int_t gain=1) const
Definition: LInfo.h:17
TH1 * fhStripRmsOverMean[kNSM][2]
Led/LedMon ratio.
Definition: LInfo.h:53
static const Int_t kNSM
Definition: LInfo.h:37
TProfile2D * GetLedHist(Int_t sm, Int_t gain=1) const
Definition: LInfo.h:23
void CreateHistograms()
Definition: LInfo.cxx:82
TH2 * GetLedDispHist(Int_t sm, Int_t gain=1) const
Definition: LInfo.h:28
Int_t fRunNo
Definition: LInfo.h:45
void FillStrip(Int_t mod, Int_t gain, Int_t strip, Double_t amp, Double_t rms)
Definition: LInfo.cxx:210
TProfile * GetStripRmsHist(Int_t sm, Int_t gain=1) const
Definition: LInfo.h:18
const char Option_t
Definition: External.C:48
TH2 * fhLedRmsOverMean[kNSM][2]
RMS over Mean for LedMon.
Definition: LInfo.h:54
TProfile2D * fhLed[kNSM][2]
Definition: LInfo.h:49
Definition: External.C:196
TH2 * GetLedOverMonHist(Int_t sm, Int_t gain=1) const
Definition: LInfo.h:26
static Int_t NRow()
Definition: LInfo.h:40