Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
//
// Make a set of histograms containing raw adcs and tdcs for
// each PMT of the HMS hodoscopes
// TDCs are after application of a common offset and
// ADCs are after subtraction of a pedestal
//
{
#define NPLANES 4
#define NSIDES 2
#define NSIGNALS 2
#define MAXBARS 16
TString SPECTROMETER="H";
TString DETECTOR="hod";
TString plane_names[NPLANES]={"1x","1y","2x","2y"};
Int_t nbars[NPLANES]={16,10,16,10};
TString sides[NSIDES]={"neg","pos"};
TString signals[NSIGNALS]={"adc","tdc"};
Double_t ADC_MIN=0.0;
Double_t ADC_MAX=4000.0;
Int_t ADC_NBINS=200;
Double_t TDC_MIN=0.0;
Double_t TDC_MAX=4000.0;
Int_t TDC_NBINS=200;
Int_t nhits[NPLANES][NSIDES][NSIGNALS];
Double_t paddles[NPLANES][NSIDES][NSIGNALS][MAXBARS];
Double_t values[NPLANES][NSIDES][NSIGNALS][MAXBARS];
TH1F* h[NPLANES*NSIDES*NSIGNALS*MAXBARS];
TFile file("hodtest.root","READ");
TFile hfile("hhodrawhists.root", "RECREATE");
file.cd();
TTree* T = (TTree*) file.Get("T");
for(Int_t ip=0;ip<NPLANES;ip++){
for(Int_t iside=0;iside<NSIDES;iside++) {
for(Int_t isignal=0;isignal<NSIGNALS;isignal++) {
TString base_name = SPECTROMETER+"."+DETECTOR+"."+
plane_names[ip]+"."+sides[iside]+signals[isignal];
TString ndata_name = "Ndata."+base_name+"pad";
TString padlist_name = base_name+"pad";
TString vallist_name = base_name+"val";
// Set branch addresses
T->SetBranchAddress(ndata_name, &nhits[ip][iside][isignal]);
T->SetBranchAddress(padlist_name, &paddles[ip][iside][isignal][0]);
T->SetBranchAddress(vallist_name, &values[ip][iside][isignal][0]);
// Create histograms
// ADC and TDC histogram for each
Int_t nbins,hmin,hmax;
if(signals[isignal]=="adc") {
nbins = ADC_NBINS;
hmin = ADC_MIN;
hmax = ADC_MAX;
} else {
nbins = TDC_NBINS;
hmin = TDC_MIN;
hmax = TDC_MAX;
}
hfile.cd();
for(Int_t ibar=0;ibar<nbars[ip];ibar++) {
Int_t hindex = ip*NSIDES*NSIGNALS*MAXBARS
+iside*NSIGNALS*MAXBARS+isignal*MAXBARS+ibar;
TString ibarname=Form("%d",ibar+1);
TString title= "hS"+plane_names[ip]+ibarname+sides[iside]
+" "+signals[isignal];
TString name="h"; name += Form("%d",hindex);
// cout << name << " " << title << " " << nbins << " " << hmin << " "
// << hmax << endl;
h[hindex]=new TH1F(name,title,nbins, hmin, hmax);
}
file.cd();
}
}
}
// Loop over the events, filling the histograms
for(Int_t iev=0, N=T->GetEntries(); iev < N; iev++) {
T->GetEntry(iev);
for(Int_t ip=0;ip<NPLANES;ip++){
for(Int_t iside=0;iside<NSIDES;iside++) {
for(Int_t isignal=0;isignal<NSIGNALS;isignal++) {
Int_t hindex_base = ip*NSIDES*NSIGNALS*MAXBARS
+iside*NSIGNALS*MAXBARS+isignal*MAXBARS;
for(Int_t ihit=0;ihit<nhits[ip][iside][isignal];ihit++) {
Int_t ibar=TMath::Nint(paddles[ip][iside][isignal][ihit])-1;
Double_t val = values[ip][iside][isignal][ihit];
Int_t hindex = hindex_base + ibar;
// cout << ip << " " << iside << " " << isignal << " " << ihit
// << " " << ibar << " " << hindex << endl;
h[hindex]->Fill(val);
}
}
}
}
// Display or save the histograms
}
hfile.cd();
hfile.Write();
}