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