Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-03-30 07:50:54

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 /// \file DetectorConstruction.cc
0027 /// \brief Implementation of the DetectorConstruction class
0028 
0029 // This example is provided by the Geant4-DNA collaboration
0030 // Any report or published results obtained using the Geant4-DNA software
0031 // and the DNA geometry given in the Geom_DNA example
0032 // shall cite the following Geant4-DNA collaboration publications:
0033 // [1] NIM B 298 (2013) 47-54
0034 // [2] Med. Phys. 37 (2010) 4692-4708
0035 // The Geant4-DNA web site is available at http://geant4-dna.org
0036 //
0037 
0038 #include "DetectorConstruction.hh"
0039 
0040 // Geant4
0041 #include "CLHEP/Units/SystemOfUnits.h"
0042 #include "ChromosomeParameterisation.hh"
0043 
0044 #include "G4Box.hh"
0045 #include "G4Ellipsoid.hh"
0046 #include "G4LogicalVolume.hh"
0047 #include "G4NistManager.hh"
0048 #include "G4Orb.hh"
0049 #include "G4PVParameterised.hh"
0050 #include "G4PVPlacement.hh"
0051 #include "G4RotationMatrix.hh"
0052 #include "G4Tubs.hh"
0053 #include "G4UnionSolid.hh"
0054 #include "G4VisAttributes.hh"
0055 #include "globals.hh"
0056 
0057 #define countof(x) (sizeof(x) / sizeof(x[0]))
0058 
0059 using namespace std;
0060 using CLHEP::degree;
0061 using CLHEP::micrometer;
0062 using CLHEP::mm;
0063 using CLHEP::nanometer;
0064 
0065 static G4VisAttributes visInvBlue(false, G4Colour(0.0, 0.0, 1.0));
0066 static G4VisAttributes visInvWhite(false, G4Colour(1.0, 1.0, 1.0));
0067 static G4VisAttributes visInvPink(false, G4Colour(1.0, 0.0, 1.0));
0068 static G4VisAttributes visInvCyan(false, G4Colour(0.0, 1.0, 1.0));
0069 static G4VisAttributes visInvRed(false, G4Colour(1.0, 0.0, 0.0));
0070 static G4VisAttributes visInvGreen(false, G4Colour(0.0, 1.0, 0.0));
0071 static G4VisAttributes visBlue(true, G4Colour(0.0, 0.0, 1.0));
0072 static G4VisAttributes visWhite(true, G4Colour(1.0, 1.0, 1.0));
0073 static G4VisAttributes visPink(true, G4Colour(1.0, 0.0, 1.0));
0074 static G4VisAttributes visCyan(true, G4Colour(0.0, 1.0, 1.0));
0075 static G4VisAttributes visRed(true, G4Colour(1.0, 0.0, 0.0));
0076 static G4VisAttributes visGreen(true, G4Colour(0.0, 1.0, 0.0));
0077 
0078 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0079 
0080 DetectorConstruction::DetectorConstruction()
0081   : G4VUserDetectorConstruction(), fBuildChromatineFiber(true), fBuildBases(false)
0082 {}
0083 
0084 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0085 
0086 DetectorConstruction::~DetectorConstruction() {}
0087 
0088 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0089 
0090 G4VPhysicalVolume* DetectorConstruction::Construct()
0091 {
0092   DefineMaterials();
0093   return ConstructDetector();
0094 }
0095 
0096 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0097 
0098 void DetectorConstruction::DefineMaterials()
0099 {
0100   // Water is defined from NIST material database
0101   G4NistManager* man = G4NistManager::Instance();
0102   man->FindOrBuildMaterial("G4_WATER");
0103 }
0104 
0105 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0106 
0107 void DetectorConstruction::LoadChromosome(const char* filename, G4VPhysicalVolume* chromBox,
0108                                           G4LogicalVolume* logicBoxros)
0109 {
0110   ChromosomeParameterisation* cp = new ChromosomeParameterisation(filename);
0111   new G4PVParameterised("box ros", logicBoxros, chromBox, kUndefined, cp->GetNumRosettes(), cp);
0112 
0113   G4cout << filename << " done" << G4endl;
0114 }
0115 
0116 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0117 
0118 G4VPhysicalVolume* DetectorConstruction::ConstructDetector()
0119 {
0120   if (fBuildBases == false && fBuildChromatineFiber == false) {
0121     G4cout << "======================================================" << G4endl;
0122     G4cout << "WARNING from DetectorConstruction::ConstructDetector:" << G4endl;
0123     G4cout << "As long as the flags fBuildBases and fBuildChromatineFiber are "
0124               "false, the output root file will be empty"
0125            << G4endl;
0126     G4cout << "This is intended for fast computation, display, testing ..." << G4endl;
0127     G4cout << "======================================================" << G4endl;
0128   }
0129 
0130   G4String name;
0131 
0132   /***************************************************************************/
0133   //                               World
0134   /***************************************************************************/
0135 
0136   DefineMaterials();
0137   G4Material* waterMaterial = G4Material::GetMaterial("G4_WATER");
0138 
0139   G4Box* solidWorld = new G4Box("world", 10.0 * mm, 10.0 * mm, 10.0 * mm);
0140   G4LogicalVolume* logicWorld = new G4LogicalVolume(solidWorld, waterMaterial, "world");
0141   G4PVPlacement* physiWorld =
0142     new G4PVPlacement(0, G4ThreeVector(), "world", logicWorld, 0, false, 0);
0143   logicWorld->SetVisAttributes(&visInvWhite);
0144 
0145   /****************************************************************************/
0146   //                             Box nucleus
0147   /****************************************************************************/
0148 
0149   G4Box* solidTin = new G4Box("tin", 13 * micrometer, 10 * micrometer, 5 * micrometer);
0150   G4LogicalVolume* logicTin = new G4LogicalVolume(solidTin, waterMaterial, "tin");
0151   G4VPhysicalVolume* physiTin =
0152     new G4PVPlacement(0, G4ThreeVector(), "tin", logicTin, physiWorld, false, 0);
0153   logicTin->SetVisAttributes(&visInvWhite);
0154 
0155   /****************************************************************************/
0156   //                            Cell nucleus
0157   /****************************************************************************/
0158 
0159   G4Ellipsoid* solidNucleus =
0160     new G4Ellipsoid("nucleus", 11.82 * micrometer, 8.52 * micrometer, 3 * micrometer, 0, 0);
0161   G4LogicalVolume* logicNucleus = new G4LogicalVolume(solidNucleus, waterMaterial, "logic nucleus");
0162   G4VPhysicalVolume* physiNucleus =
0163     new G4PVPlacement(0, G4ThreeVector(), "physi nucleus", logicNucleus, physiTin, false, 0);
0164   logicNucleus->SetVisAttributes(&visPink);
0165 
0166   /****************************************************************************/
0167   //                        Chromosomes territories
0168   /****************************************************************************/
0169   // NOTE: The only supported values for the rotation are
0170   // 0 and 90 degrees on the Y axis.
0171   G4double chromosomePositionSizeRotation[][7] = {
0172     {4.467, 2.835, 0, 1.557, 1.557, 1.557, 90},
0173     {-4.467, 2.835, 0, 1.557, 1.557, 1.557, 0},
0174     {4.423, -2.831, 0, 1.553, 1.553, 1.553, 90},
0175     {-4.423, -2.831, 0, 1.553, 1.553, 1.553, 0},
0176     {1.455, 5.63, 0, 1.455, 1.455, 1.455, 0},
0177     {-1.455, 5.63, 0, 1.455, 1.455, 1.455, 90},
0178     {1.435, 0, 1.392, 1.435, 1.435, 1.435, 0},
0179     {-1.435, 0, 1.392, 1.435, 1.435, 1.435, 90},
0180     {1.407, 0, -1.450, 1.407, 1.407, 1.407, 90},  // 5 right
0181     {-1.407, 0, -1.450, 1.407, 1.407, 1.407, 0},  // 5 left
0182     {1.380, -5.437, 0, 1.380, 1.380, 1.380, 0},
0183     {-1.380, -5.437, 0, 1.380, 1.380, 1.380, 90},
0184     {1.347, 2.782, -1.150, 1.347, 1.347, 1.347, 90},
0185     {-1.347, 2.782, -1.150, 1.347, 1.347, 1.347, 0},
0186     {1.311, -2.746, -1.220, 1.311, 1.311, 1.311, 90},
0187     {-1.311, -2.746, -1.220, 1.311, 1.311, 1.311, 0},
0188     {7.251, -2.541, 0, 1.275, 1.275, 1.275, 0},
0189     {-6.701, 0, -0.85, 1.275, 1.275, 1.275, 90},
0190     {4.148, 0, 1.278, 1.278, 1.278, 1.278, 90},  // 10 right
0191     {-4.148, 0, 1.278, 1.278, 1.278, 1.278, 0},  // 10 left
0192     {4.147, 0, -1.277, 1.277, 1.277, 1.277, 0},
0193     {-4.147, 0, -1.277, 1.277, 1.277, 1.277, 90},
0194     {8.930, 0.006, 0, 1.272, 1.272, 1.272, 90},
0195     {-7.296, 2.547, 0, 1.272, 1.272, 1.272, 90},
0196     {1.207, -2.642, 1.298, 1.207, 1.207, 1.207, 0},
0197     {-1.207, -2.642, 1.298, 1.207, 1.207, 1.207, 90},
0198     {1.176, 2.611, 1.368, 1.176, 1.176, 1.176, 0},
0199     {-1.176, 2.611, 1.368, 1.176, 1.176, 1.176, 90},
0200     {4.065, 5.547, 0, 1.155, 1.155, 1.155, 90},  // 15 right
0201     {-4.065, 5.547, 0, 1.155, 1.155, 1.155, 0},  // 15 left
0202     {6.542, 0.159, 1.116, 1.116, 1.116, 1.116, 0},
0203     {-9.092, 0, 0, 1.116, 1.116, 1.116, 0},
0204     {6.507, 0.159, -1.081, 1.081, 1.081, 1.081, 90},
0205     {-7.057, -2.356, 0, 1.081, 1.081, 1.081, 90},
0206     {3.824, -5.448, 0, 1.064, 1.064, 1.064, 90},
0207     {-3.824, -5.448, 0, 1.064, 1.064, 1.064, 0},
0208     {5.883, -5.379, 0, 0.995, 0.995, 0.995, 0},
0209     {-9.133, -2.111, 0, 0.995, 0.995, 0.995, 0},
0210     {6.215, 5.387, 0, 0.995, 0.995, 0.995, 0},  // 20 right
0211     {-6.971, -4.432, 0, 0.995, 0.995, 0.995, 90},  // 20 left
0212     {9.583, 2.177, 0, 0.899, 0.899, 0.899, 90},
0213     {-9.467, 2.03, 0, 0.899, 0.899, 0.899, 0},
0214     {9.440, -2.180, 0, 0.914, 0.914, 0.914, 90},
0215     {-6.34, 0, 1.339, 0.914, 0.914, 0.914, 0},
0216     {-6.947, 4.742, 0, 0.923, 0.923, 0.923, 90},  // Y
0217     {7.354, 2.605, 0, 1.330, 1.330, 1.330, 0}  // X
0218   };
0219 
0220   G4RotationMatrix* rotch = new G4RotationMatrix;
0221   rotch->rotateY(90 * degree);
0222 
0223   vector<G4VPhysicalVolume*> physiBox(48);
0224 
0225   for (unsigned int i = 0; i < countof(chromosomePositionSizeRotation); i++) {
0226     G4double* p = &chromosomePositionSizeRotation[i][0];
0227     G4double* size = &chromosomePositionSizeRotation[i][3];
0228     G4double rotation = chromosomePositionSizeRotation[i][6];
0229     G4ThreeVector pos(p[0] * micrometer, p[1] * micrometer, p[2] * micrometer);
0230     G4RotationMatrix* rot = rotation == 0 ? 0 : rotch;
0231 
0232     ostringstream ss;
0233     ss << "box" << (i / 2) + 1 << (i % 2 ? 'l' : 'r');
0234     name = ss.str();
0235     ss.str("");
0236     ss.clear();
0237 
0238     /*
0239      snprintf(name, countof(name), "box%d%c",
0240      (i / 2) + 1, i % 2 ? 'l' : 'r');
0241      */
0242     G4Box* solidBox =
0243       new G4Box(name, size[0] * micrometer, size[1] * micrometer, size[2] * micrometer);
0244     G4LogicalVolume* logicBox = new G4LogicalVolume(solidBox, waterMaterial, name);
0245     physiBox[i] = new G4PVPlacement(rot, pos, "chromo", logicBox, physiNucleus, false, 0);
0246     logicBox->SetVisAttributes(&visBlue);
0247   }
0248 
0249   /**************************************************************************/
0250   //                 Box containing the chromatin flowers
0251   /**************************************************************************/
0252 
0253   G4Tubs* solidBoxros = new G4Tubs("solid box ros", 0 * nanometer, 399 * nanometer, 20 * nanometer,
0254                                    0 * degree, 360 * degree);
0255   G4LogicalVolume* logicBoxros = new G4LogicalVolume(solidBoxros, waterMaterial, "box ros");
0256   logicBoxros->SetVisAttributes(&visInvBlue);
0257 
0258   // Loading flower box position for each chromosome territory
0259 
0260   for (int k = 0; k < 22; k++) {
0261     ostringstream oss;
0262     oss << "chromo" << k + 1 << ".dat";
0263     name = oss.str();
0264     oss.str("");
0265     oss.clear();
0266     // snprintf(name, countof(name), "chromo%d.dat", k + 1);
0267     LoadChromosome(name.c_str(), physiBox[k * 2], logicBoxros);
0268     LoadChromosome(name.c_str(), physiBox[k * 2 + 1], logicBoxros);
0269   }
0270 
0271   LoadChromosome("chromoY.dat", physiBox[44], logicBoxros);
0272   LoadChromosome("chromoX.dat", physiBox[45], logicBoxros);
0273 
0274   /****************************************************************************/
0275   if (fBuildChromatineFiber) {
0276     // chromatin fiber envelope
0277     G4Tubs* solidEnv = new G4Tubs("chromatin fiber", 0, 15.4 * nanometer, 80.5 * nanometer,
0278                                   0 * degree, 360 * degree);
0279     G4LogicalVolume* logicEnv = new G4LogicalVolume(solidEnv, waterMaterial, "LV chromatin fiber");
0280     logicEnv->SetVisAttributes(&visInvPink);
0281 
0282     // Chromatin fiber position
0283     for (G4int i = 0; i < 7; i++) {
0284       G4RotationMatrix* rotFiber = new G4RotationMatrix;
0285       rotFiber->rotateX(90 * degree);
0286       rotFiber->rotateY(i * 25.72 * degree);
0287       G4ThreeVector posFiber = G4ThreeVector(0, 152 * nanometer, 0);
0288       posFiber.rotateZ(i * 25.72 * degree);
0289       new G4PVPlacement(rotFiber, posFiber, logicEnv, "physi env", logicBoxros, false, 0);
0290 
0291       rotFiber = new G4RotationMatrix;
0292       rotFiber->rotateX(90 * degree);
0293       rotFiber->rotateY((7 + i) * 25.72 * degree);
0294       posFiber = G4ThreeVector(0, 152 * nanometer, 0);
0295       posFiber.rotateZ((7 + i) * 25.72 * degree);
0296       new G4PVPlacement(rotFiber, posFiber, logicEnv, "physi env", logicBoxros, false, 0);
0297 
0298       rotFiber = new G4RotationMatrix;
0299       rotFiber->rotateX(90 * degree);
0300       rotFiber->rotateY((25.72 + (i - 14) * 51.43) * degree);
0301       posFiber = G4ThreeVector(-36.5 * nanometer, 312 * nanometer, 0);
0302       posFiber.rotateZ((i - 14) * 51.43 * degree);
0303       new G4PVPlacement(rotFiber, posFiber, logicEnv, "physi env", logicBoxros, false, 0);
0304 
0305       rotFiber = new G4RotationMatrix;
0306       rotFiber->rotateX(90 * degree);
0307       rotFiber->rotateY(180 * degree);
0308       rotFiber->rotateY((i - 21) * 51.43 * degree);
0309       posFiber = G4ThreeVector(-103 * nanometer, 297 * nanometer, 0);
0310       posFiber.rotateZ((i - 21) * 51.43 * degree);
0311       new G4PVPlacement(rotFiber, posFiber, logicEnv, "physi env", logicBoxros, false, 0);
0312     }
0313 
0314     if (fBuildBases) {
0315       // Histones
0316       G4Tubs* solidHistone = new G4Tubs("solid histone", 0, 3.25 * nanometer, 2.85 * nanometer,
0317                                         0 * degree, 360 * degree);
0318       G4LogicalVolume* logicHistone =
0319         new G4LogicalVolume(solidHistone, waterMaterial, "logic histone");
0320 
0321       // Base pair
0322       G4Orb* solidBp1 = new G4Orb("blue sphere", 0.17 * nanometer);
0323       G4LogicalVolume* logicBp1 = new G4LogicalVolume(solidBp1, waterMaterial, "logic blue sphere");
0324       G4Orb* solidBp2 = new G4Orb("pink sphere", 0.17 * nanometer);
0325       G4LogicalVolume* logicBp2 = new G4LogicalVolume(solidBp2, waterMaterial, "logic pink sphere");
0326 
0327       // Phosphodiester group
0328 
0329       G4Orb* solidSugar_48em1_nm = new G4Orb("sugar", 0.48 * nanometer);
0330 
0331       G4ThreeVector posi(0.180248 * nanometer, 0.32422 * nanometer, 0.00784 * nanometer);
0332       G4UnionSolid* uniDNA =
0333         new G4UnionSolid("move", solidSugar_48em1_nm, solidSugar_48em1_nm, 0, posi);
0334 
0335       G4ThreeVector posi2(-0.128248 * nanometer, 0.41227 * nanometer, 0.03584 * nanometer);
0336       G4UnionSolid* uniDNA2 =
0337         new G4UnionSolid("move2", solidSugar_48em1_nm, solidSugar_48em1_nm, 0, posi2);
0338 
0339       /************************************************************************
0340        Phosphodiester group Position
0341        ************************************************************************/
0342 
0343       for (G4int n = 2; n < 200; n++) {
0344         G4double SP1[2][3] = {
0345           {(-0.6 * nanometer) * cos(n * 0.26), 0, (0.6 * nanometer) * sin(n * 0.26)},
0346           {(0.6 * nanometer) * cos(n * 0.26), 0, (-0.6 * nanometer) * sin(0.26 * n)}};
0347         G4double matriceSP1[3][3] = {
0348           {cos(n * 0.076), -sin(n * 0.076), 0}, {sin(n * 0.076), cos(n * 0.076), 0}, {0, 0, 1}};
0349         G4double matriceSP2[2][3];
0350 
0351         for (G4int i = 0; i < 3; i++) {
0352           G4double sumSP1 = 0;
0353           G4double sumSP2 = 0;
0354           for (G4int j = 0; j < 3; j++) {
0355             sumSP1 += matriceSP1[i][j] * SP1[0][j];
0356             sumSP2 += matriceSP1[i][j] * SP1[1][j];
0357           }
0358           matriceSP2[0][i] = sumSP1;
0359           matriceSP2[1][i] = sumSP2;
0360         }
0361 
0362         G4double heliceSP[3] = {(4.85 * nanometer) * cos(n * 0.076),
0363                                 (4.85 * nanometer) * sin(n * 0.076), (n * 0.026 * nanometer)};
0364 
0365         for (G4int i = 0; i < 3; i++) {
0366           matriceSP2[0][i] += heliceSP[i];
0367           matriceSP2[1][i] += heliceSP[i];
0368         }
0369         G4ThreeVector posSugar1(matriceSP2[0][2], matriceSP2[0][1],
0370                                 (matriceSP2[0][0]) - (4.25 * nanometer));
0371         G4ThreeVector posSugar2(matriceSP2[1][2], matriceSP2[1][1],
0372                                 (matriceSP2[1][0]) - (5.45 * nanometer));
0373 
0374         ostringstream ss;
0375         ss << "sugar_" << n;
0376         name = ss.str().c_str();
0377         ss.str("");
0378         ss.clear();
0379 
0380         //  snprintf(name, countof(name), "sugar %d", n);
0381         uniDNA = new G4UnionSolid(name, uniDNA, solidSugar_48em1_nm, 0, posSugar1);
0382 
0383         ss << "sugar_" << n;
0384         name = ss.str().c_str();
0385         ss.str("");
0386         ss.clear();
0387 
0388         //  snprintf(name, countof(name), "sugar %d", n);
0389         uniDNA2 = new G4UnionSolid(name, uniDNA2, solidSugar_48em1_nm, 0, posSugar2);
0390       }
0391       G4LogicalVolume* logicSphere3 = new G4LogicalVolume(uniDNA, waterMaterial, "logic sugar 2");
0392       G4LogicalVolume* logicSphere4 = new G4LogicalVolume(uniDNA2, waterMaterial, "logic sugar 4");
0393 
0394       /**************************************************************************
0395        Base pair Position
0396        **************************************************************************/
0397       for (G4int n = 0; n < 200; n++) {
0398         G4double bp1[2][3] = {
0399           {(-0.34 * nanometer) * cos(n * 0.26), 0, (0.34 * nanometer) * sin(n * 0.26)},
0400           {(0.34 * nanometer) * cos(n * 0.26), 0, (-0.34 * nanometer) * sin(0.26 * n)}};
0401         G4double matriceBP1[3][3] = {
0402           {cos(n * 0.076), -sin(n * 0.076), 0}, {sin(n * 0.076), cos(n * 0.076), 0}, {0, 0, 1}};
0403         G4double matriceBP2[2][3];
0404 
0405         for (G4int i = 0; i < 3; i++) {
0406           G4double sumBP1 = 0;
0407           G4double sumBP2 = 0;
0408           for (G4int j = 0; j < 3; j++) {
0409             sumBP1 += matriceBP1[i][j] * bp1[0][j];
0410             sumBP2 += matriceBP1[i][j] * bp1[1][j];
0411           }
0412           matriceBP2[0][i] = sumBP1;
0413           matriceBP2[1][i] = sumBP2;
0414         }
0415         G4double heliceBP[3] = {(4.8 * nanometer) * cos(n * 0.076),
0416                                 (4.8 * nanometer) * sin(n * 0.076), n * 0.026 * nanometer};
0417 
0418         for (G4int i = 0; i < 3; i++) {
0419           matriceBP2[0][i] += heliceBP[i];
0420           matriceBP2[1][i] += heliceBP[i];
0421         }
0422         G4ThreeVector position1(matriceBP2[0][2], matriceBP2[0][1],
0423                                 matriceBP2[0][0] - (4.25 * nanometer));
0424         G4ThreeVector position2(matriceBP2[1][2], matriceBP2[1][1],
0425                                 matriceBP2[1][0] - (5.45 * nanometer));
0426 
0427         new G4PVPlacement(0, position1, logicBp1, "physi blue sphere", logicSphere3, false, 0);
0428         new G4PVPlacement(0, position2, logicBp2, "physi pink sphere", logicSphere4, false, 0);
0429       }
0430 
0431       /****************************************************************************/
0432       //                 Initial position of different elements
0433       /****************************************************************************/
0434       // DNA and histone positions
0435       for (int j = 0; j < 90; j++) {
0436         // DNA (bp-SP)
0437         G4RotationMatrix* rotStrand1 = new G4RotationMatrix;
0438         rotStrand1->rotateZ(j * -51.43 * degree);
0439         G4ThreeVector posStrand1(-2.7 * nanometer, 9.35 * nanometer,
0440                                  (-69.9 * nanometer) + (j * 1.67 * nanometer));
0441         posStrand1.rotateZ(j * 51.43 * degree);
0442         new G4PVPlacement(rotStrand1, posStrand1, logicSphere3, "physi sugar 2", logicEnv, false,
0443                           0);
0444 
0445         G4RotationMatrix* rotStrand2 = new G4RotationMatrix;
0446         rotStrand2->rotateZ(j * -51.43 * degree);
0447         G4ThreeVector posStrand2(-2.7 * nanometer, 9.35 * nanometer,
0448                                  (-68.7 * nanometer) + (j * 1.67 * nanometer));
0449         posStrand2.rotateZ(j * 51.43 * degree);
0450         new G4PVPlacement(rotStrand2, posStrand2, logicSphere4, "physi sugar 4", logicEnv, false,
0451                           0);
0452 
0453         // histones
0454         G4RotationMatrix* rotHistone = new G4RotationMatrix;
0455         rotHistone->rotateY(90 * degree);
0456         rotHistone->rotateX(j * (-51.43 * degree));
0457         G4ThreeVector posHistone(0.0, 9.35 * nanometer, (-74.15 + j * 1.67) * nanometer);
0458         posHistone.rotateZ(j * 51.43 * degree);
0459         new G4PVPlacement(rotHistone, posHistone, logicHistone, "PV histone", logicEnv, false, 0);
0460       }
0461       /************************************************************************/
0462       //                        Visualisation colors
0463       /************************************************************************/
0464 
0465       logicBp1->SetVisAttributes(&visInvCyan);
0466       logicBp2->SetVisAttributes(&visInvPink);
0467 
0468       logicSphere3->SetVisAttributes(&visInvWhite);
0469       logicSphere4->SetVisAttributes(&visInvRed);
0470 
0471       logicHistone->SetVisAttributes(&visInvBlue);
0472     }
0473   }
0474 
0475   G4cout << "Geometry has been loaded" << G4endl;
0476   return physiWorld;
0477 }