350{
351 using namespace Acts::VectorHelpers;
352
353
354 double rmin_at_center = std::numeric_limits<double>::max();
355 double rmin_at_choke = std::numeric_limits<double>::max();
356 double rmax = std::numeric_limits<double>::lowest();
357 double zmin = std::numeric_limits<double>::max();
358 double zmax = std::numeric_limits<double>::lowest();
359 double cutout_zmin_abs = std::numeric_limits<double>::max();
360
361
362
363
364
365 for (const auto& box : boxStore) {
366 Acts::Vector3
vmin = box->min().cast<
double>();
367 Acts::Vector3
vmax = box->max().cast<
double>();
368
369 double vrmin =
perp(vmin);
370 double vrmax =
perp(vmax);
371
372 rmin_at_choke = std::min(rmin_at_choke, std::min(vrmin, vrmax));
373
374 rmax = std::max(rmax, std::max(vrmin, vrmax));
377
378 if (std::abs(
vmin.z()) < 100) {
379 rmin_at_center = std::min(vrmin, rmin_at_center);
380 }
381 if (std::abs(
vmax.z()) < 100) {
382 rmin_at_center = std::min(vrmax, rmin_at_center);
383 }
384 }
385
386 for (const auto& box : boxStore) {
387 Acts::Vector3
vmin = box->min().cast<
double>();
388 Acts::Vector3
vmax = box->max().cast<
double>();
389 double vrmin =
perp(vmin);
390 double vrmax =
perp(vmax);
391
392 if (vrmin < rmin_at_center * 0.9) {
393 cutout_zmin_abs = std::min(cutout_zmin_abs, std::abs(
vmin.z()));
394 }
395 if (vrmax < rmin_at_center * 0.9) {
396 cutout_zmin_abs = std::min(cutout_zmin_abs, std::abs(
vmax.z()));
397 }
398 }
399
401 double dz2 = cutout_zmin_abs;
402
403
404 double envZ = 5;
405 double envR = 5;
406 dz1 += envZ;
407 dz2 -= envZ;
408 rmax += envR;
409 if(rmin_at_choke > envR) rmin_at_choke -= envR;
410 rmin_at_center -= envR;
411
412
414 << " rmin at choke: " << rmin_at_choke
415 << " rmax: " << rmax << " zmin: " << zmin << " zmax: " << zmax
416 << " coutout_zmin_abs: " << cutout_zmin_abs);
417
418
419
420
421
422
423 auto idCylBds
424 = dynamic_cast<const Acts::CylinderVolumeBounds*>(&insideVolume->volumeBounds());
425 if (idCylBds == nullptr) {
426 ATH_MSG_ERROR(
"Unable to dynamic cast volume bounds to Acts::CylinderVolumeBounds");
427 throw std::runtime_error("Error casting to CylinderVolumeBounds");
428 }
429
430 double idRMax = idCylBds->get(CVBBV::eMaxR);
431 double idRMin = idCylBds->get(CVBBV::eMinR);
432 double idHlZ = idCylBds->get(CVBBV::eHalfLengthZ);
433
435
437
438 const auto&
trf = insideVolume->localToGlobalTransform(gctx);
439
440 if (!
trf.isApprox(Acts::Transform3::Identity())) {
442
443
444
445
446 Acts::RotationMatrix3 rot =
trf.rotation();
447 bool unityRot = rot.isApprox(Acts::RotationMatrix3::Identity());
448
450
451
452 const Acts::Vector3 trl =
trf.translation();
453 bool transZOnly = std::abs(1 - std::abs(Acts::Vector3::UnitZ().
dot(trl.normalized()))) < 1
e-6;
454
456 ATH_MSG_VERBOSE(
"TRL "<< trl.normalized().dot(Acts::Vector3::UnitZ()));
457
458 if(!unityRot || !transZOnly) {
459 ATH_MSG_ERROR(
"The ID appears to be shifted from the origin. I cannot handle this.");
461 throw std::runtime_error("Error building calo");
462 }
463 else {
464 ATH_MSG_VERBOSE(
"Checked: non unitarity is ONLY due to shift along z axis: that's ok");
465 double prevIdHlZ = idHlZ;
466 idHlZ += std::abs(trl.z());
467 ATH_MSG_VERBOSE(
"Modifying effective half length of ID cylinder: " << prevIdHlZ <<
" => " << idHlZ);
468 }
469 }
470
471
472 if (idRMax > rmin_at_center || idHlZ > dz2 || (idRMin > rmin_at_choke && idRMin != 0.)) {
474 ATH_MSG_ERROR(
"This can be because the ID overlaps into the calo volume");
475 ATH_MSG_ERROR(
"Or because the Calo choke radius is SMALLER than the ID inner radius");
476 ATH_MSG_ERROR(
"Currently, I can only make the choke radius smaller, I can not make it larger");
477 ATH_MSG_ERROR(
"nor can I manipulate the ID volume bounds at this point.");
478 ATH_MSG_ERROR(
"ID rMax: " << idRMax <<
" Calo rMin@center: " << rmin_at_center);
479 ATH_MSG_ERROR(
"ID hlZ: " << idHlZ <<
" Calo inner Z hl: " << dz2);
480 ATH_MSG_ERROR(
"ID rMin: " << idRMin <<
" Calo rMin@choke: " << rmin_at_choke);
482 throw std::runtime_error("Error building calo");
483 }
484
485
486
487
488 rmin_at_choke = idRMin;
489
490 std::shared_ptr<Acts::CutoutCylinderVolumeBounds> volBds = nullptr;
491 volBds = std::make_shared<Acts::CutoutCylinderVolumeBounds>(
492 rmin_at_choke, rmin_at_center, rmax, dz1, dz2);
493
495
496 return volBds;
497}
Scalar perp() const
perp method - perpendicular length
#define ATH_MSG_VERBOSE(x)
std::string toString(const Translation3D &translation, int precision=4)
GeoPrimitvesToStringConverter.
dot(G, fn, nodesToHighlight=[])