| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // This example is provided by the Geant4-DNA 26 // This example is provided by the Geant4-DNA collaboration
27 // Any report or published results obtained us 27 // Any report or published results obtained using the Geant4-DNA software
28 // shall cite the following Geant4-DNA collabo 28 // shall cite the following Geant4-DNA collaboration publication:
29 // Med. Phys. 37 (2010) 4692-4708 29 // Med. Phys. 37 (2010) 4692-4708
30 // and papers 30 // and papers
31 // M. Batmunkh et al. J Radiat Res Appl Sci 8 31 // M. Batmunkh et al. J Radiat Res Appl Sci 8 (2015) 498-507
32 // O. Belov et al. Physica Medica 32 (2016) 15 << 32 // O. Belov et al. Physica Medica 32 (2016) 1510-1520
33 // The Geant4-DNA web site is available at htt 33 // The Geant4-DNA web site is available at http://geant4-dna.org
34 // << 34 //
35 // ------------------------------------------- 35 // -------------------------------------------------------------------
36 // November 2016 36 // November 2016
37 // ------------------------------------------- 37 // -------------------------------------------------------------------
38 // 38 //
39 /// \file DetectorConstruction.cc << 39 /// \file DetectorConstruction.cc
40 /// \brief Implementation of the DetectorConst 40 /// \brief Implementation of the DetectorConstruction class
41 41
42 #include "DetectorConstruction.hh" 42 #include "DetectorConstruction.hh"
43 <<
44 #include "G4Colour.hh" <<
45 #include "G4PhysicalConstants.hh" <<
46 #include "G4ProductionCuts.hh" <<
47 #include "G4Region.hh" <<
48 #include "G4RotationMatrix.hh" <<
49 #include "G4SystemOfUnits.hh" 43 #include "G4SystemOfUnits.hh"
>> 44 #include "G4Region.hh"
>> 45 #include "G4ProductionCuts.hh"
>> 46 #include "G4PhysicalConstants.hh"
>> 47 #include "G4Colour.hh"
50 #include "G4VisAttributes.hh" 48 #include "G4VisAttributes.hh"
51 << 49 #include "G4RotationMatrix.hh"
52 #include <algorithm> << 50 #include <algorithm>
53 #include <iostream> 51 #include <iostream>
54 52
55 //....oooOO0OOooo........oooOO0OOooo........oo 53 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
56 54
57 DetectorConstruction::DetectorConstruction() 55 DetectorConstruction::DetectorConstruction()
>> 56 :G4VUserDetectorConstruction(),
>> 57 fpDefaultMaterial(0), fpWaterMaterial(0),
>> 58 fpRegion(0), fCheckOverlaps(false)
58 { 59 {
59 fNeuronLoadParamz = new NeuronLoadDataFile() << 60
60 DefineMaterials(); << 61 }
61 } <<
62 62
63 //....oooOO0OOooo........oooOO0OOooo........oo 63 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
64 64
65 DetectorConstruction::~DetectorConstruction() 65 DetectorConstruction::~DetectorConstruction()
66 { 66 {
67 delete fNeuronLoadParamz; << 67
68 } 68 }
69 69
70 //....oooOO0OOooo........oooOO0OOooo........oo 70 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
71 71
72 void DetectorConstruction::DefineMaterials() << 72 G4VPhysicalVolume* DetectorConstruction::Construct()
73 { 73 {
>> 74
>> 75 // load Neuron positions and obtain parameters!
>> 76 fNeuronLoadParamz = new NeuronLoadDataFile() ;
>> 77
>> 78 DefineMaterials();
>> 79 return ConstructDetector();
>> 80
>> 81 }
>> 82
>> 83 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 84
>> 85 void DetectorConstruction::DefineMaterials()
>> 86 {
74 // Water is defined from NIST material datab 87 // Water is defined from NIST material database
75 G4NistManager* man = G4NistManager::Instance << 88 G4NistManager * man = G4NistManager::Instance();
76 fpWaterMaterial = man->FindOrBuildMaterial(" << 89 G4Material * H2O = man->FindOrBuildMaterial("G4_WATER");
77 fpWorldMaterial = man->FindOrBuildMaterial(" << 90 // Default materials in setup.
>> 91 fpWaterMaterial = H2O;
>> 92 G4Material * Vacuum = man->FindOrBuildMaterial("G4_Galactic");
>> 93 //G4Material * Air = man->FindOrBuildMaterial("G4_AIR");
>> 94 fpDefaultMaterial = Vacuum;
78 } 95 }
79 96
80 //....oooOO0OOooo........oooOO0OOooo........oo 97 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
81 98
82 G4VPhysicalVolume* DetectorConstruction::Const << 99 G4VPhysicalVolume* DetectorConstruction::ConstructDetector()
83 { 100 {
84 G4cout << " ---- Begin of Neuron Constructio << 101 G4cout <<" ---- Begin of Neuron Construction! ---- "
85 << "\n" << 102 <<"\n"<<" =========================================================="<<G4endl;
86 << " ================================ <<
87 103
88 // ========================================= << 104 // ===============================================
89 // WORLD VOLUME - filled by default material 105 // WORLD VOLUME - filled by default material
90 // ========================================= 106 // ===============================================
91 107
92 // Dimensions of world volume are calculated << 108 // Dimensions of world volume are calculated as overall dimensions of neuron!
93 << 109
94 G4double worldSizeX; << 110 G4double worldSizeX;
95 worldSizeX = fNeuronLoadParamz->GetdiagnlLen << 111 worldSizeX = 1.0*fNeuronLoadParamz->GetdiagnlLength()*um;
96 << 112
97 if (worldSizeX <= 0.0) { << 113 if (!worldSizeX)
98 worldSizeX = 10. * cm; << 114 {
99 } << 115 worldSizeX = 10.*cm;
100 << 116 }
101 G4double worldSizeY = worldSizeX; << 117
102 G4double worldSizeZ = worldSizeX; << 118 G4double worldSizeY = worldSizeX;
103 G4cout << " Side length of word volume is ca << 119 G4double worldSizeZ = worldSizeX;
104 G4VSolid* worldS = new G4Box("World", worldS << 120 G4cout << " Side length of word volume is calculated : "
105 << 121 << worldSizeX/um <<" um"<< G4endl;
106 G4LogicalVolume* worldLV = new G4LogicalVolu << 122 G4VSolid* worldS = new G4Box("World", //its name
>> 123 worldSizeX/2, worldSizeY/2, worldSizeZ/2); //its size
>> 124
>> 125 G4LogicalVolume* worldLV = new G4LogicalVolume(worldS, //its solid
>> 126 fpDefaultMaterial, //its material
>> 127 "World"); //its name
107 128
108 // Visualization attributes 129 // Visualization attributes
109 G4VisAttributes* worldVisAtt = new G4VisAttr << 130 G4VisAttributes* worldVisAtt =
>> 131 new G4VisAttributes(G4Colour(0.5,0.5,0.5,0.1)); //Gray
>> 132 //worldVisAtt->SetForceSolid(true);
110 worldVisAtt->SetVisibility(true); 133 worldVisAtt->SetVisibility(true);
111 worldLV->SetVisAttributes(worldVisAtt); 134 worldLV->SetVisAttributes(worldVisAtt);
112 135
113 G4VPhysicalVolume* worldPV = new G4PVPlaceme << 136
114 << 137 G4VPhysicalVolume* worldPV = new G4PVPlacement(0, //no rotation
115 << 138 G4ThreeVector(), //at (0,0,0)
116 << 139 "World", //its name
117 << 140 worldLV, //its logical volume
118 << 141 0, //its mother volume
119 << 142 false, //no boolean operation
120 << 143 0, //copy number
>> 144 true); // checking overlaps forced to YES
121 145
122 // ========================================= 146 // ===============================================
123 // HOMOGENEOUS MEDIUM - LARGE SPHERE VOLUME 147 // HOMOGENEOUS MEDIUM - LARGE SPHERE VOLUME
124 // ========================================= 148 // ===============================================
125 149
126 // Radius of water sphere calculated as over 150 // Radius of water sphere calculated as overall dimensions of neuron.
127 fTotMassMedium = fNeuronLoadParamz->GetTotMa << 151 fTotMassMedium = fNeuronLoadParamz->GetTotMassMedium() ;
128 fTotSurfMedium = fNeuronLoadParamz->GetTotSu << 152 fTotSurfMedium = fNeuronLoadParamz->GetTotSurfMedium() ;
129 G4double RadiusMedium = fNeuronLoadParamz->G << 153 G4double RadiusMedium = fNeuronLoadParamz->GetdiagnlLength()*um / 2.;
130 G4cout << " Radius of homogeneous medium is << 154 G4cout << " Radius of homogeneous medium is calculated : "
131 << G4endl; << 155 << RadiusMedium/um <<" um"<< G4endl;
132 G4VSolid* mediumS = new G4Orb("Medium", Radi << 156 G4VSolid* mediumS = new G4Orb("Medium", RadiusMedium);
133 << 157 //G4VSolid* mediumS = new G4Box("Medium", RadiusMedium*um,
134 G4LogicalVolume* mediumLV = new G4LogicalVol << 158 // RadiusMedium*um,RadiusMedium*um);
135 << 159
>> 160 G4LogicalVolume* mediumLV =
>> 161 new G4LogicalVolume(mediumS, fpWaterMaterial, "Medium");
>> 162
136 // Visualization attributes 163 // Visualization attributes
137 G4VisAttributes* mediumVisAtt = new G4VisAtt << 164 G4VisAttributes* mediumVisAtt =
138 // mediumVisAtt->SetForceSolid(true); << 165 new G4VisAttributes(G4Colour(0.0,1.0,0.0,0.02)); //Green
139 // mediumVisAtt->SetForceWireframe (true); << 166 //mediumVisAtt->SetForceSolid(true);
>> 167 //mediumVisAtt->SetForceWireframe (true);
140 mediumVisAtt->SetForceLineSegmentsPerCircle( 168 mediumVisAtt->SetForceLineSegmentsPerCircle(180);
141 mediumVisAtt->SetVisibility(true); 169 mediumVisAtt->SetVisibility(true);
142 mediumLV->SetVisAttributes(mediumVisAtt); 170 mediumLV->SetVisAttributes(mediumVisAtt);
143 << 171
144 G4VPhysicalVolume* mediumPV = << 172 G4VPhysicalVolume* mediumPV = new G4PVPlacement(0,
145 new G4PVPlacement(0, G4ThreeVector(), "Med << 173 G4ThreeVector(),
>> 174 "Medium",
>> 175 mediumLV,
>> 176 worldPV,
>> 177 false,
>> 178 0,
>> 179 fCheckOverlaps);
146 180
147 // ========================================= 181 // ===============================================
148 // TARGET - BOUNDING SLICE including NEURON << 182 // TARGET - BOUNDING SLICE including NEURON
149 // ========================================= 183 // ===============================================
150 184
151 // Dimensions of bounding slice volume defin 185 // Dimensions of bounding slice volume defined as overall measure of neuron
152 << 186
153 G4double TargetSizeX = fNeuronLoadParamz->Ge << 187 G4double TargetSizeX = fNeuronLoadParamz->GetwidthB()*um;
154 G4double TargetSizeY = fNeuronLoadParamz->Ge << 188 G4double TargetSizeY = fNeuronLoadParamz->GetheightB()*um;
155 G4double TargetSizeZ = fNeuronLoadParamz->Ge << 189 G4double TargetSizeZ = fNeuronLoadParamz->GetdepthB()*um;
156 fTotMassSlice = fNeuronLoadParamz->GetTotMas << 190 fTotMassSlice = fNeuronLoadParamz->GetTotMassSlice() ;
157 G4cout << " Overall dimensions (um) of neuro << 191 G4cout << " Overall dimensions (um) of neuron morphology : " << "\n"
158 << "\n" << 192 << '\t'<< " width = " <<TargetSizeX/um<< " height = " <<
159 << '\t' << " width = " << TargetSizeX << 193 TargetSizeY/um
160 << " depth = " << TargetSizeZ / um << << 194 << " depth = " <<TargetSizeZ/um<<G4endl;
161 << 195
162 G4cout << " Volume (um3), surface (um2) and 196 G4cout << " Volume (um3), surface (um2) and mass (ug) of Bounding Slice are"
163 << " calculated : " << 197 << " calculated : " << "\n"
164 << "\n" << 198 << '\t'<<fNeuronLoadParamz->GetTotVolSlice()/std::pow(um,3)<<"; "<<'\t'
165 << '\t' << fNeuronLoadParamz->GetTotV << 199 <<fNeuronLoadParamz->GetTotSurfSlice()/(um*um)
166 << fNeuronLoadParamz->GetTotSurfSlice << 200 <<"; "<<'\t'<<fNeuronLoadParamz->GetTotMassSlice()*1e6/g<< "\n "<<G4endl;
167 << fNeuronLoadParamz->GetTotMassSlice << 201
168 << 202 G4Box* boundingS = new G4Box("BoundingSlice",
169 G4Box* boundingS = << 203 TargetSizeX/2.,TargetSizeY/2.,TargetSizeZ/2.);
170 new G4Box("BoundingSlice", TargetSizeX / 2 <<
171 204
172 G4LogicalVolume* boundingLV = new G4LogicalV << 205 G4LogicalVolume* boundingLV =
>> 206 new G4LogicalVolume(boundingS,fpWaterMaterial, "BoundingSlice");
173 207
174 // Visualization attributes with opacity! 208 // Visualization attributes with opacity!
175 G4VisAttributes* TargetVisAtt = new G4VisAtt << 209 G4VisAttributes* TargetVisAtt =
>> 210 new G4VisAttributes(G4Colour(1.0,1.0,1.0,0.1));
176 TargetVisAtt->SetForceSolid(true); 211 TargetVisAtt->SetForceSolid(true);
177 TargetVisAtt->SetVisibility(true); 212 TargetVisAtt->SetVisibility(true);
178 boundingLV->SetVisAttributes(TargetVisAtt); << 213 boundingLV->SetVisAttributes(TargetVisAtt);
179 new G4PVPlacement(0, G4ThreeVector(), "Bound << 214 //G4VPhysicalVolume* boundingPV =
180 fCheckOverlaps); << 215 new G4PVPlacement(0,
>> 216 G4ThreeVector(),
>> 217 "BoundingSlice",
>> 218 boundingLV,
>> 219 mediumPV,
>> 220 false,
>> 221 0,
>> 222 fCheckOverlaps);
181 223
182 // ========================================= 224 // ===============================================
183 // NEURON MORPHOLOGY 225 // NEURON MORPHOLOGY
184 // ========================================= 226 // ===============================================
185 227
186 G4cout << " Volume (um3), surface (um2) and << 228 G4cout << " Volume (um3), surface (um2) and mass(ug) of Neuron "
187 << "are calculated : " << 229 << "are calculated : "<< "\n"
188 << "\n" << 230 << '\t'<<fNeuronLoadParamz->GetTotVolNeuron()/std::pow(um,3)<<"; "<<'\t'
189 << '\t' << fNeuronLoadParamz->GetTotV << 231 <<fNeuronLoadParamz->GetTotSurfNeuron()/(um*um)
190 << fNeuronLoadParamz->GetTotSurfNeuro << 232 <<"; "<<'\t'<<fNeuronLoadParamz->GetTotMassNeuron()*1e6/g<<G4endl;
191 << fNeuronLoadParamz->GetTotMassNeuro << 233 fTotMassNeuron = fNeuronLoadParamz->GetTotMassNeuron() ;
192 fTotMassNeuron = fNeuronLoadParamz->GetTotMa << 234 G4cout << " Total number of compartments into Neuron : "
193 G4cout << " Total number of compartments int << 235 <<fNeuronLoadParamz->GetnbNeuroncomp()<<G4endl;
194 << G4endl; <<
195 G4cout << " Shift values (um) for Neuron tra 236 G4cout << " Shift values (um) for Neuron translation are calculated : "
196 << "\n" << 237 << "\n"
197 << '\t' << " shiftX = " << fNeuronLoa << 238 << '\t' << " shiftX = " <<fNeuronLoadParamz->GetshiftX()<< " shiftY = "
198 << " shiftY = " << fNeuronLoadParamz- << 239 << fNeuronLoadParamz->GetshiftY()
199 << " shiftZ = " << fNeuronLoadParamz- << 240 << " shiftZ = " <<fNeuronLoadParamz->GetshiftZ()<< "\n"<< G4endl;
200 << 241
201 // Soma in Violet with opacity // 0.85,0.4 << 242 // Soma in Violet with opacity // 0.85,0.44,0.84
202 fSomaColour = new G4VisAttributes; << 243 fSomaColour = new G4VisAttributes;
203 fSomaColour->SetColour(G4Colour(G4Colour(22 << 244 fSomaColour->SetColour(G4Colour(G4Colour(22/255. , 200/255. , 30/255.)));
204 fSomaColour->SetForceSolid(true); // true << 245 fSomaColour->SetForceSolid(true); // true
205 fSomaColour->SetVisibility(true); << 246 fSomaColour->SetVisibility(true);
206 << 247
207 // Dendrites in Dark-Blue << 248 // Dendrites in Dark-Blue
208 fDendColour = new G4VisAttributes; << 249 fDendColour = new G4VisAttributes;
209 fDendColour->SetColour(G4Colour(G4Colour(0.0 << 250 fDendColour->SetColour(G4Colour(G4Colour(0.0, 0.0, 0.5)));
210 fDendColour->SetForceSolid(true); << 251 fDendColour->SetForceSolid(true);
211 // fDendColour->SetVisibility(true); << 252 //fDendColour->SetVisibility(true);
212 << 253
213 // Axon in Maroon << 254 // Axon in Maroon
214 fAxonColour = new G4VisAttributes; << 255 fAxonColour = new G4VisAttributes;
215 fAxonColour->SetColour(G4Colour(G4Colour(0.5 << 256 fAxonColour->SetColour(G4Colour(G4Colour(0.5, 0.0, 0.0)));
216 fAxonColour->SetForceSolid(true); << 257 fAxonColour->SetForceSolid(true);
217 fAxonColour->SetVisibility(true); << 258 fAxonColour->SetVisibility(true);
218 << 259
219 // Spines in Dark-Green << 260 // Spines in Dark-Green
220 fSpineColour = new G4VisAttributes; << 261 fSpineColour = new G4VisAttributes;
221 fSpineColour->SetColour(G4Colour(G4Colour(0. << 262 fSpineColour->SetColour(G4Colour(G4Colour(0.0 , 100/255. , 0.0)));
222 fSpineColour->SetForceSolid(true); << 263 fSpineColour->SetForceSolid(true);
223 fSpineColour->SetVisibility(true); << 264 fSpineColour->SetVisibility(true);
224 << 265
225 // Whole neuron in semitransparent navy blue << 266 // Whole neuron in semitransparent navy blue
226 fNeuronColour = new G4VisAttributes; << 267 fNeuronColour = new G4VisAttributes;
227 fNeuronColour->SetColour(G4Colour(G4Colour(0 << 268 fNeuronColour->SetColour(G4Colour(G4Colour(0.0,0.4,0.8,0.9)));
228 fNeuronColour->SetForceSolid(true); << 269 fNeuronColour->SetForceSolid(true);
229 fNeuronColour->SetVisibility(true); << 270 fNeuronColour->SetVisibility(true);
230 << 271
231 // Placement volumes: G4examples/extended/pa << 272 // Placement volumes: G4examples/extended/parameterisations/gflash
232 273
233 // ========================================= 274 // =======================================================================
234 // Structure-1: Soma 275 // Structure-1: Soma
235 276
236 // Create Target G4Region and add logical vo 277 // Create Target G4Region and add logical volume
237 // Active Geant4-DNA processes in this regio << 278 // Active Geant4-DNA processes in this region
238 fpRegion = new G4Region("Soma"); << 279 fpRegion = new G4Region("Soma");
239 G4ProductionCuts* cuts = new G4ProductionCut 280 G4ProductionCuts* cuts = new G4ProductionCuts();
240 G4double defCut = 1 * nanometer; << 281 G4double defCut = 1*nanometer;
241 cuts->SetProductionCut(defCut, "gamma"); << 282 cuts->SetProductionCut(defCut,"gamma");
242 cuts->SetProductionCut(defCut, "e-"); << 283 cuts->SetProductionCut(defCut,"e-");
243 cuts->SetProductionCut(defCut, "e+"); << 284 cuts->SetProductionCut(defCut,"e+");
244 cuts->SetProductionCut(defCut, "proton"); << 285 cuts->SetProductionCut(defCut,"proton");
>> 286
>> 287 fnbSomacomp = fNeuronLoadParamz->GetnbSomacomp() ;
>> 288 fMassSomaTot = fNeuronLoadParamz->GetMassSomaTot() ;
>> 289 fMassSomacomp = new G4double[fnbSomacomp];
>> 290 fPosSomacomp = new G4ThreeVector[fnbSomacomp];
>> 291 if (fNeuronLoadParamz->GetnbSomacomp()==0)
>> 292 {
>> 293 G4cout <<" ---- Soma not found! ---- "<< G4endl;
>> 294 }
>> 295 else
>> 296 {
>> 297 G4cout <<" ---- Soma for construction: ---- "<< G4endl;
>> 298 G4cout << " Total number of compartments into Soma : "
>> 299 <<fNeuronLoadParamz->GetnbSomacomp()<<G4endl;
>> 300 for (G4int i=0; i<fNeuronLoadParamz->GetnbSomacomp() ; i++)
>> 301 {
>> 302 fsomaS [i] = new G4Orb("Soma", fNeuronLoadParamz->GetRadSomacomp(i)* um);
>> 303 //fsomaS [i] = new G4Ellipsoid("Soma", Ra *um, Rb *um, Rc *um, 0,0) ;
>> 304 // you can change parameters of Soma with a single compartment
>> 305 fsomaLV[i] = new G4LogicalVolume(fsomaS[i], fpWaterMaterial, "Soma");
>> 306 fsomaLV[i] ->SetVisAttributes(fSomaColour);
>> 307 fsomaPV[i] = new G4PVPlacement(
>> 308 0, // no rotation
>> 309 G4ThreeVector(
>> 310 (fNeuronLoadParamz->GetPosSomacomp(i).x()-fNeuronLoadParamz->GetshiftX())*um,
>> 311 (fNeuronLoadParamz->GetPosSomacomp(i).y()-fNeuronLoadParamz->GetshiftY())*um,
>> 312 (fNeuronLoadParamz->GetPosSomacomp(i).z()-fNeuronLoadParamz->GetshiftZ())*um),
>> 313 fsomaLV[i],
>> 314 "Soma",
>> 315 boundingLV,
>> 316 false,
>> 317 i);
>> 318 fMassSomacomp[i] = fNeuronLoadParamz->GetMassSomacomp(i) ;
>> 319 fPosSomacomp[i] = fNeuronLoadParamz->GetPosSomacomp(i) ;
245 fpRegion->SetProductionCuts(cuts); 320 fpRegion->SetProductionCuts(cuts);
246 G4ThreeVector shift(fNeuronLoadParamz->Getsh << 321 fpRegion->AddRootLogicalVolume(fsomaLV[i]);
247 fNeuronLoadParamz->Getsh <<
248 322
249 G4int n = fNeuronLoadParamz->GetnbSomacomp() << 323 }
250 if (n <= 0) { << 324 }
251 G4cout << " ---- Soma not found! ---- " << <<
252 } <<
253 else { <<
254 G4cout << " ---- Soma for construction: -- <<
255 G4cout << " Total number of compartments i <<
256 fnbSomacomp = n; <<
257 fMassSomaTot = fNeuronLoadParamz->GetMassS <<
258 fMassSomacomp.resize(n, 0.0); <<
259 fPosSomacomp.resize(n); <<
260 fsomaS.resize(n, nullptr); <<
261 fsomaLV.resize(n, nullptr); <<
262 fsomaPV.resize(n, nullptr); <<
263 <<
264 for (G4int i = 0; i < n; ++i) { <<
265 fsomaS[i] = new G4Orb("Soma", fNeuronLoa <<
266 // you can change parameters of Soma wit <<
267 G4ThreeVector pos = (fNeuronLoadParamz-> <<
268 fsomaLV[i] = new G4LogicalVolume(fsomaS[ <<
269 fsomaLV[i]->SetVisAttributes(fSomaColour <<
270 fsomaPV[i] = new G4PVPlacement(0, // no <<
271 pos, fsom <<
272 fMassSomacomp[i] = fNeuronLoadParamz->Ge <<
273 fPosSomacomp[i] = fNeuronLoadParamz->Get <<
274 fpRegion->AddRootLogicalVolume(fsomaLV[i <<
275 } <<
276 } <<
277 325
278 // ========================================= << 326 // =======================================================================
279 // Structure-2: Dendrites 327 // Structure-2: Dendrites
280 328
281 n = fNeuronLoadParamz->GetnbDendritecomp(); << 329 // Active Geant4-DNA processes in this region
282 if (n <= 0) { << 330 fpRegion = new G4Region("Dendrites");
283 G4cout << " ---- Dendrites not found! ---- << 331
>> 332 fnbDendritecomp = fNeuronLoadParamz->GetnbDendritecomp() ;
>> 333 fMassDendTot = fNeuronLoadParamz->GetMassDendTot() ;
>> 334 fMassDendcomp = new G4double[fnbDendritecomp];
>> 335 fDistADendSoma = new G4double[fnbDendritecomp];
>> 336 fDistBDendSoma = new G4double[fnbDendritecomp];
>> 337 fPosDendcomp = new G4ThreeVector[fnbDendritecomp];
>> 338 if (fNeuronLoadParamz->GetnbDendritecomp()==0)
>> 339 {
>> 340 G4cout <<" ---- Dendrites not found! ---- "<< G4endl;
284 } 341 }
285 else { << 342 else
286 fnbDendritecomp = n; << 343 {
287 G4cout << " ---- Dendrites for constructio << 344 G4cout <<" ---- Dendrites for construction: ---- "<< G4endl;
288 G4cout << " Total number of compartments i << 345 G4cout << " Total number of compartments into Dendrites : "
289 << 346 <<fNeuronLoadParamz->GetnbDendritecomp()<<G4endl;
290 // Active Geant4-DNA processes in this reg << 347 for (G4int i=1; i<fNeuronLoadParamz->GetnbDendritecomp() ; i++)
291 auto region = new G4Region("Dendrites"); << 348 {
292 region->SetProductionCuts(cuts); << 349 fdendriteS [i] = new G4Tubs( "Dendrites",
293 << 350 0., fNeuronLoadParamz->GetRadDendcomp(i)*um,
294 fMassDendTot = fNeuronLoadParamz->GetMassD << 351 fNeuronLoadParamz->GetHeightDendcomp(i)*um/2.,
295 fMassDendcomp.resize(n, 0.0); << 352 0., 2.*pi );
296 fDistADendSoma.resize(n, 0.0); << 353 fdendriteLV[i] = new G4LogicalVolume(fdendriteS[i],fpWaterMaterial,
297 fDistBDendSoma.resize(n, 0.0); << 354 "Dendrites");
298 fPosDendcomp.resize(n); << 355 fdendriteLV[i] ->SetVisAttributes(fDendColour);
299 fdendriteS.resize(n, nullptr); << 356
300 fdendriteLV.resize(n, nullptr); << 357 fdendritePV[i] = new G4PVPlacement(
301 fdendritePV.resize(n, nullptr); << 358 // rot, // rotation checking with function ComputeTransformation!
302 for (G4int i = 1; i < n; ++i) { << 359 G4Transform3D( fNeuronLoadParamz->
303 fdendriteS[i] = new G4Tubs("Dendrites", << 360 GetRotDendcomp(i), // RotationMatrix with Inverse
304 fNeuronLoadPa << 361 G4ThreeVector(
305 fdendriteLV[i] = new G4LogicalVolume(fde << 362 (fNeuronLoadParamz->GetPosDendcomp(i).x()-fNeuronLoadParamz->GetshiftX())*um,
306 fdendriteLV[i]->SetVisAttributes(fDendCo << 363 (fNeuronLoadParamz->GetPosDendcomp(i).y()-fNeuronLoadParamz->GetshiftY())*um,
307 << 364 (fNeuronLoadParamz->GetPosDendcomp(i).z()-fNeuronLoadParamz->GetshiftZ())*um)),
308 G4ThreeVector pos = (fNeuronLoadParamz-> << 365 fdendriteLV[i],
309 // rotation checking with function Compu << 366 "Dendrites",
310 // RotationMatrix with Inverse << 367 boundingLV,
311 fdendritePV[i] = new G4PVPlacement(G4Tra << 368 false,
312 fdend << 369 i);
313 fMassDendcomp[i] = fNeuronLoadParamz->Ge << 370 fMassDendcomp[i] = fNeuronLoadParamz->GetMassDendcomp(i) ;
314 fPosDendcomp[i] = fNeuronLoadParamz->Get << 371 fPosDendcomp[i] = fNeuronLoadParamz->GetPosDendcomp(i) ;
315 fDistADendSoma[i] = fNeuronLoadParamz->G << 372 fDistADendSoma[i] = fNeuronLoadParamz->GetDistADendSoma(i) ;
316 fDistBDendSoma[i] = fNeuronLoadParamz->G << 373 fDistBDendSoma[i] = fNeuronLoadParamz->GetDistBDendSoma(i) ;
317 region->AddRootLogicalVolume(fdendriteLV << 374 fpRegion->SetProductionCuts(cuts);
>> 375 fpRegion->AddRootLogicalVolume(fdendriteLV[i]);
>> 376
318 } 377 }
319 } 378 }
320 << 379
321 // ========================================= 380 // =======================================================================
322 // Structure-3: Axon << 381 // Structure-3: Axon
>> 382
>> 383 // Active Geant4-DNA processes in this region
>> 384 fpRegion = new G4Region("Axon");
323 385
324 n = fNeuronLoadParamz->GetnbAxoncomp(); << 386 fnbAxoncomp = fNeuronLoadParamz->GetnbAxoncomp() ;
325 if (n <= 0) { << 387 fMassAxonTot = fNeuronLoadParamz->GetMassAxonTot() ;
326 G4cout << " ---- Axon not found! ---- " << << 388 fMassAxoncomp = new G4double[fnbAxoncomp];
>> 389 fDistAxonsoma = new G4double[fnbAxoncomp];
>> 390 fPosAxoncomp = new G4ThreeVector[fnbAxoncomp];
>> 391 if (fNeuronLoadParamz->GetnbAxoncomp()==0)
>> 392 {
>> 393 G4cout <<" ---- Axon not found! ---- "<< G4endl;
327 } 394 }
328 else { << 395 else
329 G4cout << " ---- Axon for construction: -- << 396 {
330 G4cout << " Total number of compartments i << 397 G4cout <<" ---- Axon for construction: ---- "<< G4endl;
331 // Active Geant4-DNA processes in this reg << 398 G4cout << " Total number of compartments into Axon : "
332 auto region = new G4Region("Axon"); << 399 << fNeuronLoadParamz->GetnbAxoncomp() <<G4endl;
333 region->SetProductionCuts(cuts); << 400 for (G4int i=1; i< fNeuronLoadParamz->GetnbAxoncomp() ; i++)
334 << 401 {
335 fnbAxoncomp = n; << 402 faxonS [i] = new G4Tubs( "Axon",
336 fMassAxonTot = fNeuronLoadParamz->GetMassA << 403 0., fNeuronLoadParamz->GetRadAxoncomp(i)*um,
337 fMassAxoncomp.resize(n, 0.0); << 404 fNeuronLoadParamz->GetHeightAxoncomp(i)*um/2.,
338 fDistAxonsoma.resize(n, 0.0); << 405 0., 2.*pi );
339 fPosAxoncomp.resize(n); << 406 faxonLV[i] = new G4LogicalVolume(faxonS[i], fpWaterMaterial, "Axon");
340 faxonS.resize(n, nullptr); << 407 faxonLV[i] ->SetVisAttributes(fAxonColour);
341 faxonLV.resize(n, nullptr); << 408 faxonPV[i] = new G4PVPlacement(G4Transform3D(
342 faxonPV.resize(n, nullptr); << 409 fNeuronLoadParamz->GetRotAxoncomp(i), // RotationMatrix with Inverse
343 << 410 G4ThreeVector(
344 for (G4int i = 1; i < n; ++i) { << 411 (fNeuronLoadParamz->GetPosAxoncomp(i).x()-fNeuronLoadParamz->GetshiftX())*um,
345 faxonS[i] = new G4Tubs("Axon", 0., fNeur << 412 (fNeuronLoadParamz->GetPosAxoncomp(i).y()-fNeuronLoadParamz->GetshiftY())*um,
346 fNeuronLoadParamz << 413 (fNeuronLoadParamz->GetPosAxoncomp(i).z()-fNeuronLoadParamz->GetshiftZ())*um)),
347 faxonLV[i] = new G4LogicalVolume(faxonS[ << 414 faxonLV[i],
348 faxonLV[i]->SetVisAttributes(fAxonColour << 415 "Axon",
349 // RotationMatrix with Inverse << 416 boundingLV,
350 G4ThreeVector pos = (fNeuronLoadParamz-> << 417 false,
351 faxonPV[i] = new G4PVPlacement(G4Transfo << 418 i);
352 faxonLV[i << 419 fMassAxoncomp[i] = fNeuronLoadParamz->GetMassAxoncomp(i) ;
353 fMassAxoncomp[i] = fNeuronLoadParamz->Ge << 420 fPosAxoncomp[i] = fNeuronLoadParamz->GetPosAxoncomp(i) ;
354 fPosAxoncomp[i] = fNeuronLoadParamz->Get << 421 fDistAxonsoma[i] = fNeuronLoadParamz->GetDistAxonsoma(i) ;
355 fDistAxonsoma[i] = fNeuronLoadParamz->Ge << 422 fpRegion->SetProductionCuts(cuts);
356 region->AddRootLogicalVolume(faxonLV[i]) << 423 fpRegion->AddRootLogicalVolume(faxonLV[i]);
>> 424
357 } 425 }
358 } << 426 }
359 // ========================================= 427 // =======================================================================
360 // Structure-4: Spines << 428 // Structure-4: Spines
361 if (fNeuronLoadParamz->GetnbSpinecomp() == 0 << 429 if (fNeuronLoadParamz->GetnbSpinecomp()==0)
362 G4cout << " ---- Spines not found! ---- " << 430 {
>> 431 G4cout <<" ---- Spines not found! ---- "<< G4endl;
363 } 432 }
364 else { << 433 else
365 G4cout << " ---- Spines for construction: << 434 {
366 G4cout << " Total number of compartments i << 435 G4cout <<" ---- Spines for construction: ---- "<< G4endl;
367 << G4endl; << 436 G4cout << " Total number of compartments into Spines : "
>> 437 << fNeuronLoadParamz->GetnbSpinecomp() <<G4endl;
368 } 438 }
369 << 439
370 G4cout << "\n ---- End of Neuron Constructio << 440 G4cout <<"\n ---- End of Neuron Construction! ---- "
371 << "\n ============================== << 441 << "\n ========================================================== \n"
372 << G4endl; << 442 << G4endl;
373 << 443
374 // ========================================= 444 // =======================================================================
375 // Active Geant4-DNA processes in BoundingSl 445 // Active Geant4-DNA processes in BoundingSlice with whole neuron structure
376 // fpRegion = new G4Region("BoundingSlice"); << 446 //fpRegion = new G4Region("BoundingSlice");
377 // fpRegion->SetProductionCuts(cuts); << 447 //fpRegion->SetProductionCuts(cuts);
378 // fpRegion->AddRootLogicalVolume(boundingLV << 448 // fpRegion->AddRootLogicalVolume(boundingLV);
379 449
380 return worldPV; 450 return worldPV;
381 } 451 }
382 452