BmnRoot
Loading...
Searching...
No Matches
ToF400DigitAnalysis_period5.c
Go to the documentation of this file.
1/*
2 * To change this license header, choose License Headers in Project Properties.
3 * To change this template file, choose Tools | Templates
4 * and open the template in the editor.
5 */
6
7#include <stdio.h>
8#include <fstream>
9#include <iostream>
10
11#include <TMath.h>
12#include "TChain.h"
13#include "TH1.h"
14#include "TH1F.h"
15#include "TH1S.h"
16#include "TH2S.h"
17#include "TH2F.h"
18#include "TCanvas.h"
19#include "TFile.h"
20#include "TList.h"
21#include "TDirectory.h"
22#include "TPad.h"
23//#include "BmnEventHeader.h"
24
25const UInt_t kNST = 48;
26const Double_t kTimeBin = 24. / 1024.;
29Double_t CorrPlane7Coeff_It1[5][4] = {
30 {12.24, -3.919, 0.2624, -0.005153},
31 {-122.9, 14.84, -0.6158, 0.008724},
32 {-12.34, 0.482, 0., 0.},
33 {-12.28, 0.479, 0., 0.},
34 {-13.14, 0.5034, 0., 0.}
35};
36Double_t CorrT0Coeff_It1[4][4] = {
37 {-5.486, 1.087, -0.07936, 0.00208},
38 {10.77, -2.635, 0.1984, -0.004631},
39 {-2.163, 0.1362, 0., 0.},
40 {200, 0, 0, 0}
41};
42Int_t HitPlane2T0 = 0, HitPlane7T0 = 0, EffPlane2 = 0, EffPlane7 = 0;
43Bool_t FlagDigit[10][kNST];
44
45struct ToF400Detector {
46 Double_t sTimeL[kNST];
47 Double_t sTimeR[kNST];
48 Double_t sTime[kNST];
49 Double_t sAmpL[kNST];
50 Double_t sAmpR[kNST];
51 Double_t sAmp[kNST];
52 Bool_t sStrHit[kNST];
53 Int_t sNHits;
54};
55
56struct T0Detector {
57 Double_t sTime;
58 Double_t sAmp;
59 Int_t sNHits;
60};
61
62void
63ToF400DigitAnalysis_period5(TString file = "", Int_t nEvForRead = 0, Int_t Periud = 5) {
64 gROOT->LoadMacro("$VMCWORKDIR/macro/run/bmnloadlibs.C");
65 bmnloadlibs(); // load BmnRoot libraries
66
67 /*if (Run < 796) {
68 CatTimeBC2Low = 250.;
69 CatTimeBC2High = 450.;
70 } else {
71 CatTimeBC2Low = 390.;
72 CatTimeBC2High = 590.;
73 }//*/
74
75 TChain *eveTree = new TChain("cbmsim");
76 TString inName;
77 inName = file;
78 cout << "Open file " << inName << endl << endl;
79 eveTree->Add(inName);
80
81 TClonesArray *ToF400Digits;
82 eveTree->SetBranchAddress("TOF400", &ToF400Digits);
83
84 TClonesArray *BC2Digits;
85 eveTree->SetBranchAddress("T0", &BC2Digits);
86
87 TClonesArray *EventHeader;
88 eveTree->SetBranchAddress("EventHeader", &EventHeader);
89
90 Long64_t nEvents = eveTree->GetEntries();
91 if (nEvForRead == 0) nEvForRead = nEvents;
92 cout << "Will be readed " << nEvForRead << " events from " << nEvents << endl;
93
94 //------variables----------------------------------------------------------------------
95 char hid[50];
96 Int_t Point;
97 ofstream f_call, f_dt;
98 TF1 *f_gaus = new TF1("f_gaus", "gaus", -20, 40);
99
100 ToF400Detector Plane2;
101 ToF400Detector Plane7;
102 T0Detector T0;
103 CleanTof(&Plane2);
104 CleanTof(&Plane7);
105 CleanBC2(&T0);
106
107 TList * ListStat = new TList();
108 TList * ListAmpTofPlane2 = new TList();
109 TList * ListTimeTofPlane2 = new TList();
110 TList * ListLRTofPlane2 = new TList();
111 TList * ListAmpTofPlane7 = new TList();
112 TList * ListTimeTofPlane7 = new TList();
113 TList * ListLRTofPlane7 = new TList();
114 TList * ListLRCorrTofPlane7 = new TList();
115 TList * ListLRvsAmpLTofPlane7 = new TList();
116 TList * ListLRvsAmpRTofPlane7 = new TList();
117 TList * ListAnalysis = new TList();
118 TString name, name2;
119
120 TH1I * hAmpTofPlane2[kNST + 1];
121 TH1I * hTimeTofPlane2[kNST + 1];
122 TH1I * hLevtMinusRightPlane2[kNST + 1];
123 for (Int_t i = 0; i < kNST + 1; i++) {
124 name = Form("hAmpTofPlane2_%d", i);
125 hAmpTofPlane2[i] = new TH1I(name, name, 1024, 0, 50.);
126 ListAmpTofPlane2->Add(hAmpTofPlane2[i]);
127 name = Form("hTimeTofPlane2_%d", i);
128 hTimeTofPlane2[i] = new TH1I(name, name, 1000, 0, 1000.);
129 ListTimeTofPlane2->Add(hTimeTofPlane2[i]);
130 name = Form("hLevtMinusRightPlane2_%d", i);
131 hLevtMinusRightPlane2[i] = new TH1I(name, name, 1024, -12., 12.);
132 ListLRTofPlane2->Add(hLevtMinusRightPlane2[i]);
133 }
134 TH2I * hAmpLRPlane2 = new TH2I("hAmpLRPlane2", "hAmpLRPlane2", 1024, 0., 50., 1024, 0., 50.);
135 ListAmpTofPlane2->Add(hAmpLRPlane2);
136
137 TH1I * hAmpTofPlane7[kNST + 1];
138 TH1I * hTimeTofPlane7[kNST + 1];
139 TH1I * hLeftMinusRightPlane7[kNST + 1];
140 TH2S * hLeftRightVsAmpLPlane7[kNST + 1];
141 TH2S * hLeftRightVsAmpRPlane7[kNST + 1];
142 TH2S * hLeftRightCorelationRPlane7[kNST + 1];
143 for (Int_t i = 0; i < kNST + 1; i++) {
144 name = Form("hAmpTofPlane7_%d", i);
145 hAmpTofPlane7[i] = new TH1I(name, name, 1024, 0, 50.);
146 ListAmpTofPlane7->Add(hAmpTofPlane7[i]);
147 name = Form("hTimeTofPlane7_%d", i);
148 hTimeTofPlane7[i] = new TH1I(name, name, 1000, 0, 1000.);
149 ListTimeTofPlane7->Add(hTimeTofPlane7[i]);
150 name = Form("hLeftRightCorelationRPlane7_%d", i);
151 hLeftRightCorelationRPlane7[i] = new TH2S(name, name, 1024, 0., 50., 1024, 0., 50.);
152 ListLRCorrTofPlane7->Add(hLeftRightCorelationRPlane7[i]);
153 name = Form("hLeftMinusRightPlane7_%d", i);
154 hLeftMinusRightPlane7[i] = new TH1I(name, name, 1024, -12., 12.);
155 ListLRTofPlane7->Add(hLeftMinusRightPlane7[i]);
156 name = Form("hLeftRightVsAmpLPlane7_%d", i);
157 hLeftRightVsAmpLPlane7[i] = new TH2S(name, name, 1024, 0., 50., 1024, -12., 12.);
158 ListLRvsAmpLTofPlane7->Add(hLeftRightVsAmpLPlane7[i]);
159 name = Form("hLeftRightVsAmpRPlane7_%d", i);
160 hLeftRightVsAmpRPlane7[i] = new TH2S(name, name, 1024, 0., 50., 1024, -12., 12.);
161 ListLRvsAmpRTofPlane7->Add(hLeftRightVsAmpRPlane7[i]);
162 }
163 TH2I * hAmpLRPlane7 = new TH2I("hAmpLRPlane7", "hAmpLRPlane7", 1024, 0., 50., 1024, 0., 50.);
164 ListAmpTofPlane7->Add(hAmpLRPlane7);
165
166 TH1I *hHitStripsPerEvPlane2 = new TH1I("hHitStripsPerEvPlane2", "Hit Strips Per Ev on Plane#2", kNST, 0, kNST);
167 ListStat->Add(hHitStripsPerEvPlane2);
168 TH1I *hHitStripsPerEvPlane7 = new TH1I("hHitStripsPerEvPlane7", "Hit Strips Per Ev on Plane#7", kNST, 0, kNST);
169 ListStat->Add(hHitStripsPerEvPlane7);
170 TH2I *hLefRightCorelation_Plane2 = new TH2I("hLefRightCorelation_Plane2", " Lef-Right Corelation Plane#2", kNST, 0, kNST, kNST, 0, kNST);
171 ListStat->Add(hLefRightCorelation_Plane2);
172 TH2I *hLefRightCorelation_Plane7 = new TH2I("hLefRightCorelation_Plane7", " Lef-Right Corelation Plane#7", kNST, 0, kNST, kNST, 0, kNST);
173 ListStat->Add(hLefRightCorelation_Plane7);
174 TH2I *hPlane2Plane7Correlation = new TH2I("hPlane2Plane7Correlation", " Plane2-Plane7 Correlation", kNST, 0, kNST, kNST, 0, kNST);
175 ListStat->Add(hPlane2Plane7Correlation);
176 TH1S *hNoize = new TH1S("hNoize", "hNoize", 500, 0, 50);
177 hNoize->GetXaxis()->SetTitle("MHz");
178 ListStat->Add(hNoize);
179 Int_t NEvStrip47 = 0;
180 TH1I *hNEvStrip47 = new TH1I("hNEvStrip47", "hNEvStrip47", 20, 0, 20);
181 ListStat->Add(hNEvStrip47);
182
183 TList * ListT0 = new TList();
184 TH1I *hHitT0PerEv = new TH1I("hHitT0PerEv", "Hit T0 Per Ev", 10, 0, 10);
185 ListT0->Add(hHitT0PerEv);
186 TH1I *hAmpT0 = new TH1I("hAmpT0", "Amp T0", 1024, 0, 50);
187 ListT0->Add(hAmpT0);
188 TH1I *hTimeT0 = new TH1I("hTimeT0", "Time T0", 1000, 0, 1000);
189 ListT0->Add(hTimeT0);
190
191 TH1I * hdt[kNST + 1];
192
193 for (Int_t i = 0; i < kNST + 1; i++) {
194 name = Form("hdt_%d", i);
195 name2 = Form("dt T0 - ToF400_Plane7_str_%d", i);
196 hdt[i] = new TH1I(name, name2, 1024 + 512, -18., 18.);
197 hdt[i]->SetLineWidth(3);
198 ListAnalysis->Add(hdt[i]);
199 }
200
201 TH2I * hdt_vs_AmpT0[kNST + 1];
202 for (Int_t i = 0; i < kNST + 1; i++) {
203 name = Form("hdt_vs_AmpT0_%d", i);
204 name2 = Form("dt T0 - ToF400_Plane7_str_%d vs AMP_T0", i);
205 hdt_vs_AmpT0[i] = new TH2I(name, name2, 1024, 0, 50, 1024 + 512, -24, 48);
206 ListAnalysis->Add(hdt_vs_AmpT0[i]);
207 }
208
209 TH2I * hdt_vs_AmpToF[kNST + 1];
210 for (Int_t i = 0; i < kNST + 1; i++) {
211 name = Form("hdt_vs_AmpToF_%d", i);
212 name2 = Form("dt T0 - ToF400_Plane7_str_%d vs AMP_TOF400", i);
213 hdt_vs_AmpToF[i] = new TH2I(name, name2, 1024, 0, 50, 1024 + 512, -24, 48);
214 ListAnalysis->Add(hdt_vs_AmpToF[i]);
215 }
216
217 TH2F *h_StripHit = new TH2F("h_StripHit", "h_StripHit", kNST, 0, kNST, 2, 0, 2);
218 TH1I *h_StripHitSum = new TH1I("h_StripHitSum", "h_StripHitSum", kNST, 0, kNST);
219 TH1F *h_StripAmp = new TH1F("h_StripAmp", "h_StripAmp", kNST, 0, kNST);
220 h_StripAmp->GetYaxis()->SetLabelSize(0.07);
221 h_StripAmp->SetLineWidth(2);
222 h_StripAmp->SetLineColor(46);
223 TH1F *h_StripTime = new TH1F("h_StripTime", "h_StripTime", kNST, 0, kNST);
224 h_StripTime->GetYaxis()->SetLabelSize(0.07);
225 h_StripTime->GetYaxis()->SetRangeUser(-1., 1.);
226 h_StripTime->SetLineWidth(2);
227 h_StripTime->SetLineColor(46);
228 TH1F *h_StripAmpL = new TH1F("h_StripAmpL", "h_StripAmpL", kNST, 0, kNST);
229 h_StripAmpL->GetYaxis()->SetLabelSize(0.07);
230 h_StripAmpL->SetLineWidth(2);
231 h_StripAmpL->SetLineColor(2);
232 TH1F *h_StripTimeL = new TH1F("h_StripTimeL", "h_StripTimeL", kNST, 0, kNST);
233 h_StripTimeL->GetYaxis()->SetLabelSize(0.07);
234 h_StripTimeL->GetYaxis()->SetRangeUser(300, 450);
235 h_StripTimeL->SetLineWidth(2);
236 h_StripTimeL->SetLineColor(2);
237 TH1F *h_StripAmpR = new TH1F("h_StripAmpR", "h_StripAmpR", kNST, 0, kNST);
238 h_StripAmpR->GetYaxis()->SetLabelSize(0.07);
239 h_StripAmpR->SetLineWidth(2);
240 h_StripAmpR->SetLineColor(9);
241 TH1F *h_StripTimeR = new TH1F("h_StripTimeR", "h_StripTimeR", kNST, 0, kNST);
242 h_StripTimeR->GetYaxis()->SetLabelSize(0.07);
243 h_StripTimeR->GetYaxis()->SetRangeUser(300, 450);
244 h_StripTimeR->SetLineWidth(2);
245 h_StripTimeR->SetLineColor(9);
246 TH1F *h_StripLeftMinusRight = new TH1F("h_StripLeftMinusRight", "h_StripLeftMinusRight", kNST, 0, kNST);
247 h_StripLeftMinusRight->GetYaxis()->SetLabelSize(0.07);
248 h_StripLeftMinusRight->GetYaxis()->SetRangeUser(-2.5, 2.5);
249 h_StripLeftMinusRight->SetLineWidth(2);
250 h_StripLeftMinusRight->SetLineColor(9);
251 TH2F *h_StripLeftMinusRight_vs_Time = new TH2F("h_StripLeftMinusRight_vs_Time", "h_StripLeftMinusRight_vs_Time", 512, -5., 5., 2048, 250., 450.);
252 h_StripLeftMinusRight_vs_Time->GetYaxis()->SetLabelSize(0.07);
253 h_StripLeftMinusRight_vs_Time->GetXaxis()->SetLabelSize(0.07);
254 h_StripLeftMinusRight_vs_Time->SetMarkerStyle(8);
255
256 TH2F *h_XYTime = new TH2F("h_XYTime", "h_XYTime", kNST, 0, kNST, 64, -16, 16.);
257 h_XYTime->GetYaxis()->SetLabelSize(0.07);
258 h_XYTime->GetXaxis()->SetLabelSize(0.07);
259 //h_XYTime->SetMarkerStyle(8);
260
261 //--------Data Analysis-----------------------------------------------------------------
262
263 TCanvas *c_1 = new TCanvas("c1", "c1", 20, 20, 1500, 750);
264 TPad *Pad_StripHit = new TPad("Pad_StripHit", "Pad_StripHit", 0.01, 0.51, 0.49, 0.99, 21);
265 TPad *Pad_StripHitSum = new TPad("Pad_StripHitSum", "Pad_StripHitSum", 0.01, 0.01, 0.49, 0.49, 21);
266 TPad *Pad_StripAmp = new TPad("Pad_StripAmp", "Pad_StripAmp", 0.51, 0.51, 0.99, 0.99, 21);
267 TPad *Pad_StripTime = new TPad("Pad_StripTime", "Pad_StripTime", 0.51, 0.01, 0.99, 0.49, 21);
268
269 Pad_StripHit->Draw();
270 Pad_StripHitSum->Draw();
271 Pad_StripAmp->Draw();
272 Pad_StripAmp->SetGridy();
273 Pad_StripTime->Draw();
274 Pad_StripTime->SetGridy();
275
276 TString NameDTFile = file;
277 Point = NameDTFile.First('.');
278 NameDTFile.Replace(Point, 15, "_dt.dat");
279 cout << "Write call to " << NameDTFile.Data() << endl;
280 f_dt.open(NameDTFile.Data());
281
282 for (Int_t iEv = 0; iEv < nEvForRead; iEv++) {
283 if (iEv % 50000 == 0) cout << "EVENT: " << iEv << endl;
284 eveTree->GetEntry(iEv);
285 CleanTof(&Plane2);
286 CleanTof(&Plane7);
287 CleanBC2(&T0);
288 h_StripHit->Reset();
289 h_StripAmp->Reset();
290 h_StripTime->Reset();
291 h_StripAmpL->Reset();
292 h_StripTimeL->Reset();
293 h_StripAmpR->Reset();
294 h_StripTimeR->Reset();
295 h_StripLeftMinusRight->Reset();
296 h_StripLeftMinusRight_vs_Time->Reset();
297 h_XYTime->Reset();
298 NEvStrip47 = 0;
299 for (Int_t iDig = 0; iDig < ToF400Digits->GetEntriesFast(); ++iDig) {
300 BmnTof1Digit* digTof = (BmnTof1Digit*) ToF400Digits->At(iDig);
301 if (digTof->GetPlane() == 2) {
302 Int_t strip = digTof->GetStrip();
303 Int_t plane = digTof->GetPlane();
304 // cout << "Plane = " << digTof->GetPlane() << " Ev# " << iEv << endl;
305 if (digTof->GetSide() == 0 && Plane2.sStrHit[strip] == kFALSE) {
306 Plane2.sAmpL[strip] = digTof->GetAmplitude();
307 Plane2.sTimeL[strip] = digTof->GetTime();
308 //h_StripHit->SetBinContent(strip + 1, digTof->GetSide() + 1, digTof->GetAmplitude());
309 }
310 if (digTof->GetSide() == 1 && Plane2.sStrHit[strip] == kFALSE) {
311 Plane2.sAmpR[strip] = digTof->GetAmplitude();
312 Plane2.sTimeR[strip] = digTof->GetTime();
313 //h_StripHit->SetBinContent(strip + 1, digTof->GetSide() + 1, digTof->GetAmplitude());
314 }
315
316 if ((Plane2.sTimeR[strip] != 0 && Plane2.sTimeL[strip] != 0
317 && strip != 47
318 && TMath::Abs((Plane2.sTimeL[strip] - Plane2.sTimeR[strip]) * 0.5) <= 2.
319 && TMath::Abs((Plane2.sAmpL[strip] - Plane2.sAmpR[strip]) * 0.5) <= 1.5
320 && Plane2.sStrHit[strip] == kFALSE)) //
321 Plane2.sStrHit[strip] = kTRUE; //*/
322 }
323
324 if (digTof->GetPlane() == 7) {
325 Int_t strip = digTof->GetStrip();
326 Int_t plane = digTof->GetPlane();
327 // cout << "Plane = " << digTof->GetPlane() << " Ev# " << iEv << endl;
328 if (digTof->GetSide() == 0 && Plane7.sStrHit[strip] == kFALSE) {
329 //if (digTof->GetTime() > 300. && digTof->GetTime() < 430.) {
330 Plane7.sAmpL[strip] = digTof->GetAmplitude();
331
332 if (strip == 0) {
333 NEvStrip47++;
334 if (Plane7.sTimeL[strip] != 0)
335 hNoize->Fill(1000. / (digTof->GetTime() - Plane7.sTimeL[strip]));
336 }
337
338 Plane7.sTimeL[strip] = digTof->GetTime();
339 h_StripHit->SetBinContent(strip + 1, digTof->GetSide() + 1, digTof->GetAmplitude());
340 //if (strip == 22) cout << "L ";
341 }
342 if (digTof->GetSide() == 1 && Plane7.sStrHit[strip] == kFALSE) {
343 // if (digTof->GetTime() > 300. && digTof->GetTime() < 430.) {
344 Plane7.sAmpR[strip] = digTof->GetAmplitude();
345 Plane7.sTimeR[strip] = digTof->GetTime();
346 h_StripHit->SetBinContent(strip + 1, digTof->GetSide() + 1, digTof->GetAmplitude());
347 //if (strip == 22) cout << "R ";
348 }
349
350 if ((Plane7.sTimeR[strip] != 0 && Plane7.sTimeL[strip] != 0
351 && strip != 47
352 && TMath::Abs((Plane7.sTimeL[strip] - Plane7.sTimeR[strip]) * 0.5) <= 2.
353 && TMath::Abs((Plane7.sAmpL[strip] - Plane7.sAmpR[strip]) * 0.5) <= 1.5
354 && Plane7.sStrHit[strip] == kFALSE)) //
355 Plane7.sStrHit[strip] = kTRUE; //*/
356 }
357 }
358 hNEvStrip47->Fill(NEvStrip47);
359
360 for (Int_t i = 0; i < kNST; i++) {
361 if (Plane2.sStrHit[i] == kFALSE || i == 0 || i == 47) {
362 Plane2.sAmpL[i] = 0;
363 Plane2.sAmpR[i] = 0;
364 Plane2.sTimeL[i] = 0;
365 Plane2.sTimeR[i] = 0;
366 Plane2.sStrHit[i] = kFALSE;
367 }
368 if (Plane7.sStrHit[i] == kFALSE || i == 0 || i == 47) {
369 Plane7.sAmpL[i] = 0;
370 Plane7.sAmpR[i] = 0;
371 Plane7.sTimeL[i] = 0;
372 Plane7.sTimeR[i] = 0;
373 Plane7.sStrHit[i] = kFALSE;
374 }
375 }
376
377 for (Int_t iDig = 0; iDig < BC2Digits->GetEntriesFast(); ++iDig) {
378 BmnTrigDigit* digBC2 = (BmnTrigDigit*) BC2Digits->At(iDig);
379 //if (digBC2->GetAmp() > 6.15 && digBC2->GetAmp() < 17.14)
380 {
381 T0.sNHits++;
382 if (T0.sAmp == 0) {
383 T0.sAmp = digBC2->GetAmp();
384 T0.sTime = digBC2->GetTime();
385 hAmpT0->Fill(T0.sAmp);
386 hTimeT0->Fill(T0.sTime);
387 }
388 }
389 }
390 hHitT0PerEv->Fill(T0.sNHits);
391
392 for (Int_t i = 0; i < kNST; i++)
393 for (Int_t j = 0; j < kNST; j++) {
394 if (Plane2.sAmpL[i] != 0 && Plane2.sAmpR[j] != 0)
395 hLefRightCorelation_Plane2->Fill(i, j);
396 if (Plane7.sAmpL[i] != 0 && Plane7.sAmpR[j] != 0)
397 hLefRightCorelation_Plane7->Fill(i, j);
398 }
399
400 for (Int_t i = 0; i < kNST; i++) {
401 if (Plane2.sAmpL[i] != 0 && Plane2.sAmpR[i] != 0) {
402 Plane2.sAmp[i] = Plane2.sAmpL[i] + Plane2.sAmpR[i];
403 Plane2.sTime[i] = (Plane2.sTimeL[i] + Plane2.sTimeR[i]) * 0.5;
404 Plane2.sNHits++;
405 h_StripHitSum->Fill(i);
406 hAmpTofPlane2[i]->Fill(Plane2.sAmp[i]);
407 hTimeTofPlane2[i]->Fill(Plane2.sTime[i]);
408 hLevtMinusRightPlane2[i]->Fill((Plane2.sTimeL[i] - Plane2.sTimeR[i]) * 0.5);
409 hAmpTofPlane2[kNST]->Fill(Plane2.sAmp[i]);
410 hTimeTofPlane2[kNST]->Fill(Plane2.sTime[i]);
411 hLevtMinusRightPlane2[kNST]->Fill((Plane2.sTimeL[i] - Plane2.sTimeR[i]) * 0.5);
412 hAmpLRPlane2->Fill(Plane2.sAmpL[i], Plane2.sAmpR[i]);
413
414 /* h_StripAmp->Fill(i + 0.5, Plane2.sAmp[i]);
415 h_StripTime->Fill(i + 0.5, Plane2.sTime[i]);
416 h_StripAmpL->Fill(i + 0.5, Plane2.sAmpL[i]);
417 h_StripTimeL->Fill(i + 0.5, Plane2.sTimeL[i]);
418 h_StripAmpR->Fill(i + 0.5, Plane2.sAmpR[i]);
419 h_StripTimeR->Fill(i + 0.5, Plane2.sTimeR[i]);
420 h_StripLeftMinusRight->Fill(i + 0.5, Plane2.sTimeL[i] - Plane2.sTimeR[i]);
421 h_StripLeftMinusRight_vs_Time->Fill(Plane2.sTimeL[i] - Plane2.sTimeR[i], Plane2.sTime[i], Plane2.sAmp[i]);
422 //*/
423 }
424 //if (Plane7.sAmpL[i] != 0 && Plane7.sAmpR[i] != 0) {
425 if (Plane7.sStrHit[i] != 0) {
426 Plane7.sAmp[i] = Plane7.sAmpL[i] + Plane7.sAmpR[i];
427 Plane7.sTime[i] = (Plane7.sTimeL[i] + Plane7.sTimeR[i]) / 2.;
428 Plane7.sNHits++;
429 h_StripHitSum->Fill(i);
430 hAmpTofPlane7[i]->Fill(Plane7.sAmp[i]);
431 hTimeTofPlane7[i]->Fill(Plane7.sTime[i]);
432 hLeftMinusRightPlane7[i]->Fill((Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5);
433 hAmpTofPlane7[kNST]->Fill(Plane7.sAmp[i]);
434 hTimeTofPlane7[kNST]->Fill(Plane7.sTime[i]);
435 hLeftMinusRightPlane7[kNST]->Fill((Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5);
436 hAmpLRPlane7->Fill(Plane7.sAmpL[i], Plane7.sAmpR[i]);
437
438 hLeftRightCorelationRPlane7[i]->Fill(Plane7.sAmpL[i], Plane7.sAmpR[i]);
439 hLeftRightCorelationRPlane7[kNST]->Fill(Plane7.sAmpL[i], Plane7.sAmpR[i]);
440
441 hLeftRightVsAmpLPlane7[i]->Fill(Plane7.sAmpL[i], (Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5);
442 hLeftRightVsAmpLPlane7[kNST]->Fill(Plane7.sAmpL[i], (Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5);
443
444 hLeftRightVsAmpRPlane7[i]->Fill(Plane7.sAmpR[i], (Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5);
445 hLeftRightVsAmpRPlane7[kNST]->Fill(Plane7.sAmpR[i], (Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5);
446
447 /* h_StripAmp->Fill(i + 0.5, Plane7.sAmp[i]);
448 h_StripTime->Fill(i + 0.5, CalculateDt(&Plane7, &T0, i));
449 h_StripAmpL->Fill(i + 0.5, Plane7.sAmpL[i]);
450 h_StripTimeL->Fill(i + 0.5, Plane7.sTimeL[i]);
451 h_StripAmpR->Fill(i + 0.5, Plane7.sAmpR[i]);
452 h_StripTimeR->Fill(i + 0.5, Plane7.sTimeR[i]);
453 h_StripLeftMinusRight->Fill(i + 0.5, (Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5);
454 //h_StripLeftMinusRight_vs_Time->Fill((Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5, Plane7.sTime[i], Plane7.sAmp[i]);
455 h_XYTime->Fill(i, ((Plane7.sTimeL[i] - Plane7.sTimeR[i]) * 0.5 / 0.06));
456
457 //*/
458 }
459 }
460 hHitStripsPerEvPlane2->Fill(Plane2.sNHits);
461 hHitStripsPerEvPlane7->Fill(Plane7.sNHits);
462
463 if (Plane2.sNHits == 1 && Plane2.sAmp[5] != 0 && T0.sNHits == 1) {
464 HitPlane2T0++;
465 if (Plane7.sNHits != 0) EffPlane7++;
466 }
467 if (Plane7.sNHits == 1 && Plane7.sAmp[42] != 0 && T0.sNHits == 1) {
468 HitPlane7T0++;
469 if (Plane2.sNHits != 0) EffPlane2++;
470 }
471 if (Plane2.sNHits != 0 && Plane7.sNHits != 0 && T0.sNHits == 1) {
472 Int_t MaxStrPlane2 = TMath::LocMax(kNST, Plane2.sAmp);
473 Int_t MaxStrPlane7 = TMath::LocMax(kNST, Plane7.sAmp);
474 hPlane2Plane7Correlation->Fill(MaxStrPlane2, MaxStrPlane7);
475 }
476
477 Double_t dt = -100.;
478 Int_t MaxStr = -1;
479 if (Plane7.sNHits == 1 && T0.sNHits == 1
480 && T0.sAmp >= 17.9 && T0.sAmp < 19.63) {
481 MaxStr = TMath::LocMax(kNST, Plane7.sAmp);
482 //dt = Plane7->sTime[MaxStr] - T0->sTime;
483 dt = CalculateDt(&Plane7, &T0, MaxStr);
484 hdt[MaxStr]->Fill(dt);
485 hdt_vs_AmpT0[MaxStr]->Fill(T0.sAmp, dt);
486 hdt_vs_AmpToF[MaxStr]->Fill(Plane7.sAmpL[MaxStr], dt);
487 hdt[kNST]->Fill(dt);
488 hdt_vs_AmpT0[kNST]->Fill(T0.sAmp, dt);
489 hdt_vs_AmpToF[kNST]->Fill(Plane7.sAmpL[MaxStr], dt);
490 }
491 Int_t iEvDig = EventHeader->GetEntriesFast();
492 if (iEvDig != 1) cout << "iEvDig == " << iEvDig << endl;
493 BmnEventHeader* digEvent = (BmnEventHeader*) EventHeader->At(0);
494 f_dt << (digEvent->GetEventId()) << "\t" << MaxStr << "\t" << dt << endl;
495
496 /* if (Plane7.sNHits != 0) {
497 // if (Plane7.sNHits > 2) {
498 cout << "Event# " << digEvent->GetEventId() << endl;
499 Int_t MaxStr = TMath::LocMax(kNST, Plane7.sTime);
500 //h_StripTime->GetYaxis()->SetRangeUser(-15., +15.);
501 //h_StripTime->GetYaxis()->SetRangeUser(Plane7.sTime[MaxStr] - 15., Plane7.sTime[MaxStr] + 5.);
502 // cout << "event # " << iEv;
503 Pad_StripHit->cd();
504 h_StripHit->Draw("COLZ");
505 Pad_StripHitSum->cd();
506 //h_StripHitSum->Draw ();
507 h_StripLeftMinusRight->Draw();
508 Pad_StripAmp->cd();
509 //h_StripLeftMinusRight_vs_Time->Draw();
510 h_XYTime->Draw("COLZ");
511 //h_StripAmp->Draw ("");
512 // h_StripAmpL->Draw("SAME");
513 // h_StripAmpR->Draw("SAME");
514
515 Pad_StripTime->cd();
516 h_StripTime->Draw("");
517 // h_StripTimeL->Draw("SAME");
518 // h_StripTimeR->Draw("SAME");
519 c_1->Update();
520 getchar();
521 }//*/
522 }
523
524 f_dt.close();
525
526 TString NameCallFile = file;
527 Point = NameCallFile.First('.');
528 NameCallFile.Replace(Point, 15, "_call.dat");
529 cout << "Write call to " << NameCallFile.Data() << endl;
530 f_call.open(NameCallFile.Data());
531 //TProfile * pr_Hist_dt_FFD_SRPC200_vs_WSRPC[3];
532 //TProfile * pr_Hist_dt_FFD_SRPC400_vs_WSRPC[3];
533 for (i = 0; i < kNST + 1; i++) {
534
535 hAmpTofPlane7[i]->Fit(f_gaus, "QW", "QW", hAmpTofPlane7[i]->GetMaximumBin() * hAmpTofPlane7[i]->GetXaxis()->GetBinWidth(1) - 0.5, hAmpTofPlane7[i]->GetMaximumBin() * hAmpTofPlane7[i]->GetXaxis()->GetBinWidth(1) + 0.8);
536 sprintf(hid, "Plane7 Strip %d\t mean %f\t integral %.0f", i, f_gaus->GetParameter(1), hAmpTofPlane7[i]->GetEntries());
537 f_call << hid << endl;
538 }
539 if (HitPlane7T0 != 0 && HitPlane2T0 != 0) {
540 f_call << "Efficiency Plane2 = " << EffPlane2 << "/" << HitPlane7T0 << "=" << EffPlane2 * 100. / HitPlane7T0 << " %" << endl;
541 cout << "Efficiency Plane2 = " << EffPlane2 << "/" << HitPlane7T0 << "=" << EffPlane2 * 100. / HitPlane7T0 << " %" << endl;
542 f_call << "Efficiency Plane7 = " << EffPlane7 << "/" << HitPlane2T0 << "=" << EffPlane7 * 100. / HitPlane2T0 << " %" << endl;
543 cout << "Efficiency Plane7 = " << EffPlane7 << "/" << HitPlane2T0 << "=" << EffPlane7 * 100. / HitPlane2T0 << " %" << endl;
544 }
545 f_call << "Strip\t Mean\t Sigma" << endl;
546 for (i = 0; i < kNST + 1; i++) {
547 f_gaus->SetParameter(1, 0.);
548 f_gaus->SetParameter(2, 0.1);
549 hdt[i]->Fit(f_gaus, "QW", "QW", -1, 1);
550 f_call << i << "\t" << f_gaus->GetParameter(1) << "\t" << f_gaus->GetParameter(2) * 1000. << endl;
551 cout << "Strip " << i << "\tMean " << f_gaus->GetParameter(1) << "\tSigma" << f_gaus->GetParameter(2) * 1000. << endl;
552 }
553 f_call.close();
554
555 TString outName = file;
556 Point = outName.First('.');
557 outName.Replace(Point, 15, "_TofAn.root");
558 TFile *fileout = new TFile(outName.Data(), "RECREATE");
559
560 TDirectory * Dir;
561 Dir = fileout->mkdir("ToF");
562 Dir->cd();
563
564 TDirectory * DirPlane2;
565 DirPlane2 = Dir->mkdir("Plane2");
566 DirPlane2->cd();
567 ListTimeTofPlane2->Write();
568 ListAmpTofPlane2->Write();
569 ListLRTofPlane2->Write();
570
571 TDirectory * DirPlane7;
572 TDirectory * DirLRPlane7;
573 DirPlane7 = Dir->mkdir("Plane7");
574 DirLRPlane7 = DirPlane7->mkdir("LRCorelation");
575 DirPlane7->cd();
576 ListTimeTofPlane7->Write();
577 ListAmpTofPlane7->Write();
578 DirLRPlane7->cd();
579 ListLRCorrTofPlane7->Write();
580 ListLRTofPlane7->Write();
581 ListLRvsAmpLTofPlane7->Write();
582 ListLRvsAmpRTofPlane7->Write();
583
584 TDirectory * DirT0;
585 DirT0 = fileout->mkdir("T0");
586 DirT0->cd();
587 ListT0->Write();
588
589 TDirectory * DirAnalysis;
590 DirAnalysis = fileout->mkdir("Analysis");
591 DirAnalysis->cd();
592 ListAnalysis->Write();
593
594 fileout->cd();
595 ListStat->Write();
596 fileout->Close();
597
598}//end of macros
599
600//----Additional functions ----------------------------------------------------------------
601
602void
604 for (Int_t st = 0; st < kNST; st++) {
605 Det->sAmpL[st] = 0.;
606 Det->sAmpR[st] = 0.;
607 Det->sAmp[st] = 0.;
608 Det->sTimeL[st] = 0.;
609 Det->sTimeR[st] = 0.;
610 Det->sTime[st] = 0.;
611 Det->sStrHit[st] = kFALSE;
612 }
613 Det->sNHits = 0;
614}
615
616void
618 Det->sTime = 0;
619 Det->sAmp = 0;
620 Det->sNHits = 0;
621}
622
623Double_t CalculateDt(ToF400Detector * Det, T0Detector * T0, Int_t Str = -1) {
624
625 if (Str == -1) Str = TMath::LocMax(kNST, Plane7->sAmp);
626 Double_t dt = Det->sTime[Str] - T0->sTime;
627 //cout << " dt = " << dt << endl;
628 /* if (Det->sAmp[Str] < 20.46) CorrPlane7_It1 = 0;
629 else if (Det->sAmp[Str] >= 20.46 && Det->sAmp[Str] < 27.25) CorrPlane7_It1 = 1;
630 else if (Det->sAmp[Str] >= 27.25 && Det->sAmp[Str] < 34.81) CorrPlane7_It1 = 2;
631 else if (Det->sAmp[Str] >= 34.81 && Det->sAmp[Str] < 37.11) CorrPlane7_It1 = 3;
632 else if (Det->sAmp[Str] >= 37.11) CorrPlane7_It1 = 4;
633 dt = dt - (CorrPlane7Coeff_It1[CorrPlane7_It1][0] + CorrPlane7Coeff_It1[CorrPlane7_It1][1] * Det->sAmp[Str] +
634 CorrPlane7Coeff_It1[CorrPlane7_It1][2] * Det->sAmp[Str] * Det->sAmp[Str] +
635 CorrPlane7Coeff_It1[CorrPlane7_It1][3] * Det->sAmp[Str] * Det->sAmp[Str] * Det->sAmp[Str]); //*/
636
637 /* if (T0->sAmp >= 6.15 && T0->sAmp < 12.3) CorrT0_It1 = 0; //deutron
638 else if (T0->sAmp >= 12.3 && T0->sAmp < 17.14) CorrT0_It1 = 1; //deutron
639 else if (T0->sAmp >= 17.9 && T0->sAmp < 19.63) CorrT0_It1 = 2; //C6+
640 else CorrT0_It1 = 3; //skip
641 dt = dt - (CorrT0Coeff_It1[CorrT0_It1][0] + CorrT0Coeff_It1[CorrT0_It1][1] * T0->sAmp +
642 CorrT0Coeff_It1[CorrT0_It1][2] * T0->sAmp * T0->sAmp +
643 CorrT0Coeff_It1[CorrT0_It1][3] * T0->sAmp * T0->sAmp * T0->sAmp); //*/
644 return dt;
645}
int i
Definition P4_F32vec4.h:22
UInt_t GetEventId()
Short_t GetPlane() const
Short_t GetSide() const
Short_t GetStrip() const
Float_t GetTime() const
Float_t GetAmplitude() const
Double_t GetAmp() const
Double_t GetTime() const
Double_t CatTimeBC2Low
Double_t CorrPlane7Coeff_It1[5][4]
const UInt_t kNST
void CleanBC2(T0Detector *Det)
const Double_t kTimeBin
Double_t CatTimeBC2High
void ToF400DigitAnalysis_period5(TString file="", Int_t nEvForRead=0, Int_t Periud=5)
void CleanTof(ToF400Detector *Det)
Double_t CalculateDt(ToF400Detector *Det, T0Detector *T0, Int_t Str=-1)
Double_t CorrT0Coeff_It1[4][4]
Int_t CorrPlane7_It1
Bool_t FlagDigit[10][kNST]
const UInt_t kNST