nanobeam/src/DetectorConstruction.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 // Please cite the following paper if you use this software
0027 // Nucl.Instrum.Meth.B260:20-27, 2007
0028 //
0029 // Based on purging magnet advanced example.
0030 //
0031
0032 #include "DetectorConstruction.hh"
0033
0034 #include "G4PhysicalConstants.hh"
0035 #include "G4SystemOfUnits.hh"
0036 #include "G4NistManager.hh"
0037 #include "G4RunManager.hh"
0038
0039 // Field
0040 #include "G4Mag_UsualEqRhs.hh"
0041 #include "G4TransportationManager.hh"
0042 #include "G4ClassicalRK4.hh"
0043 #include "G4PropagatorInField.hh"
0044
0045 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0046
0047 G4ThreadLocal TabulatedField3D* DetectorConstruction::fField = 0;
0048
0049 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0050
0051 DetectorConstruction::DetectorConstruction()
0052 {
0053  fDetectorMessenger = new DetectorMessenger(this);
0054
0055  // Default values (square field, coef calculation, profile)
0056
0057  fModel=1;
0058  fG1=-11.964623;
0059  fG2=16.494652;
0060  fG3=9.866770;
0061  fG4=-6.244493;
0062  fCoef=0;
0063  fProfile=1;
0064  fGrid=0;
0065
0066 }
0067
0068 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0069
0070 DetectorConstruction::~DetectorConstruction()
0071 { delete fDetectorMessenger;}
0072
0073 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0074
0075 G4VPhysicalVolume* DetectorConstruction::Construct()
0076
0077 {
0078   DefineMaterials();
0079   return ConstructVolumes();
0080 }
0081
0082 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0083
0084 void DetectorConstruction::DefineMaterials()
0085 {
0086   G4String name, symbol;
0087   G4double density;
0088
0089   G4double z, a;
0090
0091   // Vacuum standard definition...
0092   density = universe_mean_density;
0093   G4Material* vacuum = new G4Material(name="Vacuum", z=1., a=1.01*g/mole,
0094     density);
0095
0096   // NIST
0097   G4NistManager *man=G4NistManager::Instance();
0098   man->SetVerbose(1);
0099
0100   //
0101
0102   G4cout << G4endl << *(G4Material::GetMaterialTable()) << G4endl;
0103
0104   // Default materials in setup.
0105   fDefaultMaterial = vacuum;
0106   fGridMaterial = man->FindOrBuildMaterial("G4_Ni");
0107 }
0108
0109 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0110
0111 G4VPhysicalVolume* DetectorConstruction::ConstructVolumes()
0112 {
0113
0114   fSolidWorld = new G4Box("World",          //its name
0115                12*m/2,12*m/2,22*m/2);   //its size
0116
0117
0118   fLogicWorld = new G4LogicalVolume(fSolidWorld,    //its solid
0119                     fDefaultMaterial,   //its material
0120                     "World");       //its name
0121
0122   fPhysiWorld = new G4PVPlacement(0,            //no rotation
0123                  G4ThreeVector(),   //at (0,0,0)
0124                                  "World",       //its name
0125                                  fLogicWorld,       //its logical volume
0126                                  NULL,          //its mother  volume
0127                                  false,         //no boolean operation
0128                                  0);            //copy number
0129
0130
0131   // MAGNET VOLUME
0132
0133   fSolidVol = new G4Box("Vol",              //its name
0134                10*m/2,10*m/2,9.120*m/2);    //its size
0135
0136
0137   fLogicVol = new G4LogicalVolume(fSolidVol,            //its solid
0138                   fDefaultMaterial, //its material
0139                   "Vol");       //its name
0140
0141   fPhysiVol = new G4PVPlacement(0,          //no rotation
0142                  G4ThreeVector(0,0,-4310*mm),   //at (0,0,0)
0143                                  "Vol",         //its name
0144                                  fLogicVol,     //its logical volume
0145                                  fPhysiWorld,       //its mother  volume
0146                                  false,         //no boolean operation
0147                                  0);            //copy number
0148
0149   // GRID
0150
0151   if (fGrid==1)
0152   {
0153
0154   G4cout << G4endl;
0155
0156   G4cout << " ********************** " << G4endl;
0157   G4cout << " **** GRID IN PLACE *** " << G4endl;
0158   G4cout << " ********************** " << G4endl;
0159
0160   G4double x_grid=5.0*mm;
0161   G4double y_grid=5.0*mm;
0162   G4double grid_Zpos=(250+200)*mm;      // 250+10 mm for object size of 50µm diam
0163
0164   //G4double thickness_grid=10*micrometer;
0165   G4double thickness_grid=100*micrometer;
0166
0167   G4double z_grid=thickness_grid/2.0;
0168
0169   fSolidGridVol= new G4Box("GridVolume",x_grid,y_grid,z_grid);   //its size
0170
0171   fLogicGridVol = new G4LogicalVolume(fSolidGridVol,        //its solid
0172                       fGridMaterial,            //its material
0173                       "GridVolume");        //its name
0174
0175   fPhysiGridVol = new G4PVPlacement(0,              //no rotation
0176                  G4ThreeVector(0,0,grid_Zpos),  // origin
0177                                  fLogicGridVol,         //its logical volume
0178                                  "GridVolume",          //its name
0179                                  fLogicWorld,               //its mother  volume
0180                                  false,             //no boolean operation
0181                                  0);
0182
0183   // Holes in grid
0184
0185   G4double holeSize= 9e-3*mm;
0186   G4double pix_grid=1.3e-2*mm;
0187   G4int    num_half_grid=100;
0188
0189   fSolidGridVol_Hole= new G4Box("GridHole",holeSize/2,holeSize/2,z_grid);   //its size
0190
0191   fLogicGridVol_Hole = new G4LogicalVolume(fSolidGridVol_Hole,          //its solid
0192                    fDefaultMaterial,                        //its material
0193                    "GridHole");                         //its name
0194
0195
0196   for(int i=-num_half_grid;i<num_half_grid;i++)
0197   {
0198         for (int j=-num_half_grid;j<num_half_grid;j++)
0199     {
0200
0201             G4double  x0_grid,y0_grid,z0_grid;
0202             G4int  number_index_grid;
0203
0204             x0_grid=pix_grid*i;
0205             y0_grid=pix_grid*j;
0206             z0_grid=0.0*mm;
0207
0208         number_index_grid=(i+num_half_grid)*1000+(j+num_half_grid);
0209
0210         fPhysiGridVol_Hole  = new G4PVPlacement(0,      //no rotation
0211                  G4ThreeVector(x0_grid,y0_grid,z0_grid),//origin
0212                                  fLogicGridVol_Hole,            //its logical volume
0213                      "GridHole",                //its name
0214                                  fLogicGridVol,                 //its mother  volume
0215                                  false,                 //no boolean operation
0216                                  number_index_grid);
0217     }
0218   }
0219
0220   // Grid imaging plane
0221
0222   G4double ContVolSizeXY = 1*m;
0223   G4double ImPlaneWidth = 0.001*mm;
0224
0225   fSolidControlVol_GridShadow =
0226     new G4Box
0227     ("ControlVol_GridShadow", ContVolSizeXY/2, ContVolSizeXY/2 , ImPlaneWidth/2);
0228
0229   fLogicControlVol_GridShadow =
0230     new G4LogicalVolume
0231     (fSolidControlVol_GridShadow, fDefaultMaterial, "ControlVol_GridShadow");
0232
0233   fPhysiControlVol_GridShadow =
0234     new G4PVPlacement
0235     ( 0, G4ThreeVector(0,0,(250+300)*mm), fLogicControlVol_GridShadow, "ControlVol_GridShadow",
0236       fLogicWorld, false, 0);
0237
0238
0239   } // end GRID
0240
0241   // STEP MINIMUM SIZE
0242   fLogicVol->SetUserLimits(new G4UserLimits(1*mm));
0243
0244   return fPhysiWorld;
0245 }
0246
0247 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0248
0249 void DetectorConstruction::SetG1(G4float value)
0250 {
0251   fG1 = value;
0252   G4RunManager::GetRunManager()->ReinitializeGeometry();
0253 }
0254
0255 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0256
0257 void DetectorConstruction::SetG2(G4float value)
0258 {
0259   fG2 = value;
0260   G4RunManager::GetRunManager()->ReinitializeGeometry();
0261 }
0262
0263 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0264
0265 void DetectorConstruction::SetG3(G4float value)
0266 {
0267   fG3 = value;
0268   G4RunManager::GetRunManager()->ReinitializeGeometry();
0269 }
0270
0271 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0272
0273 void DetectorConstruction::SetG4(G4float value)
0274 {
0275   fG4 = value;
0276   G4RunManager::GetRunManager()->ReinitializeGeometry();
0277 }
0278
0279 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0280
0281 void DetectorConstruction::SetModel(G4int modelChoice)
0282 {
0283   if (modelChoice==1) fModel=1;
0284   if (modelChoice==2) fModel=2;
0285   if (modelChoice==3) fModel=3;
0286   G4RunManager::GetRunManager()->ReinitializeGeometry();
0287 }
0288
0289 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0290
0291 void DetectorConstruction::SetCoef(G4int val)
0292 {
0293   fCoef=val;
0294   G4RunManager::GetRunManager()->ReinitializeGeometry();
0295 }
0296
0297 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0298
0299 G4int DetectorConstruction::GetCoef()
0300 {
0301   return fCoef;
0302 }
0303
0304 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0305
0306 void DetectorConstruction::SetProfile(G4int myProfile)
0307 {
0308   fProfile=myProfile;
0309   G4RunManager::GetRunManager()->ReinitializeGeometry();
0310 }
0311
0312 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0313
0314 void DetectorConstruction::SetGrid(G4int myGrid)
0315 {
0316   fGrid=myGrid;
0317   G4RunManager::GetRunManager()->ReinitializeGeometry();
0318 }
0319
0320 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0321
0322 void DetectorConstruction::ConstructSDandField()
0323 {
0324       fField = new TabulatedField3D(fG1, fG2, fG3, fG4, fModel);
0325
0326       //This is thread-local
0327       G4FieldManager* fFieldMgr =
0328     G4TransportationManager::GetTransportationManager()->GetFieldManager();
0329
0330       G4Mag_UsualEqRhs* fEquation = new G4Mag_UsualEqRhs (fField);
0331
0332       G4ClassicalRK4* fStepper = new G4ClassicalRK4 (fEquation);
0333
0334       G4ChordFinder* fChordFinder = new G4ChordFinder(fField,1e-9*m,fStepper);
0335
0336       fFieldMgr->SetChordFinder(fChordFinder);
0337       fFieldMgr->SetDetectorField(fField);
0338
0339       // SI: 01-07-2018 : following settings were initially set to 1e-9*m
0340       //  instead of 1e-9*m, but they now induce warnings as
0341       //  *** G4Exception : GeomNav1002
0342       //  issued by : G4PropagatorInField::ComputeStep
0343
0344       fFieldMgr->GetChordFinder()->SetDeltaChord(1e-7*m);
0345       fFieldMgr->SetDeltaIntersection(1e-7*m);
0346       fFieldMgr->SetDeltaOneStep(1e-7*m);
0347
0348       //
0349
0350       // To avoid G4MagIntegratorDriver::OneGoodStep:Stepsize underflows in Stepper
0351
0352       if (fCoef==1)
0353       {
0354         G4PropagatorInField* fPropInField =
0355           G4TransportationManager::GetTransportationManager()->GetPropagatorInField();
0356         fPropInField->SetMinimumEpsilonStep(1e-11);
0357         fPropInField->SetMaximumEpsilonStep(1e-10);
0358
0359       }
0360       else
0361       {
0362         G4PropagatorInField* fPropInField =
0363           G4TransportationManager::GetTransportationManager()->GetPropagatorInField();
0364         fPropInField->SetMinimumEpsilonStep(1e-9);
0365         fPropInField->SetMaximumEpsilonStep(1e-8);
0366       }
0367
0368 }
0369