File indexing completed on 2026-03-29 07:51:20
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029 #include "DetectorConstruction.hh"
0030
0031 #include "DetectorMessenger.hh"
0032
0033 #include "G4GeometryManager.hh"
0034 #include "G4LogicalVolume.hh"
0035 #include "G4LogicalVolumeStore.hh"
0036 #include "G4Material.hh"
0037 #include "G4NistManager.hh"
0038 #include "G4PVPlacement.hh"
0039 #include "G4PhysicalConstants.hh"
0040 #include "G4PhysicalVolumeStore.hh"
0041 #include "G4ReflectionFactory.hh"
0042 #include "G4RotationMatrix.hh"
0043 #include "G4SolidStore.hh"
0044 #include "G4SystemOfUnits.hh"
0045 #include "G4Transform3D.hh"
0046 #include "G4Trd.hh"
0047 #include "G4Tubs.hh"
0048
0049
0050
0051 DetectorConstruction::DetectorConstruction()
0052 {
0053 fMessenger = new DetectorMessenger(this);
0054 }
0055
0056
0057
0058 DetectorConstruction::~DetectorConstruction()
0059 {
0060 delete fMessenger;
0061 }
0062
0063
0064
0065 G4VPhysicalVolume* DetectorConstruction::Construct()
0066 {
0067
0068 G4NistManager* nist = G4NistManager::Instance();
0069 G4Material* material = nist->FindOrBuildMaterial("G4_AIR");
0070
0071
0072
0073 G4GeometryManager::GetInstance()->OpenGeometry();
0074 G4PhysicalVolumeStore::GetInstance()->Clean();
0075 G4LogicalVolumeStore::GetInstance()->Clean();
0076 G4SolidStore::GetInstance()->Clean();
0077 G4ReflectionFactory::Instance()->Clean();
0078
0079
0080
0081 G4double rmin = 0.;
0082 G4double rmax = 5 * cm;
0083 G4double hz = 5 * cm;
0084 G4double phiMin = 0.;
0085 G4double deltaPhi = 360 * degree;
0086
0087 auto solidWorld = new G4Tubs("World",
0088 rmin, rmax, hz, phiMin, deltaPhi);
0089
0090 fWorldVolume = new G4LogicalVolume(solidWorld,
0091 material,
0092 "World");
0093
0094 G4VPhysicalVolume* physiWorld = new G4PVPlacement(nullptr,
0095 G4ThreeVector(),
0096 fWorldVolume,
0097 "World",
0098 nullptr,
0099 false,
0100 0);
0101
0102
0103
0104 G4double dX1 = 1 * cm;
0105 G4double dX2 = 1 * cm;
0106 G4double dY1 = 1 * cm;
0107 G4double dY2 = 2 * cm;
0108 G4double dZ = 3 * cm;
0109
0110 auto solidTrd = new G4Trd("trd",
0111 dX1 / 2, dX2 / 2, dY1 / 2, dY2 / 2, dZ / 2);
0112
0113 fTrdVolume = new G4LogicalVolume(solidTrd,
0114 material,
0115 "trd");
0116
0117
0118
0119 switch (fMethod) {
0120 case kWithDirectMatrix:
0121 PlaceWithDirectMatrix();
0122 break;
0123 case kWithInverseMatrix:
0124 PlaceWithInverseMatrix();
0125 break;
0126 case kWithAxialRotations:
0127 PlaceWithAxialRotations();
0128 break;
0129 case kWithEulerAngles:
0130 PlaceWithEulerAngles();
0131 break;
0132 case kWithReflections:
0133 PlaceWithReflections();
0134 break;
0135 default:;
0136 ;
0137 }
0138
0139
0140
0141 return physiWorld;
0142 }
0143
0144
0145
0146 void DetectorConstruction::PlaceWithDirectMatrix()
0147 {
0148 G4double og = 3 * cm;
0149
0150
0151
0152 G4double phi = 30 * deg;
0153
0154 G4ThreeVector u = G4ThreeVector(0, 0, -1);
0155 G4ThreeVector v = G4ThreeVector(-std::sin(phi), std::cos(phi), 0.);
0156 G4ThreeVector w = G4ThreeVector(std::cos(phi), std::sin(phi), 0.);
0157 G4RotationMatrix rotm1 = G4RotationMatrix(u, v, w);
0158 G4cout << "\n --> phi = " << phi / deg << " deg; direct rotation matrix : ";
0159 rotm1.print(G4cout);
0160 G4ThreeVector position1 = og * w;
0161 G4Transform3D transform1 = G4Transform3D(rotm1, position1);
0162
0163 new G4PVPlacement(transform1,
0164 fTrdVolume,
0165 "Trd",
0166 fWorldVolume,
0167 false,
0168 1);
0169
0170
0171
0172 phi = phi + 90 * deg;
0173 v = G4ThreeVector(-std::sin(phi), std::cos(phi), 0.);
0174 w = G4ThreeVector(std::cos(phi), std::sin(phi), 0.);
0175 G4RotationMatrix rotm2 = G4RotationMatrix(u, v, w);
0176 G4ThreeVector position2 = og * w;
0177 G4Transform3D transform2 = G4Transform3D(rotm2, position2);
0178 new G4PVPlacement(transform2,
0179 fTrdVolume,
0180 "Trd",
0181 fWorldVolume,
0182 false,
0183 2);
0184 }
0185
0186
0187
0188 void DetectorConstruction::PlaceWithInverseMatrix()
0189 {
0190 G4double og = 3 * cm;
0191
0192
0193
0194 G4double phi = 30 * deg;
0195
0196 G4ThreeVector u = G4ThreeVector(0, 0, -1);
0197 G4ThreeVector v = G4ThreeVector(-std::sin(phi), std::cos(phi), 0.);
0198 G4ThreeVector w = G4ThreeVector(std::cos(phi), std::sin(phi), 0.);
0199 G4RotationMatrix rotm1 = G4RotationMatrix(u, v, w);
0200 auto rotm1Inv = new G4RotationMatrix(rotm1.inverse());
0201 G4cout << "\n --> phi = " << phi / deg << " deg; inverse rotation matrix : ";
0202 rotm1Inv->print(G4cout);
0203 G4ThreeVector position1 = og * w;
0204
0205 new G4PVPlacement(rotm1Inv, position1,
0206 fTrdVolume,
0207 "Trd",
0208 fWorldVolume,
0209 false,
0210 1);
0211
0212
0213
0214 phi = phi + 90 * deg;
0215 v = G4ThreeVector(-std::sin(phi), std::cos(phi), 0.);
0216 w = G4ThreeVector(std::cos(phi), std::sin(phi), 0.);
0217 G4RotationMatrix rotm2 = G4RotationMatrix(u, v, w);
0218 auto rotm2Inv = new G4RotationMatrix(rotm2.inverse());
0219 G4ThreeVector position2 = og * w;
0220
0221 new G4PVPlacement(rotm2Inv,
0222 position2,
0223 fTrdVolume,
0224 "Trd",
0225 fWorldVolume,
0226 false,
0227 2);
0228 }
0229
0230
0231
0232 void DetectorConstruction::PlaceWithAxialRotations()
0233 {
0234 G4double og = 3 * cm;
0235
0236
0237
0238 G4double phi = 30 * deg, theta = 90 * deg;
0239 G4ThreeVector rotAxis = G4ThreeVector(std::sin(theta - pi / 2), 0., std::cos(theta - pi / 2));
0240 G4RotationMatrix rotm1 = G4RotationMatrix();
0241 rotm1.rotateY(theta);
0242 rotm1.rotate(phi, rotAxis);
0243 G4cout << "\n --> direct rotation matrix : "
0244 << " theta = " << theta / deg << " deg;"
0245 << " phi = " << phi / deg << " deg;";
0246 rotm1.print(G4cout);
0247 G4ThreeVector w =
0248 G4ThreeVector(std::sin(theta) * std::cos(phi), std::sin(phi), std::cos(theta) * std::cos(phi));
0249 G4ThreeVector position1 = og * w;
0250 G4Transform3D transform1(rotm1, position1);
0251
0252 new G4PVPlacement(transform1,
0253 fTrdVolume,
0254 "Trd",
0255 fWorldVolume,
0256 false,
0257 1);
0258
0259
0260
0261 phi = phi + 90 * deg;
0262
0263
0264 auto rotm2Inv = new G4RotationMatrix();
0265 rotm2Inv->rotate(-phi, rotAxis);
0266 rotm2Inv->rotateY(-theta);
0267 w =
0268 G4ThreeVector(std::sin(theta) * std::cos(phi), std::sin(phi), std::cos(theta) * std::cos(phi));
0269 G4ThreeVector position2 = og * w;
0270
0271 new G4PVPlacement(rotm2Inv,
0272 position2,
0273 fTrdVolume,
0274 "Trd",
0275 fWorldVolume,
0276 false,
0277 2);
0278 }
0279
0280
0281
0282 void DetectorConstruction::PlaceWithEulerAngles()
0283 {
0284
0285
0286
0287
0288
0289
0290 G4double og = 3 * cm;
0291
0292
0293
0294 G4double phi = 30 * deg;
0295 G4double phi_euler = phi + pi / 2;
0296 G4double theta_euler = 90 * deg;
0297 G4double psi_euler = -90 * deg;
0298
0299 G4RotationMatrix rotm1Inv = G4RotationMatrix(phi_euler, theta_euler, psi_euler);
0300 G4RotationMatrix rotm1 = rotm1Inv.inverse();
0301
0302 G4cout << "\n --> phi = " << phi / deg << " deg; direct rotation matrix : ";
0303 rotm1.print(G4cout);
0304 G4ThreeVector w = G4ThreeVector(std::cos(phi), std::sin(phi), 0.);
0305 G4ThreeVector position1 = og * w;
0306 G4Transform3D transform1 = G4Transform3D(rotm1, position1);
0307
0308 new G4PVPlacement(transform1,
0309 fTrdVolume,
0310 "Trd",
0311 fWorldVolume,
0312 false,
0313 1);
0314
0315
0316
0317 phi = phi + 90 * deg;
0318
0319 phi_euler = phi + pi / 2;
0320 auto rotm2Inv = new G4RotationMatrix(phi_euler, theta_euler, psi_euler);
0321 w = G4ThreeVector(std::cos(phi), std::sin(phi), 0.);
0322 G4ThreeVector position2 = og * w;
0323
0324 new G4PVPlacement(rotm2Inv,
0325 position2,
0326 fTrdVolume,
0327 "Trd",
0328 fWorldVolume,
0329 false,
0330 2);
0331 }
0332
0333
0334
0335 void DetectorConstruction::PlaceWithReflections()
0336 {
0337
0338
0339
0340
0341 G4double og = 3 * cm;
0342
0343
0344
0345
0346
0347 G4double phi = 30 * deg;
0348 G4RotationMatrix rotm1;
0349
0350 rotm1.rotateZ(phi);
0351 G4ThreeVector uz = G4ThreeVector(std::cos(phi), std::sin(phi), 0);
0352 G4ThreeVector position = og * uz;
0353 G4Transform3D transform1(rotm1, position);
0354 G4Transform3D translateZ = HepGeom::Translate3D(0, 0, 3. * cm);
0355
0356 new G4PVPlacement(translateZ * transform1,
0357 fTrdVolume,
0358 "Trd",
0359 fWorldVolume,
0360 false,
0361 1);
0362
0363
0364 phi = phi + pi / 2;
0365 G4RotationMatrix rotm2;
0366
0367 rotm2.rotateZ(phi);
0368 uz = G4ThreeVector(std::cos(phi), std::sin(phi), 0.);
0369 position = og * uz;
0370 G4Transform3D transform2 = G4Transform3D(rotm2, position);
0371
0372 new G4PVPlacement(translateZ * transform2,
0373 fTrdVolume,
0374 "Trd",
0375 fWorldVolume,
0376 false,
0377 2);
0378
0379
0380
0381
0382
0383 translateZ = HepGeom::Translate3D(0, 0, -3. * cm);
0384 G4Transform3D reflect3D = HepGeom::ReflectZ3D();
0385
0386 G4ReflectionFactory::Instance()->Place(translateZ * transform1 * reflect3D,
0387 "Trd",
0388 fTrdVolume,
0389 fWorldVolume,
0390 false,
0391 3);
0392
0393
0394 G4ReflectionFactory::Instance()->Place(translateZ * transform2 * reflect3D,
0395 "Trd",
0396 fTrdVolume,
0397 fWorldVolume,
0398 false,
0399 4);
0400 }
0401
0402
0403
0404 #include "G4RunManager.hh"
0405
0406 void DetectorConstruction::SetMethod(EMethod method)
0407 {
0408 fMethod = method;
0409 G4RunManager::GetRunManager()->DefineWorldVolume(Construct());
0410 }
0411
0412