diff --git a/.github/workflows/linux-eic-shell.yml b/.github/workflows/linux-eic-shell.yml index 7375e241c..f108b7de3 100644 --- a/.github/workflows/linux-eic-shell.yml +++ b/.github/workflows/linux-eic-shell.yml @@ -103,7 +103,7 @@ jobs: if-no-files-found: error - run: | source install/setup.sh - sed -i '/\1|g; s|\(/fiber>\)|\1|g' sed -i '/- + If you use this software, please cite it using the + metadata from this file. +type: software +authors: + - given-names: EPIC Collaboration +repository-code: 'https://github.com/eic/epic' +abstract: DD4hep Geometry Description of the EPIC Experiment diff --git a/compact/ecal/barrel_interlayers.xml b/compact/ecal/barrel_interlayers.xml index a6472046f..cf28c16fb 100644 --- a/compact/ecal/barrel_interlayers.xml +++ b/compact/ecal/barrel_interlayers.xml @@ -27,7 +27,7 @@ - + @@ -39,13 +39,14 @@ + For Pb/SiFi (GlueX): X0 ~ 1.45 cm For W/SiFi (sPHENIX): X0 ~ 0.7 cm (but different fiber orientation) - - + + @@ -142,11 +143,7 @@ - - - - + @@ -154,6 +151,7 @@ @@ -170,6 +168,7 @@ diff --git a/src/BarrelCalorimeterInterlayers_geo.cpp b/src/BarrelCalorimeterInterlayers_geo.cpp index e6049d864..812fda6d5 100644 --- a/src/BarrelCalorimeterInterlayers_geo.cpp +++ b/src/BarrelCalorimeterInterlayers_geo.cpp @@ -179,6 +179,7 @@ void buildFibers(Detector& desc, SensitiveDetector& sens, Volume& s_vol, xml_com { auto [s_trd_x1, s_thick, s_length, hphi] = dimensions; double f_radius = getAttrOrDefault(x_fiber, _U(radius), 0.1 * cm); + double f_cladding_thickness = getAttrOrDefault(x_fiber, _Unicode(cladding_thickness), 0.0 * cm); double f_spacing_x = getAttrOrDefault(x_fiber, _Unicode(spacing_x), 0.122 * cm); double f_spacing_z = getAttrOrDefault(x_fiber, _Unicode(spacing_z), 0.134 * cm); std::string f_id_grid = getAttrOrDefault(x_fiber, _Unicode(identifier_grid), "grid"); @@ -195,8 +196,11 @@ void buildFibers(Detector& desc, SensitiveDetector& sens, Volume& s_vol, xml_com auto grid_div = getNdivisions(s_trd_x1, s_thick, 2.0 * cm, 2.0 * cm); // Calculate polygonal grid coordinates (vertices) auto grid_vtx = gridPoints(grid_div.first, grid_div.second, s_trd_x1, s_thick, hphi); - Tube f_tube(0, f_radius, s_length); - Volume f_vol("fiber_vol", f_tube, desc.material(x_fiber.materialStr())); + double f_radius_core = f_radius-f_cladding_thickness; + Tube f_tube_clad(f_radius_core, f_radius, s_length); + Volume f_vol_clad("fiber_vol", f_tube_clad, desc.material(x_fiber.materialStr())); + Tube f_tube_core(0, f_radius_core, s_length); + Volume f_vol_core("fiber_core_vol", f_tube_core, desc.material(x_fiber.materialStr())); vector f_id_count(grid_div.first * grid_div.second, 0); auto f_pos = fiberPositions(f_radius, f_spacing_x, f_spacing_z, s_trd_x1, s_thick, hphi); @@ -208,7 +212,8 @@ void buildFibers(Detector& desc, SensitiveDetector& sens, Volume& s_vol, xml_com } double l_pos_y = line.front().y(); // use assembly as intermediate volume container to reduce number of daughter volumes - Assembly lfibers(Form("fiber_array_line_%lu", il)); + Assembly lfibers_clad(Form("fiber_clad_array_line_%lu", il)); + Assembly lfibers_core(Form("fiber_core_array_line_%lu", il)); for (auto& p : line) { int f_grid_id = -1; int f_id = -1; @@ -236,19 +241,22 @@ void buildFibers(Detector& desc, SensitiveDetector& sens, Volume& s_vol, xml_com } if (x_fiber.isSensitive()) { - f_vol.setSensitiveDetector(sens); + f_vol_core.setSensitiveDetector(sens); } - f_vol.setAttributes(desc, x_fiber.regionStr(), x_fiber.limitsStr(), x_fiber.visStr()); + f_vol_core.setAttributes(desc, x_fiber.regionStr(), x_fiber.limitsStr(), x_fiber.visStr()); // Fiber placement // Transform3D f_tr(RotationZYX(0,0,M_PI*0.5),Position(p.x(), 0, p.y())); // PlacedVolume fiber_phv = s_vol.placeVolume(f_vol, Position(p.x(), 0., p.y())); - PlacedVolume fiber_phv = lfibers.placeVolume(f_vol, Position(p.x(), 0., 0.)); - fiber_phv.addPhysVolID(f_id_grid, f_grid_id + 1).addPhysVolID(f_id_fiber, f_id + 1); + PlacedVolume core_phv = lfibers_core.placeVolume(f_vol_core, Position(p.x(), 0., 0.)); + core_phv.addPhysVolID(f_id_grid, f_grid_id + 1).addPhysVolID(f_id_fiber, f_id + 1); + lfibers_clad.placeVolume(f_vol_clad, Position(p.x(), 0., 0.)); } - lfibers.ptr()->Voxelize(""); + lfibers_core.ptr()->Voxelize(""); + lfibers_clad.ptr()->Voxelize(""); Transform3D l_tr(RotationZYX(0, 0, M_PI * 0.5), Position(0., 0, l_pos_y)); - s_vol.placeVolume(lfibers, l_tr); + s_vol.placeVolume(lfibers_core, l_tr); + s_vol.placeVolume(lfibers_clad, l_tr); } } diff --git a/src/DRICH_geo.cpp b/src/DRICH_geo.cpp index 496828379..9cab8cf60 100644 --- a/src/DRICH_geo.cpp +++ b/src/DRICH_geo.cpp @@ -134,10 +134,10 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec default: printout(FATAL, "DRICH_geo", "UNKNOWN debugOpticsMode"); return det; - }; + } aerogelVis = sensorVis = mirrorVis; gasvolVis = vesselVis = desc.invisible(); - }; + } // sensor size rescaling for large sensor (catching all photons) double sensorRescale = 0; @@ -234,7 +234,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec vesselSolid = vesselBox; gasvolSolid = gasvolUnion; break; - }; + } // volumes Volume vesselVol(detName, vesselSolid, vesselMat); @@ -298,42 +298,37 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec // SkinSurface aerogelSkin(desc, aerogelDE, "mirror_optical_surface", aerogelSurf, aerogelVol); // aerogelSkin.isValid(); - // airgap placement and surface properties - PlacedVolume airgapPV; + // airgap and filter placement and surface properties if (!debugOptics) { + auto airgapPlacement = Translation3D(radiatorPos.x(), radiatorPos.y(), radiatorPos.z()) * // re-center to originFront RotationY(radiatorPitch) * // change polar angle Translation3D(0., 0., (aerogelThickness + airgapThickness) / 2.); // move to aerogel backplane - airgapPV = gasvolVol.placeVolume(airgapVol, airgapPlacement); + auto airgapPV = gasvolVol.placeVolume(airgapVol, airgapPlacement); DetElement airgapDE(det, "airgap_de", 0); airgapDE.setPlacement(airgapPV); - } - // filter placement and surface properties - PlacedVolume filterPV; - if (!debugOptics) { auto filterPlacement = Translation3D(0., 0., airgapThickness) * // add an air gap Translation3D(radiatorPos.x(), radiatorPos.y(), radiatorPos.z()) * // re-center to originFront RotationY(radiatorPitch) * // change polar angle Translation3D(0., 0., (aerogelThickness + filterThickness) / 2.); // move to aerogel backplane - filterPV = gasvolVol.placeVolume(filterVol, filterPlacement); + auto filterPV = gasvolVol.placeVolume(filterVol, filterPlacement); DetElement filterDE(det, "filter_de", 0); filterDE.setPlacement(filterPV); // SkinSurface filterSkin(desc, filterDE, "mirror_optical_surface", filterSurf, filterVol); // filterSkin.isValid(); - }; - // radiator z-positions (w.r.t. IP) - if (!debugOptics) { - auto airgapZpos = vesselPos.z() + airgapPV.position().z(); + // radiator z-positions (w.r.t. IP); only needed downstream if !debugOptics + double aerogelZpos = vesselPos.z() + aerogelPV.position().z(); + double airgapZpos = vesselPos.z() + airgapPV.position().z(); + double filterZpos = vesselPos.z() + filterPV.position().z(); + desc.add(Constant("DRICH_aerogel_zpos", std::to_string(aerogelZpos))); desc.add(Constant("DRICH_airgap_zpos", std::to_string(airgapZpos))); - auto filterZpos = vesselPos.z() + filterPV.position().z(); desc.add(Constant("DRICH_filter_zpos", std::to_string(filterZpos))); - }; - auto aerogelZpos = vesselPos.z() + aerogelPV.position().z(); - desc.add(Constant("DRICH_aerogel_zpos", std::to_string(aerogelZpos))); + } + // radiator material names desc.add(Constant("DRICH_aerogel_material", aerogelMat.ptr()->GetName(), "string")); desc.add(Constant("DRICH_airgap_material", airgapMat.ptr()->GetName(), "string")); @@ -427,7 +422,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec if (debugMirror) { mirrorTheta1 = 0; mirrorTheta2 = M_PI; /*mirrorPhiw=2*M_PI;*/ - }; + } // solid : create sphere at origin, with specified angular limits; // phi limits are increased to fill gaps (overlaps are cut away later) @@ -472,7 +467,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec // if debugging sphere properties, restrict number of sensors drawn if (debugSensors) { sensorSide = 2 * M_PI * sensorSphRadius / 64; - }; + } // solid and volume: single sensor module Box sensorSolid(sensorSide / 2., sensorSide / 2., sensorThickness / 2.); @@ -518,66 +513,66 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec // thetaGen loop: iterate less than "0.5 circumference / sensor size" times double nTheta = M_PI * sensorSphRadius / (sensorSide + sensorGap); - if(debugOpticsMode != 4){ + if(debugOpticsMode != 4) { for (int t = 0; t < (int)(nTheta + 0.5); t++) { - double thetaGen = t / ((double)nTheta) * M_PI; - - // phiGen loop: iterate less than "circumference at this latitude / sensor size" times - double nPhi = 2 * M_PI * sensorSphRadius * std::sin(thetaGen) / (sensorSide + sensorGap); - for (int p = 0; p < (int)(nPhi + 0.5); p++) { - double phiGen = p / ((double)nPhi) * 2 * M_PI - M_PI; // shift to [-pi,pi] - - // determine global phi and theta - // - convert {radius,thetaGen,phiGen} -> {xGen,yGen,zGen} - double xGen = sensorSphRadius * std::sin(thetaGen) * std::cos(phiGen); - double yGen = sensorSphRadius * std::sin(thetaGen) * std::sin(phiGen); - double zGen = sensorSphRadius * std::cos(thetaGen); - // - convert {xGen,yGen,zGen} -> global {x,y,z} via rotation - double x = zGen; - double y = xGen; - double z = yGen; - // - convert global {x,y,z} -> global {phi,theta} - // double phi = std::atan2(y, x); - // double theta = std::acos(z / sensorSphRadius); - - // shift global coordinates so we can apply spherical patch cuts - double zCheck = z + sensorSphCenterZ; - double xCheck = x + sensorSphCenterX; - double yCheck = y; - double rCheck = std::hypot(xCheck, yCheck); - double phiCheck = std::atan2(yCheck, xCheck); - - // patch cut - bool patchCut = std::fabs(phiCheck) < sensorSphPatchPhiw && zCheck > sensorSphPatchZmin && - rCheck > sensorSphPatchRmin && rCheck < sensorSphPatchRmax; - if (debugSensors) - patchCut = std::fabs(phiCheck) < sensorSphPatchPhiw; - if (patchCut) { - - // append sensor position to centroid calculation - // if (isec == 0) { - // sensorCentroidX += xCheck; - // sensorCentroidZ += zCheck; - // sensorCount++; - // }; - - // placement (note: transformations are in reverse order) - // - transformations operate on global coordinates; the corresponding - // generator coordinates are provided in the comments - auto sensorPlacement = + double thetaGen = t / ((double)nTheta) * M_PI; + + // phiGen loop: iterate less than "circumference at this latitude / sensor size" times + double nPhi = 2 * M_PI * sensorSphRadius * std::sin(thetaGen) / (sensorSide + sensorGap); + for (int p = 0; p < (int)(nPhi + 0.5); p++) { + double phiGen = p / ((double)nPhi) * 2 * M_PI - M_PI; // shift to [-pi,pi] + + // determine global phi and theta + // - convert {radius,thetaGen,phiGen} -> {xGen,yGen,zGen} + double xGen = sensorSphRadius * std::sin(thetaGen) * std::cos(phiGen); + double yGen = sensorSphRadius * std::sin(thetaGen) * std::sin(phiGen); + double zGen = sensorSphRadius * std::cos(thetaGen); + // - convert {xGen,yGen,zGen} -> global {x,y,z} via rotation + double x = zGen; + double y = xGen; + double z = yGen; + // - convert global {x,y,z} -> global {phi,theta} + // double phi = std::atan2(y, x); + // double theta = std::acos(z / sensorSphRadius); + + // shift global coordinates so we can apply spherical patch cuts + double zCheck = z + sensorSphCenterZ; + double xCheck = x + sensorSphCenterX; + double yCheck = y; + double rCheck = std::hypot(xCheck, yCheck); + double phiCheck = std::atan2(yCheck, xCheck); + + // patch cut + bool patchCut = std::fabs(phiCheck) < sensorSphPatchPhiw && zCheck > sensorSphPatchZmin && + rCheck > sensorSphPatchRmin && rCheck < sensorSphPatchRmax; + if (debugSensors) + patchCut = std::fabs(phiCheck) < sensorSphPatchPhiw; + if (patchCut) { + + // append sensor position to centroid calculation + // if (isec == 0) { + // sensorCentroidX += xCheck; + // sensorCentroidZ += zCheck; + // sensorCount++; + // } + + // placement (note: transformations are in reverse order) + // - transformations operate on global coordinates; the corresponding + // generator coordinates are provided in the comments + auto sensorPlacement = sectorRotation * // rotate about beam axis to sector Translation3D(sensorSphPos.x(), sensorSphPos.y(), sensorSphPos.z()) * // move sphere to reference position RotationX(phiGen) * // rotate about `zGen` RotationZ(thetaGen) * // rotate about `yGen` Translation3D(-sensorThickness / 2.0, 0., - 0.) * // pull back so sensor active surface is at spherical surface + 0.) * // pull back so sensor active surface is at spherical surface Translation3D(sensorSphRadius, 0., 0.) * // push radially to spherical surface RotationY(M_PI / 2) * // rotate sensor to be compatible with generator coords RotationZ(-M_PI / 2); // correction for readout segmentation mapping - auto sensorPV = gasvolVol.placeVolume(sensorVol, sensorPlacement); + auto sensorPV = gasvolVol.placeVolume(sensorVol, sensorPlacement); - // generate LUT for module number -> sensor position, for readout mapping tests - // if(isec==0) printf("%d %f %f\n",imod,sensorPV.position().x(),sensorPV.position().y()); + // generate LUT for module number -> sensor position, for readout mapping tests + // if(isec==0) printf("%d %f %f\n",imod,sensorPV.position().x(),sensorPV.position().y()); // properties sensorPV.addPhysVolID("sector", isec) @@ -589,15 +584,15 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec if (!debugOptics || debugOpticsMode == 3 || debugOpticsMode == 5) { SkinSurface sensorSkin(desc, sensorDE, "sensor_optical_surface_" + modsecName, sensorSurf, sensorVol); sensorSkin.isValid(); - }; + } // increment sensor module number imod++; - }; // end patch cuts - }; // end phiGen loop - }; // end thetaGen loop - }; // end sensor loop + } // end patch cuts + } // end phiGen loop + } // end thetaGen loop + } // end if(debugOpticsMode != 4) // large sensor placement (for focal point testing): if (debugOpticsMode == 4){ @@ -621,7 +616,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec SkinSurface sensorSkin(desc, sensorDE, Form("sensor_optical_surface%d", isec), sensorSurf, sensorVolScaled); sensorSkin.isValid(); - }; // end large sensor placement + } // end large sensor placement // for debugOpticsMode 5: drawing calculated focal points // of thrown beams of parallel photons @@ -659,19 +654,19 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec DetElement FPDE(det, Form("FP_de%d_%d", isec, fpnum), 0); FPDE.setPlacement(FPPV); fpnum++; - }; - }; - }; + } + } + } // calculate centroid sensor position // if (isec == 0) { // sensorCentroidX /= sensorCount; // sensorCentroidZ /= sensorCount; - // }; + // } // END SENSOR MODULE LOOP ------------------------ - }; // END SECTOR LOOP ////////////////////////// + } // END SECTOR LOOP ////////////////////////// return det; }