34 : m_detector(
Detector::getInstance())
41 : m_detector(description)
64 std::set<const TGeoNode*>& cont;
65 PvCollector(
Region r, std::set<const TGeoNode*>& c) : rg(r), cont(c) {}
66 void collect(TGeoNode* pv) {
68 for (Int_t idau = 0, ndau = pv->GetNdaughters(); idau < ndau; ++idau)
69 collect(pv->GetDaughter(idau));
71 void operator()(TGeoNode* pv) {
78 for (Int_t idau = 0, ndau = pv->GetNdaughters(); idau < ndau; ++idau)
79 (*this)(pv->GetDaughter(idau));
86 printout(ALWAYS,
"MaterialScan",
"+++ Set new scanning region to: %s [%ld placements]",
90 printout(ALWAYS,
"MaterialScan",
"+++ No region restrictions set. Back to full scanning mode.");
106 std::set<const TGeoNode*>& cont;
107 PvCollector(std::set<const TGeoNode*>& c) : cont(c) {}
108 void operator()(TGeoNode* pv) {
110 for (Int_t idau = 0, ndau = pv->GetNdaughters(); idau < ndau; ++idau) {
111 TGeoNode* daughter = pv->GetDaughter(idau);
123 printout(ALWAYS,
"MaterialScan",
"+++ Set new scanning volume to: %s [%ld placements]",
127 printout(ALWAYS,
"MaterialScan",
"+++ No subdetector restrictions set. Back to full scanning mode.");
145 std::set<const TGeoNode*>& cont;
146 PvCollector(
Material mat, std::set<const TGeoNode*>& c) : material(mat), cont(c) {}
147 void operator()(TGeoNode* pv) {
148 Volume vol = pv->GetVolume();
150 if ( mat.
ptr() == material.
ptr() ) {
153 for (Int_t idau = 0, ndau = pv->GetNdaughters(); idau < ndau; ++idau)
154 (*this)(pv->GetDaughter(idau));
161 printout(ALWAYS,
"MaterialScan",
"+++ Set new scanning material to: %s [%ld placements]",
165 printout(ALWAYS,
"MaterialScan",
"+++ No material restrictions set. Back to full scanning mode.");
171 Vector3D p0(x0, y0, z0), p1(x1, y1, z1);
181 double sum_lambda = 0;
182 double path_length = 0, total_length = 0;
183 const char* fmt1 =
" | %7s %-25s %3.0f %8.3f %8.4f %11.4f %11.4f %10.3f %8.2f %11.6f %11.6f (%7.2f,%7.2f,%7.2f)\n";
184 const char* fmt2 =
" | %7s %-25s %3.0f %8.3f %8.4f %11.6g %11.6g %10.3f %8.2f %11.6f %11.6f (%7.2f,%7.2f,%7.2f)\n";
185 const char* line =
" +--------------------------------------------------------------------------------------------------------------------------------------------------\n";
187 direction = (p1-p0).unit();
189 ::printf(
"%s + Material scan between: x_0 = (%7.2f,%7.2f,%7.2f) [cm] and x_1 = (%7.2f,%7.2f,%7.2f) [cm] : \n%s",
190 line,p0[0],p0[1],p0[2],p1[0],p1[1],p1[2],line);
191 ::printf(
" | \\ %-16s Atomic Radiation Interaction Path Integrated Integrated Material\n",
"Material");
192 ::printf(
" | Num. \\ %-16s Number/Z Mass/A Density Length Length Thickness Length X0 Lambda Endpoint/Startpoint\n",
"Name");
193 ::printf(
" | Layer \\ %-16s [g/mole] [g/cm3] [cm] [cm] [cm] [cm] [cm] [cm] ( cm, cm, cm)\n",
"");
196 for(
unsigned i=0, n=placements.size(); i<n; ++i){
197 TGeoNode* pv = placements[i].first.ptr();
198 double length = placements[i].second;
199 total_length += length;
200 end = p0 + total_length * direction;
203 ::printf(
"%p %s %s %s\n",
204 placements[i].first.ptr(),
205 placements[i].first->GetName(),
206 placements[i].first->GetVolume()->GetName(),
207 placements[i].first->GetMedium()->GetName());
211 TGeoMaterial* curr_mat = placements[i].first->GetMedium()->GetMaterial();
212 TGeoMaterial* next_mat = (i+1)<n ? placements[i+1].first->GetMedium()->GetMaterial() :
nullptr;
213 double nx0 = length / curr_mat->GetRadLen();
214 double nLambda = length / curr_mat->GetIntLen();
217 sum_lambda += nLambda;
218 path_length += length;
219 materials.emplace_back(placements[i].first->GetMedium(),length);
220 const char* fmt = curr_mat->GetRadLen() >= 1e5 ? fmt2 : fmt1;
221 std::string mname = curr_mat->GetName();
222 if ( next_mat && (i+1)<n ) {
224 mname += next_mat->GetName();
227 ::printf(fmt,
"(start)" , curr_mat->GetName(), curr_mat->GetZ(), curr_mat->GetA(),
228 curr_mat->GetDensity(), curr_mat->GetRadLen()/dd4hep::cm, curr_mat->GetIntLen()/dd4hep::cm,
230 p0[0]/dd4hep::cm, p0[1]/dd4hep::cm, p0[2]/dd4hep::cm);
234 ::printf(fmt,
"(end)", curr_mat->GetName(), curr_mat->GetZ(), curr_mat->GetA(),
235 curr_mat->GetDensity(), curr_mat->GetRadLen()/dd4hep::cm, curr_mat->GetIntLen()/dd4hep::cm,
237 p0[0]/dd4hep::cm, p0[1]/dd4hep::cm, p0[2]/dd4hep::cm);
240 ::printf(fmt, std::to_string(i+1).c_str(), mname.c_str(), next_mat->GetZ(), next_mat->GetA(),
241 next_mat->GetDensity(), next_mat->GetRadLen()/dd4hep::cm, next_mat->GetIntLen()/dd4hep::cm,
242 length/dd4hep::cm, path_length/dd4hep::cm, sum_x0, sum_lambda,
243 end[0]/dd4hep::cm, end[1]/dd4hep::cm, end[2]/dd4hep::cm);
247 const MaterialData& avg = matMgr.createAveragedMaterial(materials);
249 ::printf(fmt,
"",
"Average Material",avg.
Z(),avg.
A(),avg.
density(),
251 path_length/dd4hep::cm, path_length/dd4hep::cm,
254 end[0]/dd4hep::cm, end[1]/dd4hep::cm, end[2]/dd4hep::cm);
260 Vector3D p0(x0, y0, z0), p1(x1, y1, z1);