ATTPCROOT  0.3.0-alpha
A ROOT-based framework for analyzing data from active target detectors
AtTpcMap.cxx
Go to the documentation of this file.
1 /*********************************************************************
2  * ATTPC Mapping Class AtTpcMap.cxx *
3  * Author: Y. Ayyad *
4  * Log: 13-02-2015 17:16 JST *
5  * *
6  *********************************************************************/
7 
8 #include "AtTpcMap.h"
9 
10 #include <FairLogger.h>
11 
12 #include <Math/Point2D.h>
13 #include <Rtypes.h>
14 #include <TH2Poly.h>
15 #include <TMath.h>
16 #include <TMathBase.h>
17 
18 #include <boost/multi_array/base.hpp>
19 #include <boost/multi_array/extent_gen.hpp>
20 #include <boost/multi_array/multi_array_ref.hpp>
21 #include <boost/multi_array/subarray.hpp>
22 
23 #include <algorithm>
24 #include <cmath>
25 #include <fstream> // IWYU pragma: keep
26 #include <iostream>
27 #include <vector>
28 
29 #undef BOOST_MULTI_ARRAY_NO_GENERATORS
30 #define BOOST_MULTI_ARRAY_NO_GENERATORS
31 
32 using std::cout;
33 using std::endl;
34 using XYPoint = ROOT::Math::XYPoint;
35 
36 constexpr auto cRED = "\033[1;31m";
37 constexpr auto cYELLOW = "\033[1;33m";
38 constexpr auto cNORMAL = "\033[0m";
39 
40 AtTpcMap::AtTpcMap() : AtMap()
41 {
42  AtPadCoord.resize(boost::extents[10240][3][2]);
43  std::fill(AtPadCoord.data(), AtPadCoord.data() + AtPadCoord.num_elements(), 0);
44  std::cout << " ATTPC Map initialized " << std::endl;
45  std::cout << " ATTPC Pad Coordinates container initialized " << std::endl;
46  fNumberPads = 10240;
47 }
48 
49 AtTpcMap::~AtTpcMap() = default;
50 
51 void AtTpcMap::Dump()
52 {
53 
54  std::ofstream coordmap;
55  coordmap.open("coordmap.txt");
56 
57  int values = 0;
58  for (index i = 0; i != 10240; ++i) {
59  coordmap << i << " ";
60  for (index j = 0; j != 3; ++j) {
61  for (index k = 0; k != 2; ++k) {
62  std::cout << " ATTPC Triangular pad coordinates - Pad Index : " << i << " X(" << j << ") - Y(" << k
63  << ") :" << AtPadCoord[i][j][k] << std::endl;
64  coordmap << AtPadCoord[i][j][k] << " ";
65  } // k
66  } // j
67 
68  coordmap << std::endl;
69 
70  } // i
71 
72  coordmap.close();
73 }
74 
75 void AtTpcMap::GeneratePadPlane()
76 {
77  if (fPadPlane) {
78  LOG(error) << "Skipping generation of pad plane because it was already parsed!";
79  return;
80  }
81 
82  std::cout << " ATTPC Map : Generating the map geometry of the ATTPC " << std::endl;
83  // Local variables
84  Float_t pads_in_half_hex = 0;
85  Float_t pads_in_hex = 0;
86  Float_t row_length = 0;
87  Float_t pads_in_half_row = 0;
88  Int_t pads_out_half_hex = 0;
89  Int_t pads_in_row = 0;
90  Int_t ort = 0;
91  Float_t pad_x_off = 0;
92  Float_t pad_y_off = 0;
93  Float_t tmp_pad_x_off = 0;
94  Float_t tmp_pad_y_off = 0;
95 
96  Float_t small_z_spacing = 2 * 25.4 / 1000.;
97  Float_t small_tri_side = 184. * 25.4 / 1000.;
98  Double_t umega_radius = 10826.772 * 25.4 / 1000.;
99  Float_t beam_image_radius = 4842.52 * 25.4 / 1000.;
100  Int_t pad_index = 0;
101  Int_t pad_index_aux = 0;
102 
103  Float_t small_x_spacing = 2. * small_z_spacing / TMath::Sqrt(3.);
104  Float_t small_y_spacing = small_x_spacing * TMath::Sqrt(3.);
105  Float_t dotted_s_tri_side = 4. * small_x_spacing + small_tri_side;
106  Float_t dotted_s_tri_hi = dotted_s_tri_side * TMath::Sqrt(3.) / 2.;
107  Float_t dotted_l_tri_side = 2. * dotted_s_tri_side;
108  Float_t dotted_l_tri_hi = dotted_l_tri_side * TMath::Sqrt(3.) / 2.;
109  Float_t large_x_spacing = small_x_spacing;
110  Float_t large_y_spacing = small_y_spacing;
111  Float_t large_tri_side = dotted_l_tri_side - 4. * large_x_spacing;
112  Float_t large_tri_hi = dotted_l_tri_side * TMath::Sqrt(3.) / 2.;
113  // Float_t row_len_s = 2**TMath::Ceil(TMath::Log(beam_image_radius/dotted_s_tri_side)/TMath::Log(2.0));
114  Float_t row_len_s = pow(2, TMath::Ceil(TMath::Log(beam_image_radius / dotted_s_tri_side) / TMath::Log(2.0)));
115  Float_t row_len_l = TMath::Floor(umega_radius / dotted_l_tri_hi);
116 
117  Float_t xoff = 0.;
118  Float_t yoff = 0.;
119 
120  for (Int_t j = 0; j < row_len_l; j++) {
121  pads_in_half_hex = 0;
122  pads_in_hex = 0;
123  // row_length = TMath::Abs(sqrt(umega_radius**2 - (j*dotted_l_tri_hi + dotted_l_tri_hi/2.)**2));
124  row_length = TMath::Abs(sqrt(pow(umega_radius, 2) - pow((j * dotted_l_tri_hi + dotted_l_tri_hi / 2.), 2)));
125 
126  if (j < row_len_s / 2.) {
127 
128  pads_in_half_hex = (2 * row_len_s - 2 * j) / 4.;
129  pads_in_hex = 2 * row_len_s - 1. - 2. * j;
130 
131  } // if row_len_s
132 
133  pads_in_half_row = row_length / dotted_l_tri_side;
134  pads_out_half_hex = static_cast<Int_t>(TMath::Nint(2 * (pads_in_half_row - pads_in_half_hex)));
135  pads_in_row = 2 * pads_out_half_hex + 4 * pads_in_half_hex - 1;
136 
137  ort = 1;
138 
139  for (Int_t i = 0; i < pads_in_row; i++) {
140 
141  // set initial ort
142  if (i == 0) {
143  if (j % 2 == 0)
144  ort = -1;
145  if (((pads_in_row - 1) / 2) % 2 == 1)
146  ort = -ort;
147 
148  } // i==0
149  else
150  ort = -ort;
151 
152  pad_x_off = -(pads_in_half_hex + pads_out_half_hex / 2.) * dotted_l_tri_side + i * dotted_l_tri_side / 2. +
153  2. * large_x_spacing + xoff;
154 
155  if (i < pads_out_half_hex || i > (pads_in_hex + pads_out_half_hex - 1) || j > (row_len_s / 2. - 1)) {
156 
157  pad_y_off = j * dotted_l_tri_hi + large_y_spacing + yoff;
158  if (ort == -1)
159  pad_y_off += large_tri_hi;
160 
161  fill_coord(pad_index, pad_x_off, pad_y_off, large_tri_side, ort);
162  pad_index += 1;
163 
164  } // if
165  else {
166 
167  pad_y_off = j * dotted_l_tri_hi + large_y_spacing + yoff;
168  if (ort == -1)
169  pad_y_off = j * dotted_l_tri_hi + 2 * dotted_s_tri_hi - small_y_spacing + yoff;
170  fill_coord(pad_index, pad_x_off, pad_y_off, small_tri_side, ort);
171 
172  pad_index += 1;
173  tmp_pad_x_off = pad_x_off + dotted_s_tri_side / 2.;
174  tmp_pad_y_off = pad_y_off + ort * dotted_s_tri_hi - 2 * ort * small_y_spacing;
175  fill_coord(pad_index, tmp_pad_x_off, tmp_pad_y_off, small_tri_side, -ort);
176 
177  pad_index += 1;
178  tmp_pad_y_off = pad_y_off + ort * dotted_s_tri_hi;
179  fill_coord(pad_index, tmp_pad_x_off, tmp_pad_y_off, small_tri_side, ort);
180 
181  pad_index += 1;
182  tmp_pad_x_off = pad_x_off + dotted_s_tri_side;
183  fill_coord(pad_index, tmp_pad_x_off, pad_y_off, small_tri_side, ort);
184  pad_index += 1;
185  }
186 
187  } // pads_in_row loop
188 
189  } // row_len_l
190 
191  for (Int_t i = 0; i < pad_index; i++) {
192  for (Int_t j = 0; j < 3; j++) {
193  AtPadCoord[i + pad_index][j][0] = AtPadCoord[i][j][0];
194  AtPadCoord[i + pad_index][j][1] = -AtPadCoord[i][j][1];
195  }
196  pad_index_aux++;
197  }
198 
199  // fPadInd = pad_index + pad_index_aux;
200  std::cout << "created pads: " << pad_index + pad_index_aux << std::endl;
201  kIsParsed = true;
202 
203  if (fPadPlane != nullptr)
204  delete fPadPlane;
205  fPadPlane = new TH2Poly(); // NOLINT
206 
207  fPadPlane->SetName("ATTPC_Plane");
208  fPadPlane->SetTitle("ATTPC_Plane");
209 
210  // for (Int_t i = 0; i < fPadInd; i++) {
211  for (Int_t i = 0; i < fNumberPads; i++) {
212 
213  Double_t px[] = {AtPadCoord[i][0][0], AtPadCoord[i][1][0], AtPadCoord[i][2][0], AtPadCoord[i][0][0]};
214  Double_t py[] = {AtPadCoord[i][0][1], AtPadCoord[i][1][1], AtPadCoord[i][2][1], AtPadCoord[i][0][1]};
215  fPadPlane->AddBin(4, px, py);
216  }
217 
218  fPadPlane->ChangePartition(500, 500);
219 }
220 
221 Int_t AtTpcMap::fill_coord(int pindex, float padxoff, float padyoff, float triside, float fort)
222 {
223  AtPadCoord[pindex][0][0] = padxoff;
224  AtPadCoord[pindex][0][1] = padyoff;
225  AtPadCoord[pindex][1][0] = padxoff + triside / 2.;
226  AtPadCoord[pindex][1][1] = padyoff + fort * triside * TMath::Sqrt(3.) / 2.;
227  AtPadCoord[pindex][2][0] = padxoff + triside;
228  AtPadCoord[pindex][2][1] = padyoff;
229  return 0;
230 }
231 
232 XYPoint AtTpcMap::CalcPadCenter(Int_t PadRef)
233 {
234  if (!kIsParsed) {
235  LOG(error) << " AtTpcMap::CalcPadCenter Error : Pad plane has not been generated or parsed ";
236  return {-9999, -9999};
237  }
238  if (PadRef < 0) { // Boost multi_array crashes with a negative index
239  LOG(debug) << " AtTpcMap::CalcPadCenter Error : Pad not found";
240  return {-9999, -9999};
241  }
242 
243  Float_t x = (AtPadCoord[PadRef][0][0] + AtPadCoord[PadRef][1][0] + AtPadCoord[PadRef][2][0]) / 3.;
244  Float_t y = (AtPadCoord[PadRef][0][1] + AtPadCoord[PadRef][1][1] + AtPadCoord[PadRef][2][1]) / 3.;
245  return {x, y};
246 }
247 ClassImp(AtTpcMap)
AtMap
Definition: AtMap.h:33
cNORMAL
constexpr auto cNORMAL
Definition: AtTpcMap.cxx:38
ClassImp
ClassImp(AtFindVertex)
AtTpcMap::GeneratePadPlane
virtual void GeneratePadPlane() override
Definition: AtTpcMap.cxx:75
cRED
constexpr auto cRED
Definition: AtTpcMap.cxx:36
XYPoint
ROOT::Math::XYPoint XYPoint
Definition: AtPatternCircle2D.cxx:16
cYELLOW
constexpr auto cYELLOW
Definition: AtTpcMap.cxx:37
AtTpcMap::Dump
virtual void Dump() override
Definition: AtTpcMap.cxx:51
XYPoint
ROOT::Math::XYPoint XYPoint
Definition: AtTpcMap.cxx:34
AtMap::AtPadCoord
multiarray AtPadCoord
Definition: AtMap.h:41
AtMap::fNumberPads
UInt_t fNumberPads
Definition: AtMap.h:48
AtMap::kIsParsed
Bool_t kIsParsed
Definition: AtMap.h:43
AtTpcMap
Definition: AtTpcMap.h:19
AtTpcMap::fill_coord
Int_t fill_coord(int pindex, float padxoff, float padyoff, float triside, float fort)
Definition: AtTpcMap.cxx:221
y
const double * y
Definition: lmcurve.cxx:20
AtTpcMap::~AtTpcMap
~AtTpcMap()
AtMap::index
multiarray::index index
Definition: AtMap.h:39
AtTpcMap::AtTpcMap
AtTpcMap()
Definition: AtTpcMap.cxx:40
AtTpcMap.h
AtMap::fPadPlane
TH2Poly * fPadPlane
Definition: AtMap.h:47
AtTpcMap::CalcPadCenter
virtual ROOT::Math::XYPoint CalcPadCenter(Int_t PadRef) override
Definition: AtTpcMap.cxx:232
Generated on Sat Nov 11 2023 15:30:38 for ATTPCROOT by   doxygen 1.8.18