19 void plot(std::vector<TH2F*> Map,
const sinfo& surface_info,
const std::string&
name){
21 std::string out_name = name+
"/"+surface_info.
name+
"/"+surface_info.
name+
"_"+surface_info.
idname;
22 gSystem->Exec( Form(
"mkdir %s", (name+
"/"+surface_info.
name).c_str()) );
25 if(surface_info.
type == 2){
27 TText *vol =
new TText(.1,.95,surface_info.
name.c_str());
29 TText *
surface =
new TText(.1,.9,surface_info.
id.c_str());
31 TText *surface_z =
new TText(.1,.85,(
"Z = " +
to_string(surface_info.
pos)).c_str() );
34 TCanvas *
c1 =
new TCanvas(
"c1",
"mat_X0",1200,1200);
35 c1->SetRightMargin(0.14);
36 c1->SetTopMargin(0.14);
37 c1->SetLeftMargin(0.14);
38 c1->SetBottomMargin(0.14);
43 c1->Print( (out_name+
"_X0.pdf").c_str());
54 if(surface_info.
type == 1){
56 TText *vol =
new TText(.1,.95,surface_info.
name.c_str());
58 TText *
surface =
new TText(.1,.9,surface_info.
id.c_str());
60 TText *surface_r =
new TText(.1,.85,(
"R = " +
to_string(surface_info.
pos)).c_str() );
62 TCanvas *
c1 =
new TCanvas(
"c1",
"mat_X0",1200,1200);
63 c1->SetRightMargin(0.14);
64 c1->SetTopMargin(0.14);
65 c1->SetLeftMargin(0.14);
66 c1->SetBottomMargin(0.14);
71 c1->Print( (out_name+
"_X0.pdf").c_str());
86 const sinfo& surface_info){
93 if(surface_info.
type == 1){
94 Map_X0 =
new TH2F((
"Map_X0_"+surface_info.
idname).c_str(),(
"Map_X0_"+surface_info.
idname).c_str(),
96 Map_L0 =
new TH2F((
"Map_L0_"+surface_info.
idname).c_str(),(
"Map_L0_"+surface_info.
idname).c_str(),
98 Map_scale =
new TH2F((
"Map_scale_"+surface_info.
idname).c_str(),(
"Map_scale_"+surface_info.
idname).c_str(),
100 Map_X0->GetXaxis()->SetTitle(
"Eta");
101 Map_X0->GetYaxis()->SetTitle(
"Phi");
102 Map_X0->GetZaxis()->SetTitle(
"X0");
103 Map_L0->GetXaxis()->SetTitle(
"Eta");
104 Map_L0->GetYaxis()->SetTitle(
"Phi");
105 Map_L0->GetZaxis()->SetTitle(
"L0");
108 if(surface_info.
type == 2){
109 Map_X0 =
new TH2F((
"Map_X0_"+surface_info.
idname).c_str(),(
"Map_X0_"+surface_info.
idname).c_str(),
111 Map_L0 =
new TH2F((
"Map_L0_"+surface_info.
idname).c_str(),(
"Map_L0_"+surface_info.
idname).c_str(),
113 Map_scale =
new TH2F((
"Map_scale_"+surface_info.
idname).c_str(),(
"Map_scale_"+surface_info.
idname).c_str(),
115 Map_X0->GetXaxis()->SetTitle(
"X [mm]");
116 Map_X0->GetYaxis()->SetTitle(
"Y [mm]");
117 Map_X0->GetZaxis()->SetTitle(
"X0");
118 Map_L0->GetXaxis()->SetTitle(
"X [mm]");
119 Map_L0->GetYaxis()->SetTitle(
"Y [mm]");
120 Map_L0->GetZaxis()->SetTitle(
"L0");
122 std::vector<TH2F*> v_hist;
123 v_hist.push_back(Map_X0);
124 v_hist.push_back(Map_L0);
125 v_hist.push_back(Map_scale);
126 surface_hist = v_hist;
131 void Fill(std::map<uint64_t,std::vector<TH2F*>>& surface_hist, std::map<uint64_t,sinfo>& surface_info,
132 const std::string& input_file,
const std::string& json_surface_file,
const int& nbprocess){
135 std::map<std::string,std::string> surface_name;
137 if(json_surface_file !=
""){
138 std::ifstream lfile(json_surface_file.c_str());
144 std::map<uint64_t,float> surface_weight;
147 TFile *tfile =
new TFile(input_file.c_str());
148 TTree *
tree = (TTree*)tfile->Get(
"material-tracks");
152 std::vector<float> *mat_X0 = 0;
153 std::vector<float> *mat_L0 = 0;
154 std::vector<float> *mat_step_length = 0;
156 std::vector<uint64_t> *sur_id = 0;
157 std::vector<int32_t> *sur_type = 0;
158 std::vector<float> *sur_x = 0;
159 std::vector<float> *sur_y = 0;
160 std::vector<float> *sur_z = 0;
161 std::vector<float> *sur_range_min = 0;
162 std::vector<float> *sur_range_max = 0;
164 tree->SetBranchAddress(
"v_phi",&v_phi);
165 tree->SetBranchAddress(
"v_eta",&v_eta);
166 tree->SetBranchAddress(
"mat_X0",&mat_X0);
167 tree->SetBranchAddress(
"mat_L0",&mat_L0);
168 tree->SetBranchAddress(
"mat_step_length",&mat_step_length);
170 tree->SetBranchAddress(
"sur_id",&sur_id);
171 tree->SetBranchAddress(
"sur_type",&sur_type);
172 tree->SetBranchAddress(
"sur_x",&sur_x);
173 tree->SetBranchAddress(
"sur_y",&sur_y);
174 tree->SetBranchAddress(
"sur_z",&sur_z);
175 tree->SetBranchAddress(
"sur_range_min",&sur_range_min);
176 tree->SetBranchAddress(
"sur_range_max",&sur_range_max);
179 if(nentries > nbprocess && nbprocess != -1) nentries = nbprocess;
181 for (Long64_t i=0;i<
nentries; i++) {
182 if(i%10000==0) std::cout <<
"processed " << i <<
" events out of " << nentries << std::endl;
186 for (
auto weight_it = surface_weight.begin(); weight_it != surface_weight.end(); weight_it++){
187 weight_it->second = 0;
190 for(
int j=0; j<mat_X0->size(); j++ ){
193 if(sur_type->at(j) == -1)
continue;
195 if(surface_hist.find(sur_id->at(j))==surface_hist.end()){
197 if(sur_type->at(j) == 1){
198 pos = sqrt(sur_x->at(j)*sur_x->at(j)+sur_y->at(j)*sur_y->at(j));
200 if(sur_type->at(j) == 2){
203 surface_weight[sur_id->at(j)] = 0;
204 Initialise_info(surface_info[sur_id->at(j)], surface_name, sur_id->at(j), sur_type->at(j),
pos, sur_range_min->at(j), sur_range_max->at(j));
205 Initialise_hist(surface_hist[sur_id->at(j)], surface_info[sur_id->at(j)]);
208 surface_weight[sur_id->at(j)]++;
212 for(
int j=0; j<mat_X0->size(); j++ ){
215 if(sur_type->at(j) == -1)
continue;
217 if(sur_type->at(j) == 1){
218 surface_hist[sur_id->at(j)][0]->Fill(v_eta, v_phi, (mat_step_length->at(j)/mat_X0->at(j)));
219 surface_hist[sur_id->at(j)][1]->Fill(v_eta, v_phi, (mat_step_length->at(j)/mat_L0->at(j)));
220 surface_hist[sur_id->at(j)][2]->Fill(v_eta, v_phi, (1/surface_weight[sur_id->at(j)]));
222 if(sur_type->at(j) == 2){
223 surface_hist[sur_id->at(j)][0]->Fill(sur_x->at(j), sur_y->at(j), (mat_step_length->at(j)/mat_X0->at(j)));
224 surface_hist[sur_id->at(j)][1]->Fill(sur_x->at(j), sur_y->at(j), (mat_step_length->at(j)/mat_L0->at(j)));
225 surface_hist[sur_id->at(j)][2]->Fill(sur_x->at(j), sur_y->at(j), (1/surface_weight[sur_id->at(j)]));
230 for (
auto hist_it = surface_hist.begin(); hist_it != surface_hist.end(); hist_it++){
231 hist_it->second[0]->Divide(hist_it->second[2]);
232 hist_it->second[1]->Divide(hist_it->second[2]);
243 void Mat_map_surface_plot(std::string input_file =
"", std::string json_surface_file =
"",
int nbprocess = -1, std::string
name =
""){
245 gStyle->SetOptStat(0);
246 gStyle->SetOptTitle(0);
248 std::map<uint64_t,std::vector<TH2F*>> surface_hist;
249 std::map<uint64_t,sinfo> surface_info;
251 Fill(surface_hist, surface_info, input_file, json_surface_file, nbprocess);
252 for (
auto hist_it = surface_hist.begin(); hist_it != surface_hist.end(); hist_it++){
253 plot(hist_it->second, surface_info[hist_it->first],
name);
254 for (
auto hist : hist_it->second){
257 hist_it->second.clear();