21 #include <phparameter/PHParameters.h>
27 #include <Geant4/G4AssemblyVolume.hh>
28 #include <Geant4/G4Box.hh>
29 #include <Geant4/G4LogicalBorderSurface.hh>
30 #include <Geant4/G4LogicalVolume.hh>
31 #include <Geant4/G4Material.hh>
32 #include <Geant4/G4MaterialPropertiesTable.hh>
33 #include <Geant4/G4OpticalSurface.hh>
34 #include <Geant4/G4PVPlacement.hh>
35 #include <Geant4/G4PhysicalConstants.hh>
36 #include <Geant4/G4Polycone.hh>
37 #include <Geant4/G4Polyhedra.hh>
38 #include <Geant4/G4RotationMatrix.hh>
39 #include <Geant4/G4SurfaceProperty.hh>
40 #include <Geant4/G4SystemOfUnits.hh>
41 #include <Geant4/G4ThreeVector.hh>
42 #include <Geant4/G4Types.hh>
43 #include <Geant4/G4VPhysicalVolume.hh>
44 #include <Geant4/G4VisAttributes.hh>
45 #include <Geant4/G4ios.hh>
56 using namespace CLHEP;
66 , sensor_PV{
nullptr,
nullptr,
nullptr,
nullptr}
96 if (subsystemSetup == DetectorSetUp::kSingle_Modular)
Construct_a_mRICH(logicWorld);
97 if (subsystemSetup == DetectorSetUp::kHSector_EWall)
106 if (subsystemSetup == DetectorSetUp::kHWall_EWall)
111 if (subsystemSetup == DetectorSetUp::kHWall_Barrel)
121 mRichParameter*
parameters =
new mRichParameter();
178 std::fill(std::begin(
z), std::end(
z), (
G4double) 0 *
mm);
179 std::fill(std::begin(rinner), std::end(rinner), (
G4double) 0 * mm);
180 std::fill(std::begin(router), std::end(router), (
G4double) 0 * mm);
189 , centerThickness(0.)
206 halfXYZ[1] = halfXYZ[0];
208 G4double NumberOfGrooves = floor(grooveDensity * (eff_diameter / 2.0));
209 G4double Rmin1 = (NumberOfGrooves - 1) * (grooveWidth);
210 G4double Rmax1 = (NumberOfGrooves - 0) * (grooveWidth);
211 halfXYZ[2] = (GetSagita(Rmax1) - GetSagita(Rmin1) + centerThickness) / 2.0;
220 G4double Aspher[8] = {0, 0, 0, 0, 0, 0, 0, 0};
224 Curvature = 0.00437636761488 /
mm;
225 Aspher[0] = 4.206739256e-05 / (
mm);
226 Aspher[1] = 9.6440152e-10 / (
mm3);
227 Aspher[2] = -1.4884317e-15 / (
mm2 *
mm3);
229 else if (lens_type == 2)
231 Curvature = 0.0132 /
mm;
232 Aspher[0] = 32.0e-05 / (
mm);
233 Aspher[1] = -2.0e-7 / (
mm3);
234 Aspher[2] = 1.2e-13 / (
mm2 *
mm3);
236 else if (lens_type == 3)
238 Curvature = 0.0150 /
mm;
239 Aspher[0] = 42.0e-05 / (
mm);
240 Aspher[1] = -3.0e-7 / (
mm3);
241 Aspher[2] = 1.2e-13 / (
mm2 *
mm3);
243 else if (lens_type == 4)
245 Curvature = 0.0175 /
mm;
246 Aspher[0] = 72.0e-05 / (
mm);
247 Aspher[1] = -5.0e-7 / (
mm3);
248 Aspher[2] = 1.2e-13 / (
mm2 *
mm3);
250 else if (lens_type == 5)
252 Curvature = 1 / (
f * (
n - 1));
258 TotAspher += Aspher[
k - 1] * std::pow(r, 2 *
k);
261 G4double ArgSqrt = 1.0 - (1.0 + Conic) * std::pow(Curvature, 2) * std::pow(r, 2);
268 G4double Sagita_value = Curvature * std::pow(r, 2) / (1.0 + std::sqrt(ArgSqrt)) + TotAspher;
281 const double myPI = 4 * atan(1);
284 hollowVolume =
new BoxPar();
285 foamHolderBox =
new BoxPar();
286 foamHolderPoly =
new PolyPar();
289 glassWindow =
new BoxPar();
295 const G4double box_thicknessXYZ[4] = {0.1 *
cm, 0.1 *
cm, 0.1 *
cm, 0.1 *
cm};
298 const G4double foamHolderThicknessXYZ[3] = {0.2 *
cm, 0.2 *
cm, 0.2 *
cm};
299 const G4double agel_halfXYZ[3] = {6.325 *
cm, 6.325 *
cm, 1.50 *
cm};
304 const G4double lensHalfx = (5.25 / 2.0) * 2.54 *
cm;
305 const G4double grooveDensity = 125.0 / (2.54 *
cm);
307 fresnelLens->n = 1.49;
308 fresnelLens->f = 6.0 * 2.54 *
cm;
309 fresnelLens->eff_diameter = 15.24 *
cm;
310 fresnelLens->diameter = 2.0 * sqrt(2.0) * lensHalfx;
311 fresnelLens->centerThickness = 0.068 * 2.54 *
cm;
312 fresnelLens->grooveWidth = (
G4double) 1.0 / grooveDensity;
313 fresnelLens->Set_halfXYZ(lensHalfx, grooveDensity);
318 const G4double glassWindow_halfXYZ[3] = {5.18 / 2.0 *
cm, 5.18 / 2.0 *
cm, 0.075 *
cm};
319 const G4double phodet_halfXYZ[3] = {2.425 *
cm, 2.425 *
cm, 0.075 *
cm};
331 G4double sensor_total_halfx = 2 * glassWindow_halfXYZ[0] + sensorGap;
334 foamHolder_halfXYZ[0] = agel_halfXYZ[0] + foamHolderThicknessXYZ[0];
335 foamHolder_halfXYZ[1] = foamHolder_halfXYZ[0];
336 foamHolder_halfXYZ[2] = foamHolderThicknessXYZ[2] / 2.0;
339 acrylicBox_halfXYZ[0] =
std::max(
std::max(foamHolder_halfXYZ[0], sensor_total_halfx + readoutThickness), fresnelLens->halfXYZ[0]) + 0.1 *
cm + box_thicknessXYZ[0];
340 acrylicBox_halfXYZ[1] = acrylicBox_halfXYZ[0];
341 acrylicBox_halfXYZ[2] = (BoxDelz + 2 * foamHolder_halfXYZ[2] + 2 * agel_halfXYZ[2] +
342 lens_gap + 2 * fresnelLens->halfXYZ[2] + fresnelLens->f + 2 * glassWindow_halfXYZ[2] +
343 2 * phodet_halfXYZ[2] + (2 * readout_halfz + BoxDelz) + box_thicknessXYZ[2] + box_thicknessXYZ[3]) /
347 hollow_halfXYZ[0] = acrylicBox_halfXYZ[0] - box_thicknessXYZ[0];
348 hollow_halfXYZ[1] = hollow_halfXYZ[0];
349 hollow_halfXYZ[2] = (2 * acrylicBox_halfXYZ[2] - box_thicknessXYZ[2] - box_thicknessXYZ[3]) / 2.0;
353 G4double foamHolder_posz = -hollow_halfXYZ[2] + BoxDelz + foamHolder_halfXYZ[2];
354 G4double agel_posz = foamHolder_posz + foamHolder_halfXYZ[2] + agel_halfXYZ[2];
356 G4double lens_z = agel_posz + agel_halfXYZ[2] + fresnelLens->halfXYZ[2] + lens_gap;
358 G4double glassWindow_z = lens_z - fresnelLens->halfXYZ[2] + fresnelLens->f + glassWindow_halfXYZ[2] - 16.;
359 G4double phodet_z = glassWindow_z + glassWindow_halfXYZ[2] + phodet_halfXYZ[2];
363 readout_z[0] = glassWindow_z - glassWindow_halfXYZ[2];
364 readout_z[1] = phodet_z + phodet_halfXYZ[2];
370 holderBox->name =
"mRICH_module";
371 for (i = 0; i < 3; i++) holderBox->halfXYZ[i] = acrylicBox_halfXYZ[i];
376 holderBox->sensitivity = 0;
378 holderBox->color =
G4Colour(0.0, 0.0, 0.0);
379 holderBox->visibility =
true;
380 holderBox->wireframe =
true;
381 holderBox->surface =
false;
386 hollowVolume->name =
"HollowVolume";
387 for (i = 0; i < 3; i++) hollowVolume->halfXYZ[i] = hollow_halfXYZ[i];
388 hollowVolume->pos = hollow_pos;
391 hollowVolume->sensitivity = 0;
393 hollowVolume->color =
G4Colour(0.0, 0.0, 0.0);
394 hollowVolume->visibility =
true;
395 hollowVolume->wireframe =
true;
396 hollowVolume->surface =
false;
401 foamHolderBox->name =
"FoamHolder";
402 for (i = 0; i < 3; i++) foamHolderBox->halfXYZ[i] = foamHolder_halfXYZ[i];
405 foamHolderBox->sensitivity = 0;
407 foamHolderBox->color =
G4Colour(0.2, 0.498, 0.369);
408 foamHolderBox->visibility =
true;
409 foamHolderBox->wireframe =
true;
410 foamHolderBox->surface =
false;
415 foamHolderPoly->name =
"FoamHolder";
417 foamHolderPoly->start = 45.0 * myPI / 180.0;
418 foamHolderPoly->theta = 2 * myPI;
419 foamHolderPoly->numSide = 4;
420 foamHolderPoly->num_zLayer = 2;
422 foamHolderPoly->z[0] = agel_posz - agel_halfXYZ[2];
423 foamHolderPoly->z[1] = agel_posz + agel_halfXYZ[2];
425 foamHolderPoly->rinner[0] = agel_halfXYZ[0];
426 foamHolderPoly->rinner[1] = agel_halfXYZ[0];
428 foamHolderPoly->router[0] = foamHolderPoly->rinner[0] + foamHolderThicknessXYZ[0];
429 foamHolderPoly->router[1] = foamHolderPoly->rinner[1] + foamHolderThicknessXYZ[0];
432 foamHolderPoly->sensitivity = 0;
434 foamHolderPoly->color =
G4Colour(0.298, 0.6, 0.471);
435 foamHolderPoly->visibility =
true;
436 foamHolderPoly->wireframe =
true;
437 foamHolderPoly->surface =
false;
443 aerogel->name =
"Aerogel";
444 for (i = 0; i < 3; i++) aerogel->halfXYZ[i] = agel_halfXYZ[i];
447 aerogel->sensitivity = 0;
448 aerogel->color =
G4Colour(1.0, 0.65, 0.0);
449 aerogel->visibility =
true;
450 aerogel->wireframe =
true;
451 aerogel->surface =
false;
459 fresnelLens->name =
"FresnelLens";
462 fresnelLens->sensitivity = 0;
463 fresnelLens->color =
G4Colour(0.0, 1.0, 1.0);
464 fresnelLens->visibility =
true;
465 fresnelLens->wireframe =
true;
466 fresnelLens->surface =
false;
471 mirror->name =
"mirror";
473 mirror->start = 45.0 * myPI / 180.0;
474 mirror->theta = 2 * myPI;
476 mirror->num_zLayer = 2;
478 mirror->z[0] = lens_z + fresnelLens->halfXYZ[2] + lens_gap;
479 mirror->z[1] = glassWindow_z - glassWindow_halfXYZ[2];
481 mirror->rinner[0] = agel_halfXYZ[0];
482 mirror->rinner[1] = sensor_total_halfx;
484 mirror->router[0] = mirror->rinner[0] + mirrorThickness;
485 mirror->router[1] = mirror->rinner[1] + mirrorThickness;
488 mirror->sensitivity = 0;
490 mirror->color =
G4Colour(1.0, 1.0, 0.0);
491 mirror->visibility =
true;
492 mirror->wireframe =
true;
497 glassWindow->name =
"glassWindow";
498 for (i = 0; i < 3; i++) glassWindow->halfXYZ[i] = glassWindow_halfXYZ[i];
499 glassWindow->pos =
G4ThreeVector(glassWindow_halfXYZ[0] + sensorGap,
500 glassWindow_halfXYZ[0] + sensorGap,
503 glassWindow->sensitivity = 0;
505 glassWindow->color =
G4Colour(0.101, 0.737, 0.612);
506 glassWindow->visibility =
true;
507 glassWindow->wireframe =
true;
508 glassWindow->surface =
false;
513 sensor->name =
"sensor";
514 for (i = 0; i < 3; i++) sensor->halfXYZ[i] = phodet_halfXYZ[i];
517 sensor->sensitivity = 0;
519 sensor->color =
G4Colour(0.0, 0.0, 0.63);
520 sensor->visibility =
true;
521 sensor->wireframe =
true;
522 sensor->surface =
true;
527 readout->name =
"readout";
529 readout->start = 45.0 * myPI / 180.0;
530 readout->theta = 2 * myPI;
531 readout->numSide = 4;
532 readout->num_zLayer = 2;
534 readout->z[0] = readout_z[0];
535 readout->z[1] = readout_z[1];
537 readout->rinner[0] = sensor_total_halfx;
538 readout->rinner[1] = readout->rinner[0];
540 readout->router[0] = readout->rinner[0] + readoutThickness;
541 readout->router[1] = readout->router[0];
544 readout->sensitivity = 0;
546 readout->color =
G4Colour(1.0, 0.0, 0.0);
547 readout->visibility =
true;
548 readout->wireframe =
true;
549 readout->surface =
false;
558 glassWindow->pos.setX(x);
559 glassWindow->pos.setY(y);
572 if (componentName.compare(
"holderBox") == 0)
574 else if (componentName.compare(
"hollowVolume") == 0)
576 else if (componentName.compare(
"foamHolderBox") == 0)
577 return foamHolderBox;
578 else if (componentName.compare(
"aerogel") == 0)
580 else if (componentName.compare(
"glassWindow") == 0)
582 else if (componentName.compare(
"sensor") == 0)
585 std::cout << __FILE__ <<
"::" << __func__ <<
":: ERROR: cannot find parameter " << componentName << std::endl;
592 if (componentName.compare(
"fresnelLens") == 0)
595 std::cout << __FILE__ <<
"::" << __func__ <<
"::ERROR: cannot find parameter " << componentName << std::endl;
601 if (componentName.compare(
"foamHolderPoly") == 0)
602 return foamHolderPoly;
603 else if (componentName.compare(
"mirror") == 0)
605 else if (componentName.compare(
"readout") == 0)
608 std::cout << __FILE__ <<
"::" << __func__ <<
"::ERROR: cannot find parameter " << componentName << std::endl;
621 visAtt->SetForceWireframe(par->
wireframe);
622 visAtt->SetForceSolid(par->
surface);
623 log->SetVisAttributes(visAtt);
637 visAtt->SetForceWireframe(par->
wireframe);
638 visAtt->SetForceSolid(par->
surface);
639 log->SetVisAttributes(visAtt);
683 myST1->
AddProperty(
"RINDEX", Ephoton, RefractiveIndex, num);
684 myST1->
AddProperty(
"SPECULARLOBECONSTANT", Ephoton, SpecularLobe, num);
685 myST1->
AddProperty(
"SPECULARSPIKECONSTANT", Ephoton, SpecularSpike, num);
686 myST1->
AddProperty(
"BACKSCATTERCONSTANT", Ephoton, Backscatter, num);
700 OpticalAirMirror->SetModel(
unified);
709 G4double ICEREFLECTIVITY[NUM] = {0.95, 0.95};
712 AirMirrorMPT->
AddProperty(
"REFLECTIVITY", XX, ICEREFLECTIVITY, NUM);
713 OpticalAirMirror->SetMaterialPropertiesTable(AirMirrorMPT);
727 for (i = 0; i < 4; i++)
765 for (igroove = 0; igroove < 1000; igroove++)
775 G4double lens_poly_rmin[3] = {iRmin1, iRmin1, iRmin2};
776 G4double lens_poly_rmax[3] = {iRmax1, iRmax1, iRmax2};
778 if (iRmax1 > par->
diameter / 2.0)
break;
787 if (iRmax1 < par->halfXYZ[0])
794 phi1 = acos(par->
halfXYZ[0] / iRmax1);
795 phi2 = asin(par->
halfXYZ[0] / iRmax1);
796 deltaPhi = phi2 - phi1;
804 if (iRmin1 < par->eff_diameter / 2.0)
808 lens_poly_z[0] = par->
halfXYZ[2];
809 lens_poly_z[1] = -par->
halfXYZ[2] + dZ;
810 lens_poly_z[2] = -par->
halfXYZ[2];
815 lens_poly_z[0] = par->
halfXYZ[2];
828 for (
int i = 0; i < repeat; i++)
830 Groove_poly[i] =
new G4Polycone(par->
name.c_str(), phi1, deltaPhi, numOfLayer, lens_poly_z, lens_poly_rmin, lens_poly_rmax);
834 Groove_log[i]->SetVisAttributes(SurfaceVisAtt);
848 for (
int i_mRICH = 0; i_mRICH < NumOfModule; ++i_mRICH)
856 std::stringstream key_position_x;
857 key_position_x <<
"mRICH_wall_hside_" << i_mRICH <<
"_position_x";
860 std::stringstream key_position_y;
861 key_position_y <<
"mRICH_wall_hside_" << i_mRICH <<
"_position_y";
864 std::stringstream key_position_z;
865 key_position_z <<
"mRICH_wall_hside_" << i_mRICH <<
"_position_z";
869 std::stringstream key_rotation_theta;
870 key_rotation_theta <<
"mRICH_wall_hside_" << i_mRICH <<
"_rotation_theta";
873 std::stringstream key_rotation_phi;
874 key_rotation_phi <<
"mRICH_wall_hside_" << i_mRICH <<
"_rotation_phi";
882 if (x != 0 || y != 0)
908 for (
int i_mRICH = 0; i_mRICH < NumOfModule; ++i_mRICH)
916 std::stringstream key_position_x;
917 key_position_x <<
"mRICH_wall_eside_" << i_mRICH <<
"_position_x";
920 std::stringstream key_position_y;
921 key_position_y <<
"mRICH_wall_eside_" << i_mRICH <<
"_position_y";
924 std::stringstream key_position_z;
925 key_position_z <<
"mRICH_wall_eside_" << i_mRICH <<
"_position_z";
953 for (
int i_mRICH = 0; i_mRICH < NumOfModule; ++i_mRICH)
961 std::stringstream key_position_x;
962 key_position_x <<
"mRICH_sector_hside_" << i_mRICH <<
"_position_x";
965 std::stringstream key_position_y;
966 key_position_y <<
"mRICH_sector_hside_" << i_mRICH <<
"_position_y";
969 std::stringstream key_position_z;
970 key_position_z <<
"mRICH_sector_hside_" << i_mRICH <<
"_position_z";
972 if (i_mRICH == 10) z -= 10.;
982 for (
int i = 0; i < numSector; i++)
1009 G4double yy[8] = {61.92, 61.92, 63.09, 63.09, 63.34, 63.34, 63.00, 63.00};
1010 G4double zz[8] = {8.94, -8.94, 27.78, -27.78, 48.26, -48.26, 71.41, -71.41};
1011 G4double rotAng[8] = {-81.78, -98.22, -66.23, -113.77, -52.69, -127.31, -41.42, -138.58};
1015 for (
int i_mRICH = 0; i_mRICH < 8; ++i_mRICH)
1023 std::stringstream key_position_x;
1024 key_position_x <<
"mRICH_sector_bside_" << i_mRICH <<
"_position_x";
1027 std::stringstream key_position_y;
1028 key_position_y <<
"mRICH_sector_bside_" << i_mRICH <<
"_position_y";
1031 std::stringstream key_position_z;
1032 key_position_z <<
"mRICH_sector_bside_" << i_mRICH <<
"_position_z";
1054 for (
int i = 0; i < nSecs; i++)
1082 for (
int i_mRICH = 0; i_mRICH < NumOfModule; ++i_mRICH)
1098 else if (i_mRICH >= 8 && i_mRICH < 20)
1100 else if (i_mRICH >= 20 && i_mRICH < 44)
1106 std::stringstream key_position_x;
1107 key_position_x <<
"mRICH_wall_eside_proj_" << i_mRICH <<
"_position_x";
1110 std::stringstream key_position_y;
1111 key_position_y <<
"mRICH_wall_eside_proj_" << i_mRICH <<
"_position_y";
1114 std::stringstream key_position_z;
1115 key_position_z <<
"mRICH_wall_eside_proj_" << i_mRICH <<
"_position_z";
1119 G4double rotAngX = atan(y / shift);