source: trunk/examples/novice/gemc/src/MDetectorConstruction.cc@ 1190

Last change on this file since 1190 was 807, checked in by garnier, 17 years ago

update

File size: 10.2 KB
Line 
1// %%%%%%%%%%
2// G4 headers
3// %%%%%%%%%%
4#include "G4Box.hh"
5#include "G4GeometryManager.hh"
6#include "G4LogicalVolumeStore.hh"
7#include "G4LogicalVolume.hh"
8#include "G4PhysicalVolumeStore.hh"
9#include "G4ProductionCuts.hh"
10#include "G4PVPlacement.hh"
11#include "G4RunManager.hh"
12#include "G4SolidStore.hh"
13#include "G4VisAttributes.hh"
14#include "G4SDManager.hh"
15
16
17#include "G4OpBoundaryProcess.hh"
18
19// %%%%%%%%%%%%%
20// gemc headers
21// %%%%%%%%%%%%%
22#include "MDetectorConstruction.h"
23#include "MDetectorMessenger.h"
24
25MDetectorConstruction::MDetectorConstruction(gemc_opts Opts)
26{
27 gemcOpt = Opts;
28 // mdetectorMessenger = new MDetectorMessenger(this);
29}
30
31MDetectorConstruction::~MDetectorConstruction()
32{
33 ;
34}
35
36
37G4VPhysicalVolume* MDetectorConstruction::Construct()
38{
39 string hd_msg = gemcOpt.args["LOG_MSG"].args + " Construction: >> ";
40 double VERB = gemcOpt.args["G4P_VERBOSITY"].arg ;
41 string catch_v = gemcOpt.args["CATCH"].args;
42 string hall_mat = gemcOpt.args["HALL_MATERIAL"].args;
43 string hall_field = gemcOpt.args["HALL_FIELD"].args;
44
45
46 // Clean old geometry, if any
47 G4GeometryManager::GetInstance()->OpenGeometry();
48 G4PhysicalVolumeStore::GetInstance()->Clean();
49 G4LogicalVolumeStore::GetInstance()->Clean();
50 G4SolidStore::GetInstance()->Clean();
51
52 // Experimental hall is a 40 meters box
53 detector queen;
54 queen.name = "root";
55 queen.mother = "akasha";
56 queen.description = "mother of us all";
57 queen.pos = G4ThreeVector(0, 0, 0);
58 queen.rot = G4RotationMatrix(G4ThreeVector(1, 0, 0),
59 G4ThreeVector(0, 1, 0),
60 G4ThreeVector(0, 0, 1));
61 queen.type = "Box";
62 queen.dimensions.push_back(20*m);
63 queen.dimensions.push_back(20*m);
64 queen.dimensions.push_back(20*m);
65 queen.material = hall_mat;
66 queen.magfield = hall_field;
67 queen.visible = 0;
68 queen.ncopy = 0;
69 queen.scanned = 1;
70 queen.create_solid(gemcOpt, Hall_Map);
71 queen.create_logical_volume(MMats, gemcOpt);
72 queen.create_physical_volumes(gemcOpt, NULL);
73 HasMagfield(queen, gemcOpt);
74
75
76 (*Hall_Map)["root"] = queen;
77
78 if(VERB > 3 || catch_v == "root") cout << hd_msg << " " << queen ;
79
80 cout << hd_msg << " Building Detector from Geometry STL Map..." << endl << endl;
81
82 // ######################################################
83 // Resetting Detector Map "scanned". Propagating "exist".
84 // ######################################################
85 if(VERB > 2) cout << hd_msg << " Mapping Physical Detector..." << endl << endl;
86
87 for( map<string, detector>::iterator i = Hall_Map->begin(); i!=Hall_Map->end(); i++)
88 {
89 if(VERB > 3) cout << hd_msg << " Scanning Detector " << i->first << endl;
90
91/* // Propagating "exist" flag to child - just in case
92 if(i->second.mother != "akasha")
93 if((*Hall_Map)[i->second.mother].exist == 0 && i->second.exist == 1)
94 {
95 if(VERB > 2) cout << hd_msg << "\t" << i->second.mother << " is not activated. Its child " << i->second.name
96 << " will be disactivated as well." << endl;
97 i->second.exist = 0;
98 }*/
99 if(i->first != "root") i->second.scanned = 0;
100
101 }
102 if(VERB > 2) cout << endl;
103
104 // ########################################################################
105 // Building Solids, Logical Volumes, Physical Volumes from the detector Map
106 // ########################################################################
107
108 vector<string> relatives;
109 string mom, kid;
110 for( map<string, detector>::iterator i = Hall_Map->begin(); i!=Hall_Map->end(); i++)
111 {
112 if(i->first != "root") relatives.push_back(i->second.name);
113 while(relatives.size() > 0)
114 {
115 mom = (*Hall_Map)[relatives.back()].mother;
116 kid = (*Hall_Map)[relatives.back()].name;
117
118 // Mom doesn't exists in the map. Stopping everything.
119 if((*Hall_Map)[mom].name == "" && mom != "akasha" && (*Hall_Map)[mom].scanned == 0)
120 {
121 cout << hd_msg << " Mom <" << mom << "> does not exist for <" << kid
122 << ">. We have a No Child Left Behind policy. Exiting. " << endl << endl;
123 exit(0);
124 }
125
126 if(VERB > 3 || kid.find(catch_v) != string::npos)
127 {
128 for(int i=0; i<relatives.size()-1; i++) cout << "\t";
129 cout << hd_msg << " Checking " << kid << ", child of " << mom
130 << ", for a living ancestor. This Geneaology Depth is " << relatives.size() << "." << endl;
131 }
132
133 // Mom is built, kid not built yet. Build the kid.
134 if((*Hall_Map)[kid].scanned == 0 && (*Hall_Map)[mom].scanned == 1)
135 {
136 if(VERB > 3 || kid.find(catch_v) != string::npos)
137 {
138 for(int i=0; i<relatives.size()-1; i++) cout << "\t";
139 cout << hd_msg << " Found: " << kid
140 << " is not built yet but its mommie " << mom << " is. Building " << kid << "..." << endl;
141 }
142 (*Hall_Map)[kid].create_solid(gemcOpt, Hall_Map);
143
144 if((*Hall_Map)[kid].create_logical_volume(MMats, gemcOpt))
145 {
146 (*Hall_Map)[kid].create_physical_volumes(gemcOpt, (*Hall_Map)[mom].GetLogical());
147 IsSensitive((*Hall_Map)[kid], gemcOpt);
148 HasMagfield((*Hall_Map)[kid], gemcOpt);
149 }
150 (*Hall_Map)[kid].scanned = 1;
151 }
152
153 // if the kid still doesn't exists it means its mom doesn't exist. Need to go up one level.
154 if((*Hall_Map)[kid].scanned == 0 && mom != "akasha") relatives.push_back(mom);
155
156 // the kid has been built. Can go down one step in geneaology
157 if((*Hall_Map)[kid].scanned == 1 && relatives.size())
158 {
159 if(VERB > 3 || kid.find(catch_v) != string::npos)
160 cout << hd_msg << " " << kid << " is built." << endl << endl;
161
162 relatives.pop_back();
163 }
164 }
165 }
166
167 G4OpticalSurface* Mirror_Surface = new G4OpticalSurface("ScintSurface");
168 Mirror_Surface->SetType(dielectric_metal);
169 Mirror_Surface->SetFinish(polished);
170 Mirror_Surface->SetModel(unified);
171
172 const G4int nEntries_Mirror = 2;
173 G4double PhotonEnergy_Mirror[nEntries_Mirror] = { 2.034*eV, 4.136*eV };
174
175 G4double Reflectivity_Mirror[nEntries_Mirror] = {0.99, 0.99};
176
177 G4MaterialPropertiesTable* Mirror_MPT = new G4MaterialPropertiesTable();
178 Mirror_MPT->AddProperty("REFLECTIVITY", PhotonEnergy_Mirror, Reflectivity_Mirror, nEntries_Mirror);
179
180 Mirror_Surface->SetMaterialPropertiesTable(Mirror_MPT);
181
182// G4LogicalBorderSurface* Mirror = new G4LogicalBorderSurface("Mirror",
183// (*Hall_Map)["root"].GetPhysical(),
184// (*Hall_Map)["Mirror"].GetPhysical(), Mirror_Surface);
185
186
187 return (*Hall_Map)["root"].GetPhysical();
188}
189
190 #include "G4UserLimits.hh"
191
192void MDetectorConstruction::IsSensitive(detector detect, gemc_opts Opts)
193{
194 string hd_msg = gemcOpt.args["LOG_MSG"].args + " Sensitivity: >> ";
195 double VERB = gemcOpt.args["HIT_VERBOSITY"].arg ;
196 string catch_v = gemcOpt.args["CATCH"].args;
197
198 string sensi = detect.sensitivity;
199
200 if(sensi != "no")
201 {
202 // Sensitive_region->AddRootLogicalVolume(detect.GetLogical());
203 G4SDManager* SDman = G4SDManager::GetSDMpointer();
204 if (!SDman) throw "Can't get G4SDManager";
205
206 if(VERB > 5 || detect.name.find(catch_v) != string::npos)
207 cout << hd_msg << " " << detect.name << " is sensitive." << endl;
208
209 // The key for the SD, Region, PC maps is the same so we can only check SD
210 map<string, MSensitiveDetector*>::iterator itr = SeDe_Map.find(sensi);
211 if(itr == SeDe_Map.end())
212 {
213 if(VERB > 3 || detect.name.find(catch_v) != string::npos)
214 cout << endl << hd_msg << " Sensitive detector <" << sensi
215 << "> doesn't exist yet. Adding <" << sensi << ">. " << endl;
216
217 SeDe_Map[sensi] = new MSensitiveDetector(sensi, gemcOpt);
218
219 // Creating G4 Region, assigning Production Cut to it.
220 SeRe_Map[sensi] = new G4Region(sensi);
221 SePC_Map[sensi] = new G4ProductionCuts;
222 SePC_Map[sensi] ->SetProductionCut(SeDe_Map[sensi]->SDID.ProdThreshold);
223 SeRe_Map[sensi]->SetProductionCuts(SePC_Map[sensi]);
224
225 if(VERB > 3 || detect.name.find(catch_v) != string::npos)
226 cout << hd_msg << " Max Step for Sensitive detector <" << sensi
227 << ">: " << SeDe_Map[sensi]->SDID.MaxStep/mm << " mm." << endl
228 << hd_msg << " Production Cut for Sensitive detector <" << sensi
229 << ">: " << SeDe_Map[sensi]->SDID.ProdThreshold/mm << " mm." << endl << endl;
230
231 // Pass Detector Map Pointer to Sensitive Detector
232 SeDe_Map[sensi]->Hall_Map = Hall_Map;
233
234 SDman->AddNewDetector( SeDe_Map[sensi]);
235
236 }
237 detect.setSensitivity(SeDe_Map[sensi]);
238
239 // Setting Max Acceptable Step for this SD
240 detect.SetUserLimits(new G4UserLimits(SeDe_Map[sensi]->SDID.MaxStep));
241
242 map<string, G4Region*>::iterator itrRE = SeRe_Map.find(sensi);
243 SeRe_Map[sensi]->AddRootLogicalVolume(detect.GetLogical());
244
245
246
247
248 }
249
250}
251
252
253void MDetectorConstruction::HasMagfield(detector detect, gemc_opts Opts)
254{
255 string hd_msg = gemcOpt.args["LOG_MSG"].args + " Magnetic Field: >> ";
256 double VERB = gemcOpt.args["MGN_VERBOSITY"].arg ;
257 string catch_v = gemcOpt.args["CATCH"].args;
258
259 string magf = detect.magfield;
260
261 if(magf != "no")
262 {
263 if(VERB > 5 || detect.name.find(catch_v) != string::npos)
264 cout << hd_msg << " " << detect.name << " is inside " << magf << " magnetic field." << endl;
265
266 map<string, MagneticField>::iterator itr = FieldMap->find(magf);
267 if(itr == FieldMap->end())
268 {
269 cout << hd_msg << " Magnetic Field <" << magf
270 << "> is not defined. Exiting." << endl;
271 exit(0);
272 }
273 if(itr->second.get_MFM() == NULL)
274 {
275 cout << hd_msg << " Magnetic Field <" << magf
276 << "> doesn't exist yet. Adding <" << magf << ">. " << endl;
277 itr->second.create_MFM();
278 }
279
280 if(itr->second.get_MFM() != NULL)
281 {
282 detect.AssignMFM(itr->second.get_MFM());
283
284 if(VERB > 6 || detect.name.find(catch_v) != string::npos)
285 cout << hd_msg << " Field <" << magf << "> assigned to " << detect.name << " " << detect.GetLogical() << endl;
286
287 }
288 }
289
290}
291
292
293void MDetectorConstruction::UpdateGeometry()
294{
295 cout << "Updating geometry... " << endl;
296 G4RunManager::GetRunManager()->GeometryHasBeenModified();
297}
298
299
300
301
302
Note: See TracBrowser for help on using the repository browser.