BDSIM: src/BDSDetectorConstruction.cc Source File

00001 //  
00002 //   BDSIM, (C) 2001-2007
00003 //   
00004 //   version 0.5-dev
00005 //  
00006 //
00007 //
00008 //   Geometry construction
00009 //
00010 //
00011 //   History
00012 //     19 May 2008 by Marchioni v.0.5-dev
00013 //     18 Mar 2008 by Malton v.0.5-dev
00014 //      3 Oct 2007 by Malton v.0.4
00015 //     21 Nov 2006 by Agapov v.0.3
00016 //     28 Mar 2006 by Agapov v.0.2
00017 //     15 Dec 2005 by Agapov beta
00018 //
00019 
00020 //=================================================================
00021 
00022 #include "BDSExecOptions.hh"
00023 #include "BDSGlobalConstants.hh"
00024 #include "BDSDebug.hh"
00025 
00026 #include "BDSDetectorConstruction.hh"
00027 
00028 #include "G4UserLimits.hh"
00029 #include "G4GeometryManager.hh"
00030 #include "G4Region.hh"
00031 #include "G4ProductionCuts.hh"
00032 
00033 #include "G4Tubs.hh"
00034 #include "G4Box.hh"
00035 #include "G4LogicalVolume.hh"
00036 #include "G4VPhysicalVolume.hh"
00037 #include "G4PVPlacement.hh"
00038 #include "G4UniformMagField.hh"
00039 #include "G4TransportationManager.hh"
00040 #include "G4PropagatorInField.hh"
00041 #include "G4SDManager.hh"
00042 #include "G4RunManager.hh"
00043 
00044 #include "G4VisAttributes.hh"
00045 #include "G4Colour.hh"
00046 #include "globals.hh"
00047 #include "G4ios.hh"
00048 #include <iostream>
00049 #include <list>
00050 #include <map>
00051 #include "BDSAcceleratorComponent.hh"
00052 
00053 #include "G4Navigator.hh"
00054 #include "G4UniformMagField.hh"
00055 
00056 #include "G4Material.hh"
00057 #include "BDSEnergyCounterSD.hh"
00058 
00059 // elements
00060 #include "BDSDrift.hh"
00061 #include "BDSPCLDrift.hh"
00062 #include "BDSSectorBend.hh"
00063 #include "BDSRBend.hh"
00064 #include "BDSKicker.hh"
00065 #include "BDSQuadrupole.hh"
00066 #include "BDSSextupole.hh"
00067 #include "BDSOctupole.hh"
00068 #include "BDSDecapole.hh"
00069 #include "BDSTMultipole.hh"
00070 #include "BDSRfCavity.hh"
00071 #include "BDSSolenoid.hh"
00072 #include "BDSSampler.hh"
00073 #include "BDSSamplerCylinder.hh"
00074 #include "BDSDump.hh"
00075 #include "BDSLaserWire.hh"
00076 #include "BDSLWCalorimeter.hh"
00077 #include "BDSMuSpoiler.hh"
00078 #include "BDSTransform3D.hh"
00079 #include "BDSElement.hh"
00080 #include "BDSComponentOffset.hh"
00081 #include "BDSCollimator.hh"
00082 // output interface
00083 #include "BDSOutput.hh"
00084 #include "BDSComponentFactory.hh"
00085 
00086 #include "G4MagneticField.hh"
00087 
00088 // GMAD interface
00089 #include "ggmad.hh"
00090 
00091 #include "G4VSampler.hh"
00092 #include "G4GeometrySampler.hh"
00093 #include "G4IStore.hh"
00094  
00095 using namespace std;
00096 
00097 //====================================
00098 
00099 //typedef list<BDSAcceleratorComponent*>  BDSBeamline;
00100 
00101 typedef list<BDSEnergyCounterSD*>  ECList;
00102 ECList* theECList;
00103 
00104 //extern BDSGlobalConstants* BDSGlobalConstants::Instance();
00105 extern G4int gflash;
00106 extern G4double gflashemax;
00107 extern G4double gflashemin;
00108 
00109 //--------------------------
00110 // SYNCHROTRON RAD ***
00111 G4double BDSLocalRadiusOfCurvature=DBL_MAX;// Used in Mean Free Path calc.
00112 //--------------------------
00113 
00114 G4Material* aMaterial;
00115 extern G4double NumSpoilerRadLen;
00116 
00117 typedef std::map<G4String,int> LogVolCountMap;
00118 LogVolCountMap* LogVolCount;
00119 
00120 typedef std::map<G4String,G4LogicalVolume*> LogVolMap;
00121 LogVolMap* LogVol;
00122 
00123 G4RotationMatrix* RotY90=new G4RotationMatrix();
00124 G4RotationMatrix* RotYM90=new G4RotationMatrix();
00125 
00126 G4RotationMatrix* RotX90=new G4RotationMatrix();
00127 G4RotationMatrix* RotXM90=new G4RotationMatrix();
00128 
00129 G4RotationMatrix* RotYM90X90=new G4RotationMatrix();
00130 G4RotationMatrix* RotYM90XM90=new G4RotationMatrix();
00131 
00132 //G4Navigator* StepperNavigator;
00133 //G4Navigator* QuadNavigator;
00134 
00135 //=========================================
00136 extern BDSOutput* bdsOutput;
00137 extern G4bool verbose;
00138 extern G4bool outline;
00139 
00140 #ifdef DEBUG
00141 bool debug = true;
00142 #else
00143 bool debug = false;
00144 #endif
00145 
00146 //=================================================================
00147 
00148 
00149 
00150 BDSDetectorConstruction::BDSDetectorConstruction():
00151   itsGeometrySampler(NULL),precisionRegion(NULL),gasRegion(NULL),
00152   solidWorld(NULL),logicWorld(NULL),physiWorld(NULL),
00153   magField(NULL),BDSUserLimits(NULL),BDSSensitiveDetector(NULL),
00154   itsIStore(NULL)
00155 {  
00156   verbose    = BDSExecOptions::Instance()->GetVerbose();
00157   outline    = BDSExecOptions::Instance()->GetOutline();
00158   gflash     = BDSExecOptions::Instance()->GetGFlash();
00159   gflashemax = BDSExecOptions::Instance()->GetGFlashEMax();
00160   gflashemin = BDSExecOptions::Instance()->GetGFlashEMin();
00161   
00162   //initialize world size vector
00163   itsWorldSize.push_back(0);
00164   itsWorldSize.push_back(1);
00165   itsWorldSize.push_back(2);
00166 
00167 // create commands for interactive definition of the beamline  
00168   G4double pi_ov_2 = asin(1.);
00169 
00170   RotY90->rotateY(pi_ov_2);
00171 
00172   RotYM90->rotateY(-pi_ov_2);
00173 
00174   RotX90->rotateX(pi_ov_2);
00175 
00176   RotXM90->rotateX(-pi_ov_2);
00177 
00178   RotYM90X90->rotateY(-pi_ov_2);
00179   RotYM90X90->rotateX( pi_ov_2);
00180 
00181   RotYM90XM90->rotateY(-pi_ov_2);
00182   RotYM90XM90->rotateX(-pi_ov_2);
00183 
00184   // GlashStuff                                                                                                                                                         
00185   theParticleBounds  = new GFlashParticleBounds();              // Energy Cuts to kill particles                                                                
00186   theParticleBounds->SetMaxEneToParametrise(*G4Electron::ElectronDefinition(),gflashemax*GeV);
00187   theParticleBounds->SetMinEneToParametrise(*G4Electron::ElectronDefinition(),gflashemin*GeV);
00188   theParticleBounds->SetEneToKill(*G4Electron::ElectronDefinition(),BDSGlobalConstants::Instance()->GetThresholdCutCharged());
00189   
00190   theParticleBounds->SetMaxEneToParametrise(*G4Positron::PositronDefinition(),gflashemax*GeV);
00191   theParticleBounds->SetMinEneToParametrise(*G4Positron::PositronDefinition(),gflashemin*GeV);
00192   theParticleBounds->SetEneToKill(*G4Positron::PositronDefinition(),BDSGlobalConstants::Instance()->GetThresholdCutCharged());
00193 
00194   theParticleBoundsVac  = new GFlashParticleBounds();              // Energy Cuts to kill particles                                                                
00195   theParticleBoundsVac->SetMaxEneToParametrise(*G4Electron::ElectronDefinition(),0*GeV);
00196   theParticleBoundsVac->SetMaxEneToParametrise(*G4Positron::PositronDefinition(),0*GeV);
00197 
00198   G4cout << __METHOD_NAME__ << "theParticleBounds - min E - electron: " 
00199          << theParticleBounds->GetMinEneToParametrise(*G4Electron::ElectronDefinition())/GeV<< " GeV" << G4endl;
00200   G4cout << __METHOD_NAME__ << "theParticleBounds - max E - electron: " 
00201          << theParticleBounds->GetMaxEneToParametrise(*G4Electron::ElectronDefinition())/GeV<< G4endl;
00202   G4cout << __METHOD_NAME__ << "theParticleBounds - kill E - electron: " 
00203          << theParticleBounds->GetEneToKill(*G4Electron::ElectronDefinition())/GeV<< G4endl;
00204   G4cout << __METHOD_NAME__ << "theParticleBounds - min E - positron: " 
00205          << theParticleBounds->GetMinEneToParametrise(*G4Positron::PositronDefinition())/GeV<< G4endl;
00206   G4cout << __METHOD_NAME__ << "theParticleBounds - max E - positron: " 
00207          << theParticleBounds->GetMaxEneToParametrise(*G4Positron::PositronDefinition())/GeV<< G4endl;
00208   G4cout << __METHOD_NAME__ << "theParticleBounds - kill E - positron: " 
00209          << theParticleBounds->GetEneToKill(*G4Positron::PositronDefinition())/GeV<< G4endl;
00210 
00211   theHitMaker          = new GFlashHitMaker();                    // Makes the EnergieSpots 
00212 
00213 }
00214 
00215 //=================================================================
00216 
00217 G4VPhysicalVolume* BDSDetectorConstruction::Construct()
00218 {
00219   theECList=new ECList;
00220 
00221   LogVolCount=new LogVolCountMap();
00222 
00223   LogVol=new LogVolMap();
00224 
00225   gasRegion = new G4Region("gasRegion");
00226   G4ProductionCuts* theGasProductionCuts = new G4ProductionCuts();
00227   theGasProductionCuts->SetProductionCut(1*m,G4ProductionCuts::GetIndex("gamma"));
00228   theGasProductionCuts->SetProductionCut(1*m,G4ProductionCuts::GetIndex("e-"));
00229   theGasProductionCuts->SetProductionCut(1*m,G4ProductionCuts::GetIndex("e+"));
00230   gasRegion->SetProductionCuts(theGasProductionCuts);
00231 
00232 
00233   if (verbose || debug) G4cout << "-->starting BDS construction \n"<<G4endl;
00234 
00235   return ConstructBDS(beamline_list);
00236 
00237 
00238 
00239 
00240 
00241 
00242 
00243 }
00244 
00245 
00246 G4VPhysicalVolume* BDSDetectorConstruction::ConstructBDS(list<struct Element>& beamline_list)
00247 {
00248   //
00249   // set default output formats:
00250   //
00251   G4cout.precision(10);
00252   
00253   
00254   //
00255   // convert the parsed atom list to list of Geant4 G4Elements
00256   //
00257   list<struct Element>::iterator it;
00258 
00259 
00260   if (verbose || debug) G4cout << "parsing the atom list..."<< G4endl;
00261   for(it = atom_list.begin();it!=atom_list.end();it++)
00262   {
00263 #ifdef DEBUG
00264     G4cout << "---->adding Atom, "
00265            << "name= " << (*it).name << " "
00266            << "symbol= " << (*it).symbol << " "
00267            << "Z= " << (*it).Z << " "
00268            << "A= " << (*it).A << "g/mole "
00269            << G4endl;
00270 #endif
00271 
00272     BDSMaterials::Instance()->AddElement((*it).name,(*it).symbol,(*it).Z,(*it).A);
00273   }
00274   if (verbose || debug) G4cout << "size of atom list: "<< atom_list.size() << G4endl;
00275 
00276 
00277   //
00278   // convert the parsed material list to list of Geant4 G4Materials
00279   //
00280   if (verbose || debug) G4cout << "parsing the material list..."<< G4endl;
00281   for(it = material_list.begin();it!=material_list.end();it++)
00282   {
00283     if((*it).Z != 0) {
00284 #ifdef DEBUG  
00285       G4cout << "---->adding Material, "
00286              << "name= "<< (*it).name << " "
00287              << "Z= " << (*it).Z << " "
00288              << "A= " << (*it).A << "g/mole "
00289              << "density= "<< (*it).density << "g/cm3 "
00290              << G4endl;
00291 #endif
00292       BDSMaterials::Instance()->AddMaterial((*it).name,(*it).Z,(*it).A,(*it).density);
00293     }
00294     else if((*it).components.size() != 0){
00295 
00296       G4State itsState;
00297       if ((*it).state=="solid") itsState = kStateSolid;
00298       else
00299       if ((*it).state=="liquid") itsState = kStateLiquid;
00300       else
00301       if ((*it).state=="gas") itsState = kStateGas;
00302       else {
00303         G4cout << "Unknown material state "<< (*it).state 
00304                << ", setting it to default (solid)"
00305                << G4endl;
00306         (*it).state="solid";
00307         itsState = kStateSolid;
00308       }
00309 
00310       if((*it).componentsWeights.size()==(*it).components.size()) {
00311 
00312 #ifdef DEBUG 
00313         G4cout << "---->adding Material, "
00314                << "name= "<< (*it).name << " "
00315                << "density= "<< (*it).density << "g/cm3 "
00316                << "state= " << (*it).state << " "
00317                << "T= " << (*it).temper << "K "
00318                << "P= " << (*it).pressure << "atm "
00319                << "ncomponents= " << (*it).components.size() << " "
00320                << G4endl;
00321 #endif
00322 
00323         BDSMaterials::Instance()->AddMaterial((G4String)(*it).name,
00324                                   (G4double)(*it).density,
00325                                   (G4State)itsState,
00326                                   (G4double)(*it).temper,
00327                                   (G4double)(*it).pressure,
00328                                   (std::list<const char*>)(*it).components,
00329                                   (std::list<G4int>)(*it).componentsWeights);
00330       }
00331       else if((*it).componentsFractions.size()==(*it).components.size()) {
00332 
00333 #ifdef DEBUG 
00334         G4cout << "---->adding Material, "
00335         << "name= "<< (*it).name << " "
00336         << "density= "<< (*it).density << "g/cm3 "
00337         << "state= " << (*it).state << " "
00338         << "T= " << (*it).temper << "K "
00339         << "P= " << (*it).pressure << "atm "
00340         << "ncomponents= " << (*it).components.size() << " "
00341         << G4endl;
00342 #endif
00343         BDSMaterials::Instance()->AddMaterial((*it).name,
00344                                   (*it).density,
00345                                   itsState,
00346                                   (*it).temper,
00347                                   (*it).pressure,
00348                                   (*it).components,
00349                                   (*it).componentsFractions);
00350       }
00351       else {
00352         G4Exception("Badly defined material - number of components is not equal to number of weights or mass fractions!", "-1", FatalErrorInArgument, "");
00353         exit(1);
00354       }
00355     }
00356     else {
00357       G4Exception("Badly defined material - need more information!", "-1", FatalErrorInArgument, "");
00358       exit(1);
00359     }
00360   }
00361   if (verbose || debug) G4cout << "size of material list: "<< material_list.size() << G4endl;
00362 
00363 
00364   //
00365   // set global magnetic field first
00366   //
00367   SetMagField(0.0); // necessary to set a global field; so chose zero
00368 
00369   
00370   
00371   // convert the parsed element list to list of BDS elements
00372   //
00373   BDSComponentFactory* theComponentFactory = new BDSComponentFactory();
00374 
00375   if (verbose || debug) G4cout << "parsing the beamline element list..."<< G4endl;
00376   for(it = beamline_list.begin();it!=beamline_list.end();it++){
00377     G4cout << "BDSDetectorConstruction creating component..." << G4endl;
00378     BDSAcceleratorComponent* temp = theComponentFactory->createComponent(it, beamline_list);
00379     G4cout << "pushing onto back of beamline..." << G4endl;
00380     if(temp){
00381       BDSBeamline::Instance()->addComponent(temp);
00382       //For the outline file...
00383       BDSBeamline::Instance()->currentItem()->SetK1((*it).k1);
00384       BDSBeamline::Instance()->currentItem()->SetK2((*it).k2);
00385       BDSBeamline::Instance()->currentItem()->SetK3((*it).k3);
00386     }
00387     G4cout << "done." << G4endl;
00388   }
00389   delete theComponentFactory;
00390   theComponentFactory = NULL;
00391 
00392   if (verbose || debug) G4cout << "size of beamline element list: "<< beamline_list.size() << G4endl;
00393   if (verbose || debug) G4cout << "size of theBeamline: "<< BDSBeamline::Instance()->size() << G4endl;
00394     
00395   //
00396   // construct the component list
00397   //
00398 
00399   if (verbose || debug) G4cout << "now constructing geometry" << G4endl;
00400   
00401   list<BDSAcceleratorComponent*>::const_iterator iBeam;
00402   
00403   G4ThreeVector rtot = G4ThreeVector(0.,0.,0.);   // world dimension
00404   G4ThreeVector rlast = G4ThreeVector(0.,0.,0.);  // edge of last element coordinates
00405   G4ThreeVector rextentmax = G4ThreeVector(0.,0.,0.);  // extent
00406   G4ThreeVector rextentmin = G4ThreeVector(0.,0.,0.);  // extent
00407   G4ThreeVector rmin = G4ThreeVector(0.,0.,0.);
00408   G4ThreeVector rmax = G4ThreeVector(0.,0.,0.);
00409 
00410   G4ThreeVector localX = G4ThreeVector(1,0,0);  // local coordinate axis
00411   G4ThreeVector localY = G4ThreeVector(0,1,0);
00412   G4ThreeVector localZ = G4ThreeVector(0,0,1);
00413 
00414   G4double s_tot = 0; // position along the beamline
00415 
00416   // define geometry scope
00417   for(BDSBeamline::Instance()->first();!BDSBeamline::Instance()->isDone();BDSBeamline::Instance()->next())
00418     {
00419 
00420 #ifdef DEBUG 
00421       G4cout << BDSBeamline::Instance()->currentItem()->GetName() << "  "
00422              << BDSBeamline::Instance()->currentItem()->GetLength() << "  "
00423              << BDSBeamline::Instance()->currentItem()->GetAngle() << "  "
00424              << G4endl;
00425 #endif
00426 
00427       BDSBeamline::Instance()->currentItem()->SetSPos(s_tot+BDSBeamline::Instance()->currentItem()->GetArcLength()/2.0);
00428 
00429       // advance coordinates , but not for cylindrical sampler
00430       if(( BDSBeamline::Instance()->currentItem()->GetType() != "csampler") || ( BDSBeamline::Instance()->currentItem()->GetLength() <= BDSGlobalConstants::Instance()->GetSamplerLength() ) )
00431         {
00432           s_tot+= BDSBeamline::Instance()->currentItem()->GetArcLength();
00433 
00434           G4double angle=BDSBeamline::Instance()->currentItem()->GetAngle();
00435           if(!angle && BDSBeamline::Instance()->currentItem()->GetType()=="transform3d")
00436             angle=BDSBeamline::Instance()->currentItem()->GetPhi();
00437           G4double theta=BDSBeamline::Instance()->currentItem()->GetTheta();
00438           G4double psi=BDSBeamline::Instance()->currentItem()->GetPsi();
00439 
00440           // define new coordinate system local frame     
00441           localX.rotate(psi,localZ);
00442           localY.rotate(psi,localZ);
00443           localZ.rotate(psi,localZ);
00444 
00445           localX.rotate(angle,localY);
00446           localY.rotate(angle,localY);
00447           localZ.rotate(angle,localY);
00448 
00449           localX.rotate(theta,localX);
00450           localY.rotate(theta,localX);
00451           localZ.rotate(theta,localX);
00452 
00453           // advance the coordinate system translation
00454           rtot += localZ * BDSBeamline::Instance()->currentItem()->GetZLength();
00455 #ifdef DEBUG
00456           G4cout << BDSBeamline::Instance()->currentItem()->GetType() << " " << rtot << G4endl;
00457 #endif
00458         }
00459 
00460 
00461       rextentmax = rtot + localX*BDSBeamline::Instance()->currentItem()->GetXLength() + localY*BDSBeamline::Instance()->currentItem()->GetYLength();
00462       rextentmin = rtot - localX*BDSBeamline::Instance()->currentItem()->GetXLength() - localY*BDSBeamline::Instance()->currentItem()->GetYLength();
00463 
00464       //      rextentmax(0) = rtot(0) + (localX*BDSBeamline::Instance()->currentItem()->GetXLength())(0) + (localX*BDSBeamline::Instance()->currentItem()->GetYLength())(0);
00465       //      rextentmax(1) = rtot(1) + (localY*BDSBeamline::Instance()->currentItem()->GetXLength())(1) + (localY*BDSBeamline::Instance()->currentItem()->GetYLength())(1);
00466       //      rextentmax(2) = rtot(2) + (localZ*BDSBeamline::Instance()->currentItem()->GetXLength())(2) + (localZ*BDSBeamline::Instance()->currentItem()->GetYLength())(2);
00467 
00468       //      rextentmin(0) = rtot(0) - localX*BDSBeamline::Instance()->currentItem()->GetXLength() - localX*BDSBeamline::Instance()->currentItem()->GetYLength();
00469       //      rextentmin(1) = rtot(1) - localY*BDSBeamline::Instance()->currentItem()->GetXLength() - localY*BDSBeamline::Instance()->currentItem()->GetYLength();
00470       //      rextentmin(2) = rtot(2) - localZ*BDSBeamline::Instance()->currentItem()->GetXLength() - localZ*BDSBeamline::Instance()->currentItem()->GetYLength();
00471 
00472 
00473         for(int i=0; i<3; i++){
00474           rmax(i) = std::max(rextentmax(i),rmax(i));
00475           rmin(i) = std::min(rextentmin(i),rmin(i));
00476         }
00477 
00478     }
00479     
00480   BDSGlobalConstants::Instance()->SetTotalS(s_tot);
00481   
00482   // -----------------------------------
00483   
00484   // determine the world size
00485   
00486   G4String LocalName="World";
00487   
00488   SetWorldSizeX(1.5*( 2*std::max(fabs( rmin(0) ),fabs( rmax(0) ) ) ) *mm);
00489   SetWorldSizeY(1.5*( 2*std::max(fabs(rmin(1)),fabs(rmax(1))) ) *mm);
00490   SetWorldSizeZ(1.5*(fabs(rmin(2)) + fabs(rmax(2))) *mm);
00491 
00492   if(verbose || debug)
00493     {
00494       
00495       G4cout<<"minX="<<rmin(0)/m<<" m"<<" maxX="<<rmax(0)/m<<" m"<<G4endl;
00496       G4cout<<"minY="<<rmin(1)/m<<" m"<<" maxY="<<rmax(1)/m<<" m"<<G4endl;
00497       G4cout<<"minZ="<<rmin(2)/m<<" m"<<" maxZ="<<rmax(2)/m<<" m"<<G4endl;
00498 
00499       G4cout<<"itsWorldSizeX = "<<GetWorldSizeX()/m<<G4endl;
00500       G4cout<<"itsWorldSizeY = "<<GetWorldSizeY()/m<<G4endl;
00501       G4cout<<"itsWorldSizeZ = "<<GetWorldSizeZ()/m<<G4endl;
00502       
00503       G4cout<<"box size="<<BDSGlobalConstants::Instance()->GetComponentBoxSize()/m<<" m"<<G4endl;
00504       G4cout<<"s_tot="<<s_tot/m<<" m"<<G4endl;
00505     }
00506 
00507   bdsOutput->zMax=s_tot;
00508   bdsOutput->transMax=std::max(GetWorldSizeX(), GetWorldSizeY());
00509 
00510   solidWorld = new G4Box("World", GetWorldSizeX(), GetWorldSizeY(), GetWorldSizeZ());
00511     
00512   logicWorld = new G4LogicalVolume(solidWorld,         //its solid
00513                                    BDSMaterials::Instance()->GetMaterial(BDSGlobalConstants::Instance()->GetVacuumMaterial()), //its material
00514                                    "World");           //its name
00515   
00516   logicWorld->SetVisAttributes (G4VisAttributes::Invisible);
00517   // set world volume visibility for debugging
00518 #ifdef DEBUG 
00519   logicWorld->SetVisAttributes(new G4VisAttributes(true));      
00520 #endif
00521         
00522   // set default max step length (only for particles which have the
00523   // G4StepLimiter process enabled)
00524 #ifndef NOUSERLIMITS
00525   G4UserLimits* WorldUserLimits =new G4UserLimits();
00526   WorldUserLimits->SetMaxAllowedStep(10*m);
00527   WorldUserLimits->SetUserMinEkine(BDSGlobalConstants::Instance()->GetThresholdCutCharged());
00528   WorldUserLimits->SetUserMaxTime(BDSGlobalConstants::Instance()->GetMaxTime());
00529   logicWorld->SetUserLimits(WorldUserLimits);
00530 #endif
00531 
00532 
00533   G4cout<<"Creating regions..."<<G4endl;
00534   
00535   precisionRegion = new G4Region("precisionRegion");
00536    
00537   G4ProductionCuts* theProductionCuts = new G4ProductionCuts();
00538   
00539   if(BDSGlobalConstants::Instance()->GetProdCutPhotonsP()>0)
00540     theProductionCuts->SetProductionCut(BDSGlobalConstants::Instance()->GetProdCutPhotonsP(),G4ProductionCuts::GetIndex("gamma"));
00541 
00542   if(BDSGlobalConstants::Instance()->GetProdCutElectronsP()>0)
00543     theProductionCuts->SetProductionCut(BDSGlobalConstants::Instance()->GetProdCutElectronsP(),G4ProductionCuts::GetIndex("e-"));
00544 
00545   if(BDSGlobalConstants::Instance()->GetProdCutPositronsP()>0)
00546     theProductionCuts->SetProductionCut(BDSGlobalConstants::Instance()->GetProdCutPositronsP(),G4ProductionCuts::GetIndex("e+"));
00547   
00548   precisionRegion->SetProductionCuts(theProductionCuts);
00549 #ifndef NOUSERLIMITS
00550   precisionRegion->SetUserLimits(WorldUserLimits);
00551 #endif
00552   // world
00553 
00554   physiWorld = new G4PVPlacement((G4RotationMatrix*)0,          // no rotation
00555                                  (G4ThreeVector)0,             // at (0,0,0)
00556                                  logicWorld,    // its logical volume
00557                                  LocalName,     // its name
00558                                  NULL,          // its mother  volume
00559                                  false,         // no boolean operation
00560                                  0,             // copy number
00561                                  BDSGlobalConstants::Instance()->GetCheckOverlaps());           // overlap checking
00562 
00563 
00564 
00565   // sensitive detectors
00566 
00567   G4SDManager* SDman = G4SDManager::GetSDMpointer();
00568 
00569   G4bool use_graphics=true;
00570   G4ThreeVector TargetPos;
00571 
00572   // set default output formats:
00573   G4cout.precision(15);
00574   
00575 #ifdef DEBUG 
00576   G4cout<<"total length="<<s_tot/m<<"m"<<G4endl;
00577 #endif
00578   
00579   // reset counters:
00580   for(BDSBeamline::Instance()->first();!BDSBeamline::Instance()->isDone();BDSBeamline::Instance()->next()){
00581 
00582     // zero length components have no logical volumes
00583     if(BDSBeamline::Instance()->currentItem()->GetLength()!=0.)
00584       {
00585         G4String logVolName = BDSBeamline::Instance()->currentItem()->GetMarkerLogicalVolume()->GetName();
00586         (*LogVolCount)[logVolName]=1;
00587       }
00588   }
00589   
00590   if (verbose || debug) G4cout<<"starting placement procedure "<<G4endl;
00591   
00592   rtot = G4ThreeVector(0.,0.,0.);
00593   localX = G4ThreeVector(1.,0.,0.); 
00594   localY = G4ThreeVector(0.,1.,0.);
00595   localZ = G4ThreeVector(0.,0.,1.);
00596   
00597   G4RotationMatrix globalRotation;
00598   
00599   for(BDSBeamline::Instance()->first();!BDSBeamline::Instance()->isDone();BDSBeamline::Instance()->next())
00600     { 
00601       //BDSBeamline::Instance()->currentItem()->SetZLower(rtot.z());
00602        
00603       G4double angle=BDSBeamline::Instance()->currentItem()->GetAngle();
00604       G4double theta=BDSBeamline::Instance()->currentItem()->GetTheta();
00605       G4double psi = BDSBeamline::Instance()->currentItem()->GetPsi();
00606       G4double tilt = BDSBeamline::Instance()->currentItem()->GetTilt();
00607       G4double phi = BDSBeamline::Instance()->currentItem()->GetPhi();
00608       G4double length = BDSBeamline::Instance()->currentItem()->GetZLength();
00609 
00610       if( BDSBeamline::Instance()->currentItem()->GetType() == "transform3d")
00611         {
00612 
00613 #ifdef DEBUG 
00614           G4cout<<"transform3d : "<<phi<<" "<<theta<<" "<<psi<<G4endl;
00615 #endif
00616 
00617 
00618           rtot(0) += BDSBeamline::Instance()->currentItem()->GetXOffset(); 
00619           rtot(1) += BDSBeamline::Instance()->currentItem()->GetYOffset(); 
00620           rtot(2) += BDSBeamline::Instance()->currentItem()->GetZOffset(); 
00621 
00622           rlast(0) += BDSBeamline::Instance()->currentItem()->GetXOffset();
00623           rlast(1) += BDSBeamline::Instance()->currentItem()->GetYOffset(); 
00624           rlast(2) += BDSBeamline::Instance()->currentItem()->GetZOffset(); 
00625 
00626           globalRotation.rotate(psi,localZ);
00627           localX.rotate(psi,localZ);
00628           localY.rotate(psi,localZ);
00629 
00630           globalRotation.rotate(phi,localY);
00631           localX.rotate(phi,localY);
00632           localZ.rotate(phi,localY);
00633           
00634           globalRotation.rotate(theta,localX);
00635           localY.rotate(theta,localX);
00636           localZ.rotate(theta,localX);
00637                           
00638           continue;
00639         }
00640       
00641       // rotation matrix for component placement
00642       G4RotationMatrix *rotateComponent = new G4RotationMatrix;
00643 
00644       // tilted bends influence reference frame, otherwise just local tilt
00645       if(BDSBeamline::Instance()->currentItem()->GetType() == "sbend" || BDSBeamline::Instance()->currentItem()->GetType() == "rbend" )
00646         {
00647           globalRotation.rotate(tilt,localZ);
00648           localX.rotate(tilt,localZ);
00649           localY.rotate(tilt,localZ);
00650         }
00651       else 
00652         rotateComponent->rotateZ(tilt);
00653     
00654       // define center of bended elements from the previos coordinate frame
00655       G4ThreeVector zHalfAngle = localZ; 
00656 
00657       if( BDSBeamline::Instance()->currentItem()->GetType() == "sbend" || BDSBeamline::Instance()->currentItem()->GetType() == "rbend"  )
00658         zHalfAngle.rotate(angle/2,localY);
00659 
00660 #ifdef DEBUG
00661       G4cout<<"zHalfAngle="<<zHalfAngle<<G4endl;
00662       G4cout<<"localZ="<<localZ<<G4endl;
00663       G4cout<<"localX="<<localX<<G4endl;
00664       G4cout<<"localY="<<localY<<G4endl;
00665       G4cout<<"rlast="<<rlast<<G4endl;
00666 #endif
00667       
00668       // target position
00669       TargetPos = rlast + zHalfAngle *  ( length/2 + BDSGlobalConstants::Instance()->GetLengthSafety()/2 ) ;
00670 
00671 #ifdef DEBUG 
00672       G4cout<<"TargetPos="<<TargetPos<<G4endl;
00673 #endif
00674 
00675       // advance the coordinates, but not for cylindrical samplers 
00676       if( ( ( BDSBeamline::Instance()->currentItem()->GetType() != "csampler") || ( length <= BDSGlobalConstants::Instance()->GetSamplerLength() ) )  && ( BDSBeamline::Instance()->currentItem()->GetType()!=_ELEMENT ))
00677         {
00678 #ifdef DEBUG 
00679           G4cout << BDSBeamline::Instance()->currentItem()->GetType() << " "
00680                  << BDSBeamline::Instance()->currentItem()->GetName() << " "
00681                  << G4endl;
00682 #endif
00683           rtot = rlast + zHalfAngle * ( length/2 + BDSGlobalConstants::Instance()->GetLengthSafety()/2 );
00684           rlast = rtot + zHalfAngle * ( length/2 + BDSGlobalConstants::Instance()->GetLengthSafety()/2 );
00685         }
00686 
00687       // rotate to the previous reference frame
00688       rotateComponent->transform(globalRotation);
00689 
00690       rotateComponent->invert();
00691 
00692       // recompute global rotation
00693       // define new coordinate system local frame         
00694  
00695       // bends transform the coordinate system
00696       if( BDSBeamline::Instance()->currentItem()->GetType() == "sbend" || BDSBeamline::Instance()->currentItem()->GetType() == "rbend")
00697         {
00698           globalRotation.rotate(angle,localY);
00699           localX.rotate(angle,localY);
00700           localZ.rotate(angle,localY);
00701           
00702           
00703           globalRotation.rotate(theta,localX);
00704           localY.rotate(theta,localX);
00705           localZ.rotate(theta,localX);
00706           
00707           // bend trapezoids defined along z-axis
00708           rotateComponent->rotateY(-twopi/4-angle/2);                                           
00709         } else {
00710         if (BDSBeamline::Instance()->currentItem()->GetMarkerLogicalVolume()->GetSolid()->GetName().contains("trapezoid") ){
00711           rotateComponent->rotateY(-twopi/4); //Drift trapezoids defined along z axis 
00712         }
00713       }
00714     
00715                                         
00716 
00717       // zero length components not placed (transform3d)
00718       if(length!=0.){
00719         G4LogicalVolume* LocalLogVol=BDSBeamline::Instance()->currentItem()->GetMarkerLogicalVolume();
00720         
00721         G4String LogVolName=LocalLogVol->GetName();
00722         //--test
00723         G4VisAttributes* VisAtt1 = new G4VisAttributes(G4Colour(0.0, 1.0, 0.0));
00724         VisAtt1->SetForceSolid(true);  
00725         VisAtt1->SetVisibility(false);
00726         LocalLogVol->SetVisAttributes(VisAtt1);
00727         //------------
00728         int nCopy=(*LogVolCount)[LogVolName]-1;
00729         (*LogVolCount)[LogVolName]++;
00730 
00731 
00732         // add the wolume to one of the regions
00733         if(BDSBeamline::Instance()->currentItem()->GetPrecisionRegion())
00734           {
00735 #ifdef DEBUG
00736             G4cout<<"ELEMENT IS IN PRECISION REGION: "<<BDSBeamline::Instance()->currentItem()->GetPrecisionRegion()<< G4endl;
00737 #endif
00738             LocalLogVol->SetRegion(precisionRegion);
00739             precisionRegion->AddRootLogicalVolume(LocalLogVol);
00740           }
00741 
00742         
00743 #ifdef DEBUG
00744         G4cout<<"SETTING UP SENSITIVE VOLUME..."<< G4endl;
00745 #endif  
00746         // set up the sensitive volumes for energy counting:
00747         BDSBeamline::Instance()->currentItem()->SetCopyNumber(nCopy);
00748         G4LogicalVolume* SensVol=BDSBeamline::Instance()->currentItem()->GetSensitiveVolume();
00749 
00750         if(SensVol)
00751           {
00752             BDSEnergyCounterSD* ECounter=new BDSEnergyCounterSD(LogVolName);
00753             BDSBeamline::Instance()->currentItem()->SetBDSEnergyCounter(ECounter);
00754             SensVol->SetSensitiveDetector(ECounter);
00755             SDman->AddNewDetector(ECounter);
00756             theECList->push_back(ECounter);
00757           }
00758 
00759 
00760 #ifdef DEBUG
00761         G4cout<<"SETTING UP MULTIPLE SENSITIVE VOLUMES..."<< G4endl;
00762 #endif  
00763         vector<G4LogicalVolume*> MultipleSensVols = BDSBeamline::Instance()->currentItem()->GetMultipleSensitiveVolumes();
00764         if( ( BDSBeamline::Instance()->currentItem()->GetType()!="sampler" && BDSBeamline::Instance()->currentItem()->GetType()!="csampler" )
00765             && MultipleSensVols.size()>0)
00766           {
00767             for(G4int i=0; i<(G4int)MultipleSensVols.size(); i++)
00768               {
00769                 BDSEnergyCounterSD* ECounter=
00770                   new BDSEnergyCounterSD(LogVolName+BDSGlobalConstants::Instance()->StringFromInt(i));
00771                 BDSBeamline::Instance()->currentItem()->SetBDSEnergyCounter(ECounter);
00772                 MultipleSensVols.at(i)->SetSensitiveDetector(ECounter);
00773                 SDman->AddNewDetector(ECounter);
00774                 theECList->push_back(ECounter);      
00775                 
00776                 if(gflash){
00777                   if((MultipleSensVols.at(i)->GetRegion() != precisionRegion) && (BDSBeamline::Instance()->currentItem()->GetType()==_ELEMENT)){//If not in the precision region....
00778                     //              if(MultipleSensVols[i]->GetMaterial()->GetState()!=kStateGas){ //If the region material state is not gas, associate with a parameterisation
00779                     G4cout << "...adding " << MultipleSensVols[i]->GetName() << " to gFlashRegion" << G4endl;
00780                     /**********************************************                                                                                                                       
00781                      * Initialise shower model                                                                                                                                          
00782                      ***********************************************/
00783                     G4String rname = "gFlashRegion_" + MultipleSensVols[i]->GetName();
00784                     gFlashRegion.push_back(new G4Region(rname.c_str()));
00785                     G4String mname = "fastShowerModel" + rname;
00786                     G4cout << "...making parameterisation..." << G4endl;
00787                     theFastShowerModel.push_back(new BDSShowerModel(mname.c_str(),gFlashRegion.back()));
00788                     theParameterisation.push_back(new GFlashHomoShowerParameterisation(BDSMaterials::Instance()->GetMaterial(MultipleSensVols[i]->GetMaterial()->GetName().c_str()))); 
00789                     theFastShowerModel.back()->SetParameterisation(*theParameterisation.back());
00790                     theFastShowerModel.back()->SetParticleBounds(*theParticleBounds) ;
00791                     theFastShowerModel.back()->SetHitMaker(*theHitMaker);
00792                     if(MultipleSensVols[i]->GetMaterial()->GetState()!=kStateGas){ //If the region material state is not gas, associate with a parameterisation
00793                       theFastShowerModel.back()->SetFlagParamType(1);//Turn on the parameterisation for e-m showers starting in sensitive material and fitting in the current stack.
00794                       theFastShowerModel.back()->SetFlagParticleContainment(1);//Turn on containment
00795                     } else {
00796                       theFastShowerModel.back()->SetFlagParamType(0);//Turn on the parameterisation for e-m showers starting in sensitive material and fitting in the current stack.
00797                       theFastShowerModel.back()->SetFlagParticleContainment(0);//Turn on containment
00798 
00799                     }
00800                     MultipleSensVols[i]->SetRegion(gFlashRegion.back());
00801                     gFlashRegion.back()->AddRootLogicalVolume(MultipleSensVols[i]);
00802                     //              gFlashRegion.back()->SetUserLimits(new G4UserLimits(BDSBeamline::Instance()->currentItem()->GetLength()/10.0));
00803                     //              MultipleSensVols[i]->SetUserLimits(new G4UserLimits(BDSBeamline::Instance()->currentItem()->GetLength()/10.0));
00804                   }               
00805                 }
00806               }
00807 
00808           }
00809 
00810 
00811             //Loop through again, unsetting gas regions
00812 
00813         //      if(MultipleSensVols.size()>0){
00814         //        for(G4int i=0; i<(G4int)MultipleSensVols.size(); i++){
00815         //          if(MultipleSensVols[i]->GetMaterial()->GetState()==kStateGas){ //If the region material state is not gas, associate with a parameterisation
00816         //            MultipleSensVols[i]->GetRegion()->RemoveRootLogicalVolume(MultipleSensVols[i]);
00817               //              MultipleSensVols[i]->SetRegion(NULL);
00818 
00819               //              MultipleSensVols[i]->SetRegion(gasRegion);
00820               //              MultipleSensVols[i]->SetUserLimits(new G4UserLimits(BDSBeamline::Instance()->currentItem()->GetLength()/10.0));
00821         //      }
00822         //      }
00823         //      }
00824 
00825         
00826         
00827   
00828         if(BDSBeamline::Instance()->currentItem()->GetType()=="sampler") {
00829           LocalName=BDSBeamline::Instance()->currentItem()->GetName()+"_phys";
00830           bdsOutput->SampName.push_back(LocalName + "_" + BDSGlobalConstants::Instance()->StringFromInt(nCopy+1));
00831         } 
00832         else if(BDSBeamline::Instance()->currentItem()->GetType()=="csampler") {
00833           LocalName=BDSBeamline::Instance()->currentItem()->GetName()+"_phys";
00834           bdsOutput->CSampName.push_back(LocalName + "_" + BDSGlobalConstants::Instance()->StringFromInt(nCopy+1));
00835         } else {
00836           //it would be nice to set correctly names also for other elements...
00837           //but need to count them!
00838           LocalName=BDSBeamline::Instance()->currentItem()->GetName()+"_phys";
00839         }
00840 
00841         /*
00842         //for torus sbend
00843         if(BDSBeamline::Instance()->currentItem()->GetType() == "sbend") {
00844 
00845           G4double rho = length/fabs(angle);
00846 
00847           G4RotationMatrix* Rot=new G4RotationMatrix();
00848           Rot->rotateX(-pi/2 * rad);
00849           //Rot->rotateZ(pi * rad);
00850           //Rot->rotateZ(- ( pi/2 - fabs(angle)/2 ) * rad);
00851           TargetPos -= zHalfAngle *  ( length/2 + BDSGlobalConstants::Instance()->GetLengthSafety()/2 ) ;
00852           TargetPos+=G4ThreeVector(-rho,0,0);
00853           //TargetPos=G4ThreeVector(0,0,rho);
00854           //if(angle < 0)
00855           //{
00856           //  Rot->rotateZ(pi);
00857           //  TargetPos=G4ThreeVector(rho,0,0);
00858           //}
00859           rotateComponent=Rot;
00860         }
00861         */
00862 
00863 #ifdef DEBUG
00864         G4cout<<"ALIGNING COMPONENT..."<< G4endl;
00865 #endif  
00866         // Align Component - most cases does nothing. 
00867         // Currently only used for BDSElement   
00868         // For other elements stores global position and rotation,
00869         // needed for BDSOutline
00870         BDSBeamline::Instance()->currentItem()->AlignComponent(//TargetPos,
00871                                  rlast,
00872                                  rotateComponent,
00873                                  globalRotation,
00874                                  rtot,
00875                                  rlast,
00876                                  localX,
00877                                  localY,
00878                                  localZ);
00879 
00880 #ifdef DEBUG
00881         G4cout<<"Placing PHYSICAL COMPONENT..."<< G4endl;
00882 #endif  
00883 
00884         G4cout << "BDSDetectorConstruction : rotateComponent = " << *rotateComponent << G4endl;
00885         G4cout << "BDSDetectorConstruction : TargetPos        = " << TargetPos << G4endl;
00886         G4PVPlacement* PhysiComponentPlace = 
00887           new G4PVPlacement(
00888                             rotateComponent,  // its rotation
00889                             TargetPos,        // its position
00890                             LocalName,        // its name
00891                             LocalLogVol,      // its logical volume
00892                             physiWorld,       // its mother  volume
00893                             false,            // no boolean operation
00894                             nCopy,            // copy number
00895                             BDSGlobalConstants::Instance()->GetCheckOverlaps());              //overlap checking
00896 
00897           fPhysicalVolumeVector.push_back(PhysiComponentPlace);
00898           vector<G4VPhysicalVolume*> MultiplePhysicalVolumes = BDSBeamline::Instance()->currentItem()->GetMultiplePhysicalVolumes();
00899           for (unsigned int i=0;i<MultiplePhysicalVolumes.size(); i++) fPhysicalVolumeVector.push_back(MultiplePhysicalVolumes.at(i));
00900                                             
00901 
00902           
00903 #ifdef DEBUG 
00904         G4cout << "Volume name: " << LocalName << G4endl;
00905 #endif
00906         if(BDSGlobalConstants::Instance()->GetRefVolume()+"_phys"==LocalName && 
00907            BDSGlobalConstants::Instance()->GetRefCopyNo()==nCopy){
00908 #ifdef DEBUG 
00909           G4cout << "Setting new transform" <<G4endl;
00910 #endif
00911           G4AffineTransform tf(globalRotation,TargetPos-G4ThreeVector(0,0,length/2));
00912           BDSGlobalConstants::Instance()->SetRefTransform(tf);
00913         }
00914 
00915         BDSBeamline::Instance()->currentItem()->PrepareField(PhysiComponentPlace);
00916 
00917         if(use_graphics)
00918           {
00919                     BDSBeamline::Instance()->currentItem()->GetVisAttributes()->SetVisibility(true);
00920 #ifdef DEBUG
00921                     BDSBeamline::Instance()->currentItem()->GetMarkerLogicalVolume()->
00922                       SetVisAttributes(BDSBeamline::Instance()->currentItem()->GetVisAttributes());
00923 #endif
00924           }
00925         
00926       }
00927     }
00928 
00929   // construct tunnel
00930   for(it = beamline_list.begin();it!=beamline_list.end();it++)
00931     {
00932       
00933       if((*it).type==_TUNNEL ) {
00934         G4cout<<"BUILDING TUNNEL : "<<(*it).l<<"  "<<(*it).name<<G4endl;
00935         
00936         G4String gFormat="",  GFile="";
00937         G4String geometry = (*it).geometryFile;
00938 
00939         // get geometry format and file
00940         G4int pos = geometry.find(":");
00941         
00942         if(pos<0) { 
00943           G4cerr<<"WARNING: invalid geometry reference format : "<<geometry<<endl;
00944           gFormat="none";
00945         }
00946         
00947         else {
00948           gFormat = geometry.substr(0,pos);
00949           GFile = geometry.substr(pos+1,geometry.length() - pos); 
00950         }
00951         
00952         G4cout<<"placing components\n: geometry format - "<<gFormat<<G4endl<<
00953           "file - "<<GFile<<G4endl;
00954         
00955         if(gFormat=="gmad") {
00956          
00957           GGmadDriver ggmad(GFile);
00958           ggmad.Construct(logicWorld);
00959           
00960         } else  G4cerr<<"Tunnel won't be build! "<<endl;
00961       }
00962       
00963     }
00964   
00965   // free the parser list
00966   beamline_list.clear();
00967 
00968   if(verbose || debug) G4cout<<"end placement, size="<<BDSBeamline::Instance()->size()<<G4endl;
00969   
00970   if(verbose || debug) G4cout<<"Detector Construction done"<<G4endl; 
00971 
00972 #ifdef DEBUG 
00973     G4cout << *(G4Material::GetMaterialTable()) << G4endl;
00974 #endif
00975 
00976     if(verbose || debug) G4cout<<"Finished listing materials, returning physiWorld"<<G4endl; 
00977 
00978     return physiWorld;
00979 }
00980 
00981 //=================================================================
00982  
00983 void BDSDetectorConstruction::SetMagField(const G4double fieldValue){
00984   
00985   G4FieldManager* fieldMgr =
00986     G4TransportationManager::GetTransportationManager()->GetFieldManager();
00987   magField = new G4UniformMagField(G4ThreeVector(0.,fieldValue,0.));  
00988   fieldMgr->SetDetectorField(magField);
00989   fieldMgr->CreateChordFinder(magField);
00990 }
00991 
00992 //=================================================================
00993   
00994 void BDSDetectorConstruction::UpdateGeometry()
00995 {
00996   G4RunManager::GetRunManager()->DefineWorldVolume(ConstructBDS(beamline_list));
00997 }
00998 
00999 //=================================================================
01000 BDSDetectorConstruction::~BDSDetectorConstruction()
01001 { 
01002   LogVolCount->clear();
01003   delete LogVolCount;
01004 
01005   LogVolMap::iterator iter;
01006   for(iter=LogVol->begin();iter!=LogVol->end();iter++){
01007     delete (*iter).second;
01008   }
01009   LogVol->clear();
01010   delete LogVol;
01011 
01012   delete precisionRegion;
01013   gFlashRegion.clear();
01014 }
01015 
01016 //=================================================================
01017 
01018 G4double* BDSDetectorConstruction::GetWorldSize(){
01019   int s=3;
01020   G4double* ret = new G4double[s];
01021   for(int i=0; i<s; i++){
01022     ret[i]=itsWorldSize[i];
01023   }
01024   return ret;
01025 }
01026 
01027 G4double BDSDetectorConstruction::GetWorldSizeX(){
01028   return itsWorldSize[0];
01029 }
01030 
01031 G4double BDSDetectorConstruction::GetWorldSizeY(){
01032   return itsWorldSize[1];
01033 }
01034 
01035 G4double BDSDetectorConstruction::GetWorldSizeZ(){
01036   return itsWorldSize[2];
01037 }
01038 
01039 void BDSDetectorConstruction::SetWorldSize(G4double* val){
01040   int sExpected = 3;
01041   int s=sizeof(val)/sizeof(val[0]);
01042   if(s!=sExpected){
01043     std::cerr << "Error: BDSDetectorConstruction::SetWorldSize(G4double*) expects an array of size " << sExpected << ". Exiting." << std::endl;
01044     exit(1);
01045   }
01046   for(int i=0; i<s; i++){
01047     itsWorldSize[i]=val[i];
01048   }
01049 }
01050 
01051 void BDSDetectorConstruction::SetWorldSizeX(G4double val){
01052   itsWorldSize[0]=val;
01053 }
01054 
01055 void BDSDetectorConstruction::SetWorldSizeY(G4double val){
01056   itsWorldSize[1]=val;
01057 }
01058 
01059 void BDSDetectorConstruction::SetWorldSizeZ(G4double val){
01060   itsWorldSize[2]=val;
01061 }