source: trunk/source/g3tog4/src/G3Division.cc@ 1116

Last change on this file since 1116 was 965, checked in by garnier, 17 years ago

update g3tog4
File size: 25.2 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26//
27// $Id: G3Division.cc,v 1.17 2006/06/29 18:12:53 gunter Exp $
28// GEANT4 tag $Name: geant4-09-02-ref-02 $
29//
30// by I.Hrivnacova, V.Berejnoi 13.10.99
31
32#include "G3Division.hh"
33#include "G3VolTableEntry.hh"
34#include "G3toG4MakeSolid.hh"
35#include "G4Para.hh"
36#include "G3Pos.hh"
37#include "G4LogicalVolume.hh"
38#include "G4VPhysicalVolume.hh"
39#include "G4PVPlacement.hh"
40#include "G4PVReplica.hh"
41#ifndef G3G4_NO_REFLECTION
42#include "G4ReflectionFactory.hh"
43#endif
44
45G3VolTableEntry* G4CreateVTE(G4String vname, G4String shape, G4int nmed,
46 G4double Rpar[], G4int npar);
47
48G3Division::G3Division(G3DivType type, G3VolTableEntry* vte,
49 G3VolTableEntry* mvte, G4int nofDivisions,
50 G4int iaxis, G4int nmed, G4double c0, G4double step)
51 : fType(type),
52 fVTE(vte),
53 fMVTE(mvte),
54 fNofDivisions(nofDivisions),
55 fIAxis(iaxis),
56 fNmed(nmed),
57 fC0(c0),
58 fStep(step),
59 fLowRange(0.),
60 fHighRange(0.),
61 fWidth(0.),
62 fOffset(0.),
63 fAxis(kXAxis)
64{
65 fVTE->SetHasNegPars(true);
66}
67
68G3Division::G3Division(G3VolTableEntry* vte, G3VolTableEntry* mvte,
69 const G3Division& division)
70 : fVTE(vte),
71 fMVTE(mvte)
72{
73 // only "input" parameters are copied from division
74 fType = division.fType;
75 fNofDivisions = division.fNofDivisions;
76 fIAxis = division.fIAxis;
77 fNmed = division.fNmed;
78 fC0 = division.fC0;
79 fStep = division.fStep;
80
81 // other parameters are set as in standard constructor
82 fLowRange = 0.;
83 fHighRange = 0.;
84 fWidth = 0.;
85 fOffset = 0.;
86 fAxis = kXAxis;
87 fVTE->SetHasNegPars(true);
88}
89
90G3Division::~G3Division()
91{}
92
93// public methods
94
95void G3Division::UpdateVTE()
96{
97 if (fVTE->HasNegPars() && !(fMVTE->HasNegPars())) {
98
99 // set nmed from mother
100 if (fNmed == 0) fNmed = fMVTE->GetNmed();
101 fVTE->SetNmed(fNmed);
102
103 SetRangeAndAxis();
104
105 // create envelope (if necessary)
106 // and solid
107 G3VolTableEntry* envVTE = 0;
108 if (fType == kDvn) envVTE = Dvn();
109 else if (fType == kDvn2) envVTE = Dvn2();
110 else if (fType == kDvt) envVTE = Dvt();
111 else if (fType == kDvt2) envVTE = Dvt2();
112
113 if (envVTE) {
114 // reset mother <-> daughter
115 fMVTE->ReplaceDaughter(fVTE, envVTE);
116 fVTE->ReplaceMother(fMVTE, envVTE);
117 envVTE->AddDaughter(fVTE);
118 envVTE->AddMother(fMVTE);
119
120 // replace mother with envelope
121 fMVTE = envVTE;
122 }
123 }
124}
125
126void G3Division::CreatePVReplica()
127{
128 G4String name = fVTE->GetName();
129 G4LogicalVolume* lv = fVTE->GetLV();
130 G4LogicalVolume* mlv = fMVTE->GetLV();
131
132 G4String shape = fMVTE->GetShape();
133 if (shape == "PARA") {
134 // The para volume cannot be replicated using G4PVReplica.
135 // (Replicating a volume along a cartesian axis means "slicing" it
136 // with slices -perpendicular- to that axis.)
137
138 // position the replicated elements
139 for (G4int i=0; i<fNofDivisions; i++) {
140 G4ThreeVector position = G4ThreeVector();
141 position[fIAxis-1] = fLowRange + fWidth/2. + i*fWidth;
142 if (position.y()!=0.)
143 position.setX(position.y()*((G4Para*)lv->GetSolid())->GetTanAlpha());
144
145 #ifndef G3G4_NO_REFLECTION
146 G4ReflectionFactory::Instance()
147 ->Place(G4Translate3D(position), name, lv, mlv, 0, i);
148
149 #else
150 new G4PVPlacement(0, position, lv, name, mlv, 0, i);
151
152 #endif
153 }
154
155 // G4PVReplica cannot be created
156 return;
157 }
158
159 #ifdef G3G4DEBUG
160 G4cout << "Create G4PVReplica name " << name << " logical volume name "
161 << lv->GetName() << " mother logical volme name "
162 << mlv->GetName() << " axis " << fAxis << " ndivisions "
163 << fNofDivisions << " width " << fWidth << " Offset "
164 << fOffset << G4endl;
165 #endif
166
167 #ifndef G3G4_NO_REFLECTION
168 G4ReflectionFactory::Instance()
169 ->Replicate(name, lv, mlv, fAxis, fNofDivisions, fWidth, fOffset);
170
171 #else
172 new G4PVReplica(name, lv, mlv, fAxis, fNofDivisions, fWidth, fOffset);
173
174 #endif
175}
176
177// private methods
178
179void G3Division::Exception(G4String where, G4String what) {
180 G4Exception("G3Division::" + where + " for " + what + " is not implemented");
181}
182
183void G3Division::SetRangeAndAxis()
184// set fHighRange, fLowRange, fAxis
185{
186 G4String shape = fMVTE->GetShape();
187 G4double *Rpar = fMVTE->GetRpar();
188
189 switch (fIAxis) {
190 case 1: fAxis = kXAxis;
191 break;
192 case 2: fAxis = kYAxis;
193 break;
194 case 3: fAxis = kZAxis;
195 break;
196 default: G4Exception("G3Division: wrong iaxis defenition");
197 }
198
199 if ( shape == "BOX" ) {
200 fHighRange = Rpar[fIAxis-1]*cm;
201 fLowRange = -fHighRange;
202 }
203 else if ( shape == "TRD1" ) {
204 if (fIAxis == 1){
205 fHighRange = std::max(Rpar[0]*cm, Rpar[1]*cm);
206 }
207 else if( fIAxis == 2) {
208 fHighRange = Rpar[2]*cm;
209 }
210 else if( fIAxis == 3) {
211 fHighRange = Rpar[3]*cm;
212 }
213 fLowRange = - fHighRange;
214 }
215 else if ( shape == "TRD2" ) {
216 if (fIAxis == 1){
217 fHighRange = std::max(Rpar[0]*cm, Rpar[1]*cm);
218 }
219 else if( fIAxis == 2) {
220 fHighRange = std::max(Rpar[2]*cm, Rpar[3]*cm);
221 }
222 else if( fIAxis == 3) {
223 fHighRange = Rpar[4]*cm;
224 }
225 }
226 else if ( shape == "TRAP" ) {
227 if ( fIAxis == 3 ) fHighRange = Rpar[0]*cm;
228 else fHighRange = 0.;
229 fLowRange = -fHighRange;
230 }
231 else if ( shape == "TUBE" ) {
232 if (fIAxis == 1){
233 fHighRange = Rpar[1]*cm;
234 fLowRange = Rpar[0]*cm;
235 fAxis = kRho;
236 }
237 else if( fIAxis == 2) {
238 fHighRange = 360.*deg;
239 fLowRange = 0.;
240 fAxis = kPhi;
241 }
242 else if( fIAxis == 3) {
243 fHighRange = Rpar[2]*cm;
244 fLowRange = -fHighRange;
245 }
246 }
247 else if ( shape == "TUBS" ) {
248 if (fIAxis == 1){
249 fHighRange = Rpar[1]*cm;
250 fLowRange = Rpar[0]*cm;
251 fAxis = kRho;
252 }
253 else if( fIAxis == 2) {
254
255 fLowRange = Rpar[3]*deg;
256 fHighRange = Rpar[4]*deg - fLowRange;
257 if ( Rpar[4]*deg <= fLowRange )fHighRange = fHighRange + 360.*deg;
258 fHighRange = fHighRange + fLowRange;
259 fAxis = kPhi;
260 }
261 else if( fIAxis == 3) {
262 fHighRange = Rpar[2]*cm;
263 fLowRange = -fHighRange;
264 }
265 }
266 else if ( shape == "CONE" ) {
267 if (fIAxis == 1){
268 fHighRange = std::max(Rpar[2]*cm,Rpar[4]*cm);
269 fLowRange = std::max(Rpar[1]*cm,Rpar[3]*cm);
270 fAxis = kRho;
271 }
272 else if( fIAxis == 2) {
273
274 fLowRange = 0.;
275 fHighRange = 360.*deg;
276 fAxis = kPhi;
277 }
278 else if( fIAxis == 3) {
279 fHighRange = Rpar[0]*cm;
280 fLowRange = -fHighRange;
281 }
282 }
283 else if ( shape == "CONS" ) {
284 if (fIAxis == 1){
285 fHighRange = std::max(Rpar[2]*cm,Rpar[4]*cm);
286 fLowRange = std::max(Rpar[1]*cm,Rpar[3]*cm);
287 fAxis = kRho;
288 }
289 else if( fIAxis == 2) {
290
291 fLowRange = Rpar[5]*deg;
292 fHighRange = Rpar[6]*deg - fLowRange;
293 if ( Rpar[6]*deg <= fLowRange )fHighRange = fHighRange + 360.*deg;
294 fHighRange = fHighRange + fLowRange;
295 fAxis = kPhi;
296 }
297 else if( fIAxis == 3) {
298 fHighRange = Rpar[2]*cm;
299 fLowRange = -fHighRange;
300 }
301 }
302 else if ( shape == "SPHE" ) {
303 if (fIAxis == 1){
304 fHighRange = Rpar[1]*cm;
305 fLowRange = Rpar[0]*cm;
306 fAxis = kRho;
307 }
308 else if( fIAxis == 2) {
309 fLowRange = std::min(Rpar[2]*deg,Rpar[3]*deg);
310 fHighRange = std::max(Rpar[2]*deg,Rpar[3]*deg);
311 fAxis = kPhi;
312 }
313 else if( fIAxis == 3) {
314 fLowRange = std::min(Rpar[4]*deg,Rpar[5]*deg);
315 fHighRange = std::max(Rpar[4]*deg,Rpar[5]*deg);
316 fAxis = kPhi; // ??????
317 }
318 }
319 else if ( shape == "PARA" ) {
320 fHighRange = Rpar[fIAxis-1]*cm;
321 fLowRange = -fHighRange;
322 }
323 else if ( shape == "PGON" ) {
324 G4int i;
325 G4int nz = G4int(Rpar[3]);
326
327 G4double pPhi1 = Rpar[0]*deg;
328 G4double dPhi = Rpar[1]*deg;
329
330 G4double *DzArray = new G4double[nz];
331 G4double *Rmax = new G4double[nz];
332 G4double *Rmin = new G4double[nz];
333 G4double rangehi[3], rangelo[3];
334 rangehi[0] = -kInfinity ;
335 rangelo[0] = kInfinity ;
336 rangehi[2] = -kInfinity ;
337 rangelo[2] = kInfinity ;
338
339 for(i=0; i<nz; i++)
340 {
341 G4int i4=3*i+4;
342 G4int i5=i4+1;
343 G4int i6=i4+2;
344
345 DzArray[i] = Rpar[i4]*cm;
346 Rmin[i] = Rpar[i5]*cm;
347 Rmax[i] = Rpar[i6]*cm;
348 rangelo[0] = std::min(rangelo[0], Rmin[i]);
349 rangehi[0] = std::max(rangehi[0], Rmax[i]);
350 rangelo[2] = std::min(rangelo[2], DzArray[i]);
351 rangehi[2] = std::max(rangehi[2], DzArray[i]);
352 }
353 for (i=0;i<nz;i++){
354 assert(Rmin[i]>=0 && Rmax[i]>=Rmin[i]);
355 }
356 rangehi[1] = pPhi1 + dPhi;
357 rangelo[1] = pPhi1;
358 fHighRange = rangehi[fIAxis-1];
359 fLowRange = rangelo[fIAxis-1];
360 if (fIAxis == 1)fAxis = kRho;
361 else if (fIAxis == 2)fAxis = kPhi;
362 else if (fIAxis == 3)fAxis = kZAxis;
363
364 delete [] DzArray;
365 delete [] Rmin;
366 delete [] Rmax;
367
368 }
369 else if ( shape == "PCON" ) {
370
371 G4int i;
372 G4double pPhi1 = Rpar[0]*deg;
373 G4double dPhi = Rpar[1]*deg;
374 G4int nz = G4int(Rpar[2]);
375
376 G4double *DzArray = new G4double[nz];
377 G4double *Rmax = new G4double[nz];
378 G4double *Rmin = new G4double[nz];
379 G4double rangehi[3],rangelo[3];
380
381 rangehi[0] = -kInfinity ;
382 rangelo[0] = kInfinity ;
383 rangehi[2] = -kInfinity ;
384 rangelo[2] = kInfinity ;
385
386 for(i=0; i<nz; i++){
387 G4int i4=3*i+3;
388 G4int i5=i4+1;
389 G4int i6=i4+2;
390
391 DzArray[i] = Rpar[i4]*cm;
392 Rmin[i] = Rpar[i5]*cm;
393 Rmax[i] = Rpar[i6]*cm;
394 rangelo[0] = std::min(rangelo[0], Rmin[i]);
395 rangehi[0] = std::max(rangehi[0], Rmax[i]);
396 rangelo[2] = std::min(rangelo[2], DzArray[i]);
397 rangehi[2] = std::max(rangehi[2], DzArray[i]);
398 }
399 for (i=0;i<nz;i++){
400 assert(Rmin[i]>=0 && Rmax[i]>=Rmin[i]);
401 }
402 rangehi[1] = pPhi1 + dPhi;
403 rangelo[1] = pPhi1;
404 fHighRange = rangehi[fIAxis-1];
405 fLowRange = rangelo[fIAxis-1];
406 if (fIAxis == 1)fAxis = kRho;
407 else if (fIAxis == 2)fAxis = kPhi;
408 else if (fIAxis == 3)fAxis = kZAxis;
409
410
411 delete [] DzArray;
412 delete [] Rmin;
413 delete [] Rmax;
414 }
415 else if ( shape == "ELTU" || shape == "HYPE" || shape == "GTRA" ||
416 shape == "CTUB") {
417 Exception("SetRangeAndAxis", shape);
418 }
419 else {
420 Exception("SetRangeAndAxis", "Unknown shape" + shape);
421 }
422
423 // verbose
424 #ifdef G3G4DEBUG
425 G4cout << "Shape " << shape << " SetRangeAndAxis: "
426 << fLowRange << " " << fHighRange << " " << fAxis << G4endl;
427 #endif
428}
429
430G3VolTableEntry* G3Division::CreateEnvelope(G4String shape, G4double hi,
431 G4double lo, G4double par[], G4int npar)
432// create new VTE with G3Pos corresponding to the
433// envelope of divided volume
434{
435 // verbose
436 // G4cout << " G3Division::CreateEnvelope " << "fIAaxis= " << fIAxis
437 // << " hi= " << hi
438 // << " lo= " << lo
439 // << G4endl;
440
441 G4double *Rpar = new G4double[npar+2];
442 for (G4int i=0; i<npar; ++i){ Rpar[i] = par[i];}
443 G4double pos[3] = {0.,0.,0.};
444
445 if ( shape == "BOX" ) {
446 Rpar[fIAxis-1] = (hi - lo)/2./cm;
447 pos [fIAxis-1] = (hi + lo)/2.;
448 }
449 else if ( shape == "TRD1" ) {
450 if ( fIAxis == 1 || fIAxis == 2 ) {
451 Exception("CreateEnvelope","TRD1-x,y");
452 }
453 else if ( fIAxis == 3 ) {
454 // x = x1 + (c-z1)(x2 -x1)/(z2-z1)
455 G4double tn, x1, z1;
456 tn = (Rpar[1] - Rpar[0])/(2.* Rpar[3]);
457 x1 = Rpar[0]; z1 = -Rpar[3];
458 Rpar[0] = x1 + tn * (lo/cm - z1);
459 Rpar[1] = x1 + tn * (hi/cm - z1);
460 Rpar[3] = (hi - lo)/2./cm;
461 pos[2] = (hi + lo)/2.;
462 }
463 }
464 else if ( shape == "TRD2" ) {
465 if ( fIAxis == 1 || fIAxis == 2) {
466 Exception("CreateEnvelope","TRD2-x,y");
467 }
468 else if ( fIAxis == 3 ) {
469 // x = x1 + (c-z1)(x2 -x1)/(z2-z1)
470 // y = y1 + (c-z1)(y2 -y1)/(z2-z1)
471 G4double tn1, tn2, x1, y1, z1;
472 tn1 = (Rpar[1] - Rpar[0])/(2.* Rpar[4]);
473 tn2 = (Rpar[3] - Rpar[2])/(2.* Rpar[4]);
474 x1 = Rpar[0]; y1 = Rpar[2]; z1 = -Rpar[3];
475 Rpar[0] = x1 + tn1 * (lo/cm - z1);
476 Rpar[1] = x1 + tn1 * (hi/cm - z1);
477 Rpar[2] = y1 + tn2 * (lo/cm - z1);
478 Rpar[3] = y1 + tn2 * (hi/cm - z1);
479 Rpar[4] = (hi - lo)/2./cm;
480 pos[2] = (hi + lo)/2.;
481 }
482 }
483 else if ( shape == "TRAP" ) {
484 Exception("CreateEnvelope","TRAP-x,y,z");
485 }
486 else if ( shape == "TUBE" ) {
487 if ( fIAxis == 1 ) {
488 Rpar[0] = lo/cm;
489 Rpar[1] = hi/cm;
490 }
491 else if ( fIAxis == 2 ) {
492 Rpar[3] = lo/deg;
493 Rpar[4] = hi/deg;
494 npar = npar + 2;
495 shape = "TUBS";
496 }
497 else if ( fIAxis == 3 ) {
498 Rpar[2] = (hi - lo)/2./cm;
499 pos [2] = (hi + lo)/2.;
500 }
501 }
502 else if ( shape == "TUBS" ) {
503 if ( fIAxis == 1 ) {
504 Rpar[0] = lo/cm;
505 Rpar[1] = hi/cm;
506 }
507 else if ( fIAxis == 2 ) {
508 Rpar[3] = lo/deg;
509 Rpar[4] = hi/deg;
510 }
511 else if ( fIAxis == 3 ) {
512 Rpar[2] = (hi - lo)/2./cm;
513 pos [2] = (hi + lo)/2.;
514 }
515 }
516 else if ( shape == "CONE" ) {
517 if ( fIAxis == 1) {
518 Exception("CreateEnvelope","CONE-x,z");
519 }
520 else if ( fIAxis == 2 ) {
521 Rpar[5] = lo/deg;
522 Rpar[6] = hi/deg;
523 npar = npar + 2;
524 shape = "CONS";
525 }
526 else if ( fIAxis == 3 ) {
527 G4double tn1, tn2, rmin, rmax, z1;
528 tn1 = (Rpar[3] - Rpar[1])/(2.* Rpar[0]);
529 tn2 = (Rpar[4] - Rpar[2])/(2.* Rpar[0]);
530 rmin = Rpar[1]; rmax = Rpar[2]; z1 = -Rpar[0];
531 Rpar[1] = rmin + tn1 * (lo/cm - z1);
532 Rpar[3] = rmin + tn1 * (hi/cm - z1);
533 Rpar[2] = rmax + tn2 * (lo/cm - z1);
534 Rpar[4] = rmax + tn2 * (hi/cm - z1);
535 Rpar[0] = (hi - lo)/2./cm;
536 pos[2] = (hi + lo)/2.;
537 }
538 }
539 else if ( shape == "CONS" ) {
540 if ( fIAxis == 1 ) {
541 Exception("CreateEnvelope","CONS-x");
542 }
543 else if ( fIAxis == 2 ) {
544 Rpar[5] = lo/deg;
545 Rpar[6] = hi/deg;
546 }
547 else if ( fIAxis == 3 ) {
548 G4double tn1, tn2, rmin, rmax, z1;
549 tn1 = (Rpar[3] - Rpar[1])/(2.* Rpar[0]);
550 tn2 = (Rpar[4] - Rpar[2])/(2.* Rpar[0]);
551 rmin = Rpar[1]; rmax = Rpar[2]; z1 = -Rpar[0];
552 Rpar[1] = rmin + tn1 * (lo/cm - z1);
553 Rpar[3] = rmin + tn1 * (hi/cm - z1);
554 Rpar[2] = rmax + tn2 * (lo/cm - z1);
555 Rpar[4] = rmax + tn2 * (hi/cm - z1);
556 Rpar[0] = (hi - lo)/2./cm;
557 pos[2] = (hi + lo)/2.;
558 }
559 }
560 else if ( shape == "SPHE" ) {
561 Exception("CreateEnvelope","SPHE-x,y,z");
562 }
563 else if ( shape == "PARA" ) {
564 Exception("CreateEnvelope","PARA-x,y,z");
565 }
566 else if ( shape == "PGON" ) {
567 if ( fIAxis == 2) {
568 Rpar[0] = lo/deg;
569 Rpar[1] = hi/deg;
570 // rotm = ???
571 }
572 else {
573 Exception("CreateEnvelope","PGON-x,z");
574 }
575 }
576 else if ( shape == "PCON" ) {
577 if ( fIAxis == 2) {
578 Rpar[0] = lo/deg;
579 Rpar[1] = hi/deg;
580 // rotm = ???
581 }
582 else {
583 Exception("CreateEnvelope","PCON-x,z");
584 }
585 }
586 else {
587 Exception("CreateEnvelope", "Unknown shape" + shape);
588 }
589
590 // create new VTE corresponding to envelope
591 G4String envName = fVTE->GetName() + "_ENV";
592 G3VolTableEntry* envVTE
593 = G4CreateVTE(envName, shape, fNmed, Rpar, npar);
594
595 // create a G3Pos object and add it to envVTE
596 G4String motherName = fMVTE->GetMasterClone()->GetName();
597 G4ThreeVector* offset = new G4ThreeVector(pos[0],pos[1],pos[2]);
598 G4String only = "ONLY";
599 G3Pos* aG3Pos = new G3Pos(motherName, 1, offset, 0, only);
600 envVTE->AddG3Pos(aG3Pos);
601
602 delete [] Rpar;
603
604 return envVTE;
605}
606
607void G3Division::CreateSolid(G4String shape, G4double par[], G4int npar)
608// create the solid corresponding to divided volume
609// and set the fOffset for replica
610{
611 G4double *Rpar = new G4double[npar+2];
612 for (G4int i=0; i<npar; ++i){ Rpar[i] = par[i];}
613
614 // verbose
615 // G4cout << "G3Division::CreateSolid volume before: "
616 // << fVTE->GetName() << " " << shape << G4endl;
617 // G4cout << " npar,Rpar: " << npar;
618 // for (G4int ii = 0; ii < npar; ++ii) G4cout << " " << Rpar[ii];
619 // G4cout << G4endl;
620
621 if ( shape == "BOX" ) {
622 if ( fIAxis == 1 ) Rpar[0] = fWidth/2./cm;
623 else if ( fIAxis == 2 ) Rpar[1] = fWidth/2./cm;
624 else if ( fIAxis == 3 ) Rpar[2] = fWidth/2./cm;
625 }
626 else if ( shape == "TRD1" ) {
627 if ( fIAxis == 1 || fIAxis == 2 ) {
628 Exception("CreateSolid", "TRD1-x,y");
629 }
630 else if ( fIAxis == 3 ) {
631 Rpar[3] = fWidth/2./cm;
632 }
633 }
634 else if ( shape == "TRD2" ) {
635 if ( fIAxis == 1 || fIAxis == 2 ) {
636 Exception("CreateSolid", "TRD2-x,y");
637 }
638 else if ( fIAxis == 3 ) {
639 Rpar[4] = fWidth/2./cm;
640 }
641 }
642 else if ( shape == "TRAP" ) {
643 if ( fIAxis == 1 || fIAxis == 2) {
644 Exception("CreateSolid", "TRAP-x,y");
645 }
646 else if ( fIAxis == 3 ) {
647 Rpar[0] = fWidth/2./cm;
648 }
649 }
650 else if ( shape == "TUBE" ) {
651 if ( fIAxis == 1 ) {
652 Rpar[1] = Rpar[0] + fWidth/cm;
653 fOffset = Rpar[0]*cm;
654 }
655 else if ( fIAxis == 2 ) {
656 Rpar[3] = 0.;
657 Rpar[4] = fWidth/deg;
658 shape = "TUBS";
659 npar = npar + 2;
660 }
661 else if ( fIAxis == 3 ) {
662 Rpar[2] = fWidth/2./cm;
663 }
664 }
665 else if ( shape == "TUBS" ) {
666 if ( fIAxis == 1 ) {
667 Rpar[1] = Rpar[0] + fWidth/cm;
668 fOffset = Rpar[0]*cm;
669 }
670 else if ( fIAxis == 2 ) {
671 fOffset = Rpar[3]*deg;
672 Rpar[3] = 0.;
673 Rpar[4] = fWidth/deg;
674 }
675 else if ( fIAxis == 3 ) {
676 Rpar[2] = fWidth/2./cm;
677 }
678 }
679 else if ( shape == "CONE" ) {
680 if ( fIAxis == 1 ) {
681 Exception("CreateSolid", "CONE-x");
682 }
683 else if ( fIAxis == 2 ) {
684 Rpar[5] = 0.;
685 Rpar[6] = fWidth/deg;
686 shape = "CONS";
687 npar = npar + 2;
688 }
689 else if ( fIAxis == 3 ) {
690 Rpar[0] = fWidth/2./cm;
691 }
692 }
693 else if ( shape == "CONS" ) {
694 if ( fIAxis == 1 ) {
695 Exception("CreateSolid", "CONS-x");
696 }
697 else if ( fIAxis == 2 ) {
698 fOffset = Rpar[5]*deg;
699 Rpar[5] = 0.;
700 Rpar[6] = fWidth/deg;
701 }
702 else if ( fIAxis == 3 ) {
703 Rpar[0] = fWidth/2./cm;
704 }
705 }
706 else if (shape == "PARA") {
707 if ( fIAxis == 1 ) {
708 Rpar[0] = fWidth/2./cm;
709 }
710 else if ( Rpar[4] == 0. && Rpar[5] == 0. ) {
711 // only special case for axis 2,3 is supported
712 if ( fIAxis == 2 ) {
713 Rpar[1] = fWidth/2./cm;
714 }
715 else if ( fIAxis == 3) {
716 Rpar[2] = fWidth/2./cm;
717 }
718 }
719 else
720 Exception("CreateSolid", shape);
721 }
722 else if (shape == "SPHE") {
723 Exception("CreateSolid", shape);
724 }
725 else if ( shape == "PGON" ) {
726 if ( fIAxis == 2 ) {
727 fOffset = Rpar[0]*deg;
728 Rpar[0] = 0.;
729 Rpar[1] = fWidth/deg;
730 Rpar[2] = 1.;
731 }
732 else
733 Exception("CreateSolid", shape);
734 }
735 else if ( shape == "PCON" ) {
736 if ( fIAxis == 2 ) {
737 fOffset = Rpar[0]*deg;
738 Rpar[0] = 0.;
739 Rpar[1] = fWidth/deg;
740 }
741 else {
742 Exception("CreateSolid", shape);
743 }
744 }
745 else {
746 Exception("CreateSolid", "Unknown shape" + shape);
747 }
748
749 // create solid and set it to fVTE
750 G4bool hasNegPars;
751 G4bool deferred;
752 G4bool okAxis[3];
753 G4VSolid* solid
754 = G3toG4MakeSolid(fVTE->GetName(), shape, Rpar, npar, hasNegPars, deferred, okAxis);
755
756 if (hasNegPars) {
757 G4String name = fVTE->GetName();
758 G4Exception("CreateSolid VTE " + name + " has negative parameters.");
759 }
760
761 // update vte
762 fVTE->SetSolid(solid);
763 fVTE->SetNRpar(npar, Rpar);
764 fVTE->SetHasNegPars(hasNegPars);
765
766 // verbose
767 // G4cout << "G3Division::CreateSolid volume after: "
768 // << fVTE->GetName() << " " << shape << G4endl;
769 // G4cout << " npar,Rpar: " << npar;
770 // for (G4int iii = 0; iii < npar; ++iii) G4cout << " " << Rpar[iii];
771 // G4cout << G4endl;
772}
773
774
775G3VolTableEntry* G3Division::Dvn()
776{
777 // no envelope need to be created
778
779 // get parameters from mother
780 G4String shape = fMVTE->GetShape();
781 G4double* Rpar = fMVTE->GetRpar();
782 G4int npar = fMVTE->GetNpar();
783
784 // set width for replica and create solid
785 fWidth = (fHighRange - fLowRange)/fNofDivisions;
786 CreateSolid(shape, Rpar, npar);
787
788 return 0;
789}
790
791G3VolTableEntry* G3Division::Dvn2()
792{
793 // to be defined as const of this class
794 G4double Rmin = 0.0001*cm;
795
796 G4String shape = fMVTE->GetShape();
797 G4double* Rpar = fMVTE->GetRpar();
798 G4int npar = fMVTE->GetNpar();
799
800 G4double c0 = fC0;
801 if (fAxis == kPhi) c0 = c0*deg;
802 else c0 = c0*cm;
803
804 // create envelope (if needed)
805 G3VolTableEntry* envVTE = 0;
806 if( std::abs(c0 - fLowRange) > Rmin) {
807 envVTE = CreateEnvelope(shape, fHighRange, c0, Rpar, npar);
808 Rpar = envVTE->GetRpar();
809 npar = envVTE->GetNpar();
810 }
811
812 // set width for replica and create solid
813 fWidth = (fHighRange - c0)/fNofDivisions;
814 CreateSolid(shape, Rpar, npar);
815
816 return envVTE;
817}
818
819G3VolTableEntry* G3Division::Dvt()
820{
821 // to be defined as const of this class
822 G4double Rmin = 0.0001*cm;
823
824 // get parameters from mother
825 G4String shape = fMVTE->GetShape();
826 G4double* Rpar = fMVTE->GetRpar();
827 G4int npar = fMVTE->GetNpar();
828
829 // calculate the number of divisions
830 G4int ndvmx = fNofDivisions;
831 G4double step = fStep;
832
833 if (fAxis == kPhi) step = step*deg;
834 else step = step*cm;
835
836 G4int ndiv = G4int((fHighRange - fLowRange + Rmin)/step);
837 // to be added warning
838 if (ndvmx > 255) ndvmx = 255;
839 if (ndiv > ndvmx && ndvmx > 0 ) ndiv = ndvmx;
840
841 // create envVTE (if needed)
842 G3VolTableEntry* envVTE = 0;
843 G4double delta = std::abs((fHighRange - fLowRange) - ndiv*step);
844 if (delta > Rmin) {
845 envVTE
846 = CreateEnvelope(shape, fHighRange-delta/2., fLowRange+delta/2.,
847 Rpar, npar);
848 Rpar = envVTE->GetRpar();
849 npar = envVTE->GetNpar();
850 }
851
852 // set width for replica and create solid
853 fWidth = step;
854 fNofDivisions = ndiv;
855 CreateSolid(shape, Rpar, npar);
856
857 return envVTE;
858}
859
860G3VolTableEntry* G3Division::Dvt2()
861{
862 // to be defined as const of this class
863 G4double Rmin = 0.0001*cm;
864
865 // get parameters from mother
866 G4String shape = fMVTE->GetShape();
867 G4double* Rpar = fMVTE->GetRpar();
868 G4int npar = fMVTE->GetNpar();
869
870 // calculate the number of divisions
871 G4int ndvmx = fNofDivisions;
872 G4double step = fStep;
873 G4double c0 = fC0;
874
875 if(fAxis == kPhi){
876 step = step*deg;
877 c0 = c0*deg;
878 }
879 else {
880 step = step*cm;
881 c0 = c0*cm;
882 }
883
884 G4int ndiv = G4int((fHighRange - c0 + Rmin)/step);
885 // to be added warning
886 if (ndvmx > 255) ndvmx = 255;
887 if (ndiv > ndvmx && ndvmx > 0 ) ndiv = ndvmx;
888
889 // create envelope (if needed)
890 G3VolTableEntry* envVTE = 0;
891 G4double delta = std::abs((fHighRange - c0) - ndiv*step);
892 if (std::abs(c0 - fLowRange) > Rmin) {
893 envVTE
894 = CreateEnvelope(shape, fHighRange-delta/2., c0+delta/2., Rpar, npar);
895 Rpar = envVTE->GetRpar();
896 npar = envVTE->GetNpar();
897 }
898
899 // set with for replica and create solid
900 fWidth = step;
901 fNofDivisions = ndiv;
902 CreateSolid(shape, Rpar, npar);
903
904 return envVTE;
905}
Note: See TracBrowser for help on using the repository browser.