MagneticField/src/FieldBuilder.cc

0001 #include "MagneticField/Engine/interface/MagneticField.h"
0002 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0003
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005
0006 #include "SimG4Core/MagneticField/interface/CMSFieldManager.h"
0007 #include "SimG4Core/MagneticField/interface/Field.h"
0008 #include "SimG4Core/MagneticField/interface/FieldBuilder.h"
0009 #include "SimG4Core/MagneticField/interface/FieldStepper.h"
0010 #include "SimG4Core/MagneticField/interface/MonopoleEquation.h"
0011
0012 #include <CLHEP/Units/SystemOfUnits.h>
0013 #include "G4ChordFinder.hh"
0014 #include "G4ClassicalRK4.hh"
0015 #include "G4FieldManager.hh"
0016 #include "G4LogicalVolumeStore.hh"
0017 #include "G4Mag_UsualEqRhs.hh"
0018 #include "G4TMagFieldEquation.hh"
0019 #include "CMSTMagFieldEquation.h"
0020 #include "G4PropagatorInField.hh"
0021
0022 using namespace sim;
0023
0024 FieldBuilder::FieldBuilder(const MagneticField *f, const edm::ParameterSet &p) : theTopVolume(nullptr), thePSet(p) {
0025   theDelta = p.getParameter<double>("delta") * CLHEP::mm;
0026   theField = new Field(f, theDelta);
0027   theFieldEquation = nullptr;
0028 }
0029
0030 FieldBuilder::~FieldBuilder() {}
0031
0032 void FieldBuilder::build(CMSFieldManager *fM, G4PropagatorInField *fP) {
0033   edm::ParameterSet thePSetForGMFM = thePSet.getParameter<edm::ParameterSet>("ConfGlobalMFM");
0034   std::string volName = thePSetForGMFM.getParameter<std::string>("Volume");
0035   edm::ParameterSet volPSet = thePSetForGMFM.getParameter<edm::ParameterSet>(volName);
0036
0037   configureForVolume(volName, volPSet, fM, fP);
0038
0039   edm::LogVerbatim("SimG4CoreMagneticField") << " FieldBuilder::build: Global magnetic field is used";
0040 }
0041
0042 void FieldBuilder::configureForVolume(const std::string &volName,
0043                                       edm::ParameterSet &volPSet,
0044                                       CMSFieldManager *fM,
0045                                       G4PropagatorInField *fP) {
0046   G4LogicalVolumeStore *theStore = G4LogicalVolumeStore::GetInstance();
0047   for (auto vol : *theStore) {
0048     if ((std::string)vol->GetName() == volName) {
0049       theTopVolume = vol;
0050       break;
0051     }
0052   }
0053
0054   std::string fieldType = volPSet.getParameter<std::string>("Type");
0055   std::string stepper = volPSet.getParameter<std::string>("Stepper");
0056
0057   edm::ParameterSet stpPSet = volPSet.getParameter<edm::ParameterSet>("StepperParam");
0058   double minStep = stpPSet.getParameter<double>("MinStep") * CLHEP::mm;
0059
0060   if (stepper == "CMSTDormandPrince45") {
0061     theFieldEquation = new CMSTMagFieldEquation<sim::Field>(theField);
0062   } else if (stepper == "G4TDormandPrince45") {
0063     theFieldEquation = new G4TMagFieldEquation<sim::Field>(theField);
0064   } else {
0065     theFieldEquation = new G4Mag_UsualEqRhs(theField);
0066   }
0067
0068   FieldStepper *dStepper = new FieldStepper(theFieldEquation, theDelta, stepper);
0069   G4ChordFinder *cf = new G4ChordFinder(theField, minStep, dStepper);
0070
0071   MonopoleEquation *monopoleEquation = new MonopoleEquation(theField);
0072   G4MagIntegratorStepper *mStepper = new G4ClassicalRK4(monopoleEquation, 8);
0073   G4ChordFinder *cfmon = new G4ChordFinder(theField, minStep, mStepper);
0074
0075   fM->InitialiseForVolume(stpPSet, theField, cf, cfmon, volName, fieldType, stepper, theDelta, fP);
0076 }