loon: GeoMediumMap.cxx Source File

00001 
00002 //
00003 // GeoMediumMap
00004 //
00005 // GeoMediumMap is a helper class to build mediums and associations
00006 // between detector components and mediums.
00007 //
00008 // Author:  S. Kasahara 11/07
00010 
00011 #include <cassert>
00012 #include <stdlib.h>
00013 
00014 #include <TGeoManager.h>
00015 #include <TList.h>
00016 
00017 #include "Util/UtilString.h"
00018 #include "Util/UtilMCFlag.h"
00019 using namespace UtilMCFlag;
00020 #include "Conventions/MinosMaterial.h"
00021 #include "GeoGeometry/GeoMediumMap.h"
00022 #include "MessageService/MsgService.h"
00023 #include "UgliGeometry/UgliLoanPool.h"
00024 #include "GeoGeometry/GeoGeometry.h"
00025 #include "GeoGeometry/GeoMedium.h"
00026 #include "Registry/Registry.h"
00027 
00028 ClassImp(GeoMediumMap)
00029 
00030 CVSID("$Id: GeoMediumMap.cxx,v 1.14 2010/07/28 21:53:45 rhatcher Exp $");
00031 
00032 //_____________________________________________________________________________
00033 GeoMediumMap::GeoMediumMap(GeoGeometry* geo,const VldContext& vldc) : 
00034                            fGeoGeometry(geo),fVldContext(vldc) {
00035   // Normal constructor. Protected. Invoked by GeoGeometry creator.
00036 
00037   MSG("Geo",Msg::kDebug) << "GeoMediumMap normal ctor @ " << this 
00038                          << " w/vldc " << vldc << " fGeoGeometry "
00039                          << fGeoGeometry << endl;
00040 
00041   UpdateGlobalManager();
00042   
00043   BuildMediums();
00044   BuildMediumMap();
00045 
00046 }
00047 
00048 //_____________________________________________________________________________
00049 void GeoMediumMap::UpdateGlobalManager() const {
00050 
00051   if ( fGeoGeometry ) fGeoGeometry -> UpdateGlobalManager();
00052   else gGeoManager = 0;
00053   
00054 }
00055 
00056 //_____________________________________________________________________________
00057 GeoMedium* GeoMediumMap::GetMedium(Geo::EDetComponent detcomp) const {
00058   // Get medium corresponding to detector component
00059 
00060   MediumMapConstItr citr = fMediumMap.find(detcomp);
00061   if ( citr == fMediumMap.end() ) {
00062     MSG("Geo",Msg::kWarning) << "GeoMediumMap::GetMedium unable to find "
00063                              << "medium for detcomponent " 
00064                              << Geo::AsString(detcomp) << "." << endl;
00065     return (GeoMedium*)0;
00066   }
00067   
00068   return const_cast<GeoMedium*>(citr->second);
00069   
00070 }
00071 
00072 //_____________________________________________________________________________
00073 Geo::EFldStrength GeoMediumMap::GetFldStrength(Geo::EDetComponent detcomp) 
00074                                                                     const {
00075   // Get field strength corresponding to detector component
00076 
00077   FldStrengthMapConstItr citr = fFldStrengthMap.find(detcomp);
00078   if ( citr == fFldStrengthMap.end() ) {
00079     MSG("Geo",Msg::kWarning) << "GeoMediumMap::GetFldStrength unable to find "
00080                              << "field strength for detcomponent " 
00081                              << Geo::AsString(detcomp) << "." << endl;
00082     return Geo::kUnknownFldStrength;
00083   }
00084   
00085   return citr->second;
00086   
00087 }
00088 
00089 //_____________________________________________________________________________
00090 void GeoMediumMap::BuildMediums() {
00091    // Private. Build this geometry's materials.
00092 
00093    UpdateGlobalManager();
00094 
00095    MSG("Geo",Msg::kDebug) << "GeoMediumMap::BuildMediums " << endl;
00096    assert(gGeoManager);
00097 
00098 //-----------List of Materials, Mixtures and Tracking Media--------------
00099    // Arguments for Medium definitions: 
00100    //            medName,medId,matId,
00101    //            isVol  (Sensitive volume flag, set 0 if not sensitive)
00102    //            iField (User defined magnetic field, set
00103    //                    0, if no magnetic field
00104    //                    1, tracking performed with GRKUTA 
00105    //                       (use for strongly inhomogeneous field)
00106    //                    2, tracking performed with GHELIX
00107    //                       (use for inhomogeneous field)
00108    //                    3, tracking performed with GHELX3
00109    //                       (use for uniform magnetic field along z axis)
00110    //            fieldM (Max field value (in kGauss)),
00111    //            tMaxFd (Max angle due to field deflection in one step (deg)),
00112    //            steMax (Max step allowed (cm)),
00113    //            deeMax (Max fractional energy loss in a step),
00114    //            epsil  (Boundary crossing precision (cm)),
00115    //            stMin  (Minimum step due to energy loss, multiple scattering,
00116    //                    Cerenkov or B-fields (cm))
00117 
00118 
00119    // Default max field values for field activated mediums 
00120    Float_t fieldM_iron = 20.; // kGauss, in iron collar or steel plane
00121    Float_t fieldM      = 6.;  // kGauss, outside of iron collar or steel plane 
00122 
00123    // Sensitive volume flags
00124    Int_t isSens    = 1;
00125    Int_t isNotSens = 0;
00126    
00127    // Default boundary crossing precision
00128    Int_t medId = 1; // id assigned to first medium, incr with each new medium
00129    Int_t matId = 1; // id assigned to first material, incr with each new mat
00130    
00131    
00132    // **********  ALUMINIUM **************
00133    //MinosMaterial::StdMaterial_t matAl = MinosMaterial::eAluminium;
00134    TGeoMaterial *mat1 = new TGeoMaterial("ALUMINIUM",26.98, 13, 2.7); 
00135    mat1->SetUniqueID(matId);
00136    medId = CreateMediumsforMaterial("ALUM",medId,matId,isNotSens,fieldM);
00137 
00138    medId++;
00139    matId++;
00140    
00141 
00142    // ************** IRON ***************
00143    //MinosMaterial::StdMaterial_t matFe = MinosMaterial::eIron;
00144    TGeoMaterial *mat2 = new TGeoMaterial("IRON",55.85, 26, 7.87);
00145    mat2->SetUniqueID(matId);
00146    medId = CreateMediumsforMaterial("IRON",medId, matId,isNotSens,fieldM_iron);
00147 
00148    medId++;
00149    matId++;
00150   
00151 
00152    // ************** AIR ****************
00153    TGeoMixture *mat3 = new TGeoMixture("AIR",2,0.1205000E-02);
00154    mat3->SetUniqueID(matId);
00155    mat3->DefineElement(0,14.01,7,0.7);
00156    mat3->DefineElement(1,16.00,8,0.3);
00157    medId = CreateMediumsforMaterial("AIR",medId, matId,isNotSens,fieldM);
00158 
00159    medId++;
00160    matId++;
00161 
00162 
00163    // ************** CONCRETE ****************
00164    //MinosMaterial::StdMaterial_t matConcrete = MinosMaterial::eConcrete;
00165    TGeoMixture *mat4 = new TGeoMixture("CONCRETE",6, 2.5); 
00166    mat4->SetUniqueID(matId);
00167    mat4->DefineElement(0,16,8,0.53);
00168    mat4->DefineElement(1,28.09,14,0.33);
00169    mat4->DefineElement(2,40.078,20,0.6300000E-01);
00170    mat4->DefineElement(3,22.99,11,0.1500000E-01);
00171    mat4->DefineElement(4,55.85,26,0.2000000E-01);
00172    mat4->DefineElement(5,26.98,13,0.4200000E-01);
00173    medId = CreateMediumsforMaterial("CONCRETE",medId, matId,isNotSens,fieldM);
00174 
00175    medId++;
00176    matId++;
00177    
00178 
00179    // *************** ROCK *****************
00180    // ROCK once used for far detector rock. Density should be closer to 2.85.
00181    // Replaced with GREENSTONE on 2007/09/13.
00182    TGeoMixture *mat5 = new TGeoMixture("ROCK",6, 2.5); 
00183    mat5->SetUniqueID(matId);
00184    mat5->DefineElement(0,16,8,0.53);  // oxygen, 53%
00185    mat5->DefineElement(1,28.09,14,0.33); // silicon 33%
00186    mat5->DefineElement(2,40.078,20,0.6300000E-01); // calcium, 6.3%
00187    mat5->DefineElement(3,22.99,11,0.1500000E-01);  // sodium, 1.5%
00188    mat5->DefineElement(4,55.85,26,0.2000000E-01);  // iron, 2%
00189    mat5->DefineElement(5,26.98,13,0.4200000E-01);  // aluminum, 4.2%
00190    medId = CreateMediumsforMaterial("ROCK",medId, matId,isNotSens,fieldM);
00191 
00192    medId++;
00193    matId++;
00194    
00195 
00196    // *************** PSTYRENE SCINT  *****************
00197    //   MinosMaterial::StdMaterial_t matScint=MinosMaterial::ePolystyreneMinos;
00198    // Pure pstyrene is C6H5CH.
00199    TGeoMixture *mat6 = new TGeoMixture("PSTYRENE SCINT",2, 1.0457); 
00200    mat6->SetUniqueID(matId);
00201    mat6->DefineElement(0,1.00794,1,0.7742105E-01);
00202    mat6->DefineElement(1,12.011,6,0.9225789);
00203    // Sensitve medium pstyrene scint, used for scint strips
00204    medId = CreateMediumsforMaterial("PSTYRENE SCINT",medId,matId,isSens,
00205                                      fieldM);
00206 
00207    medId++;
00208    matId++;
00209 
00210 
00211 
00212    // *************** COEX TIO2 PSTYRENE  *****************
00213    //MinosMaterial::StdMaterial_t matTiO2=MinosMaterial::ePolystyreneMinos;
00214    // Pure pstyrene is C6H5CH. Coex Pstyrene is  10% TiO2 by weight
00215    TGeoMixture *mat7 = new TGeoMixture("COEX TIO2 PSTYRENE",4, 1.0457); 
00216    mat7->SetUniqueID(matId);
00217    mat7->DefineElement(0,1.00794,1,0.6967895E-01);
00218    mat7->DefineElement(1,12.011,6,0.8303211);
00219    mat7->DefineElement(2,47.93,22,0.5996587E-01);
00220    mat7->DefineElement(3,15.9994,8,0.4003413E-01);
00221    // Coex TiO2 used for white coating on surface of scint strips
00222    medId = CreateMediumsforMaterial("COEX TIO2 PSTYRENE",medId,matId,
00223                                      isNotSens,fieldM);
00224 
00225    medId++;
00226    matId++;
00227    
00228 
00229    // ******************** FARCOIL **********************
00230    TGeoMixture *mat8 = new TGeoMixture("FARCOIL",5,   2.34881    );
00231    mat8->SetUniqueID(matId);
00232    mat8->DefineElement(0,1.00794,1,0.1281415E-01);
00233    mat8->DefineElement(1,12.0107,6,0.6021570E-01);
00234    mat8->DefineElement(2,14.00674,7,0.1974151E-03);
00235    mat8->DefineElement(3,15.9994,8,0.2134690E-01);
00236    mat8->DefineElement(4,63.546,29,0.9054258);
00237    medId = CreateMediumsforMaterial("FARCOIL",medId, matId,isNotSens,fieldM);
00238 
00239    medId++;
00240    matId++;
00241    
00242 
00243    // ******************* DOLOMITE ********************
00244    // DOLOMITE once used for near det rock. Density should be closer to 2.5?
00245    // Replaced with NEARBEDROCK on 2007/09/13.
00246    TGeoMixture *mat9 = new TGeoMixture("DOLOMITE",5,   2.850000    );
00247    mat9->SetUniqueID(matId);
00248    mat9->DefineElement(0,40.078,20,0.2161123); // calcium 21.6%
00249    mat9->DefineElement(1,24.305,12,0.1310597); // magnesium 13.1%
00250    mat9->DefineElement(2,12.011,6,0.1295336);  // carbon 13.0%
00251    mat9->DefineElement(3,15.9994,8,0.5226618); // oxygen 52.3%
00252    mat9->DefineElement(4,1.0079,1,0.6326329E-03);  // hydrogen 0.06%
00253    medId = CreateMediumsforMaterial("DOLOMITE",medId, matId,isNotSens,fieldM);
00254 
00255    medId++;
00256    matId++;
00257 
00258    
00259    // ******************** FEPL **********************
00260    // FEPL MINOS Steel plane mixture  
00261    //MinosMaterial::StdMaterial_t matSteel=MinosMaterial::eIronMinos;
00262    Double_t densFEPL = 7.847;
00263    if ( fVldContext.GetDetector() == Detector::kCalDet ) densFEPL = 7.85;
00264    TGeoMixture *mat10 = new TGeoMixture("FEPL",17, densFEPL);
00265    mat10->SetUniqueID(matId);
00266    mat10->DefineElement(0,54,26,0.5556792E-01);
00267    mat10->DefineElement(1,56,26,0.9046049);
00268    mat10->DefineElement(2,57,26,0.2126433E-01);
00269    mat10->DefineElement(3,58,26,0.2879539E-02);
00270    mat10->DefineElement(4,54.94,25,0.1016047E-01);
00271    mat10->DefineElement(5,28.09,14,0.2917362E-02);
00272    mat10->DefineElement(6,12.01,6,0.1810777E-02);
00273    mat10->DefineElement(7,58.7,28,0.2112573E-03);
00274    mat10->DefineElement(8,52,24,0.1810777E-03);
00275    mat10->DefineElement(9,30.97,15,0.1508981E-03);
00276    mat10->DefineElement(10,32.06,16,0.8047896E-04);
00277    mat10->DefineElement(11,63.55,29,0.7041910E-04);
00278    mat10->DefineElement(12,26.98,13,0.3017961E-04);
00279    mat10->DefineElement(13,14.01,7,0.3017961E-04);
00280    mat10->DefineElement(14,118.69,50,0.2011974E-04);
00281    mat10->DefineElement(15,47.9,22,0.1005987E-04);
00282    mat10->DefineElement(16,50.94,23,0.1005987E-04);
00283    medId = CreateMediumsforMaterial("FEPL",medId, matId,isNotSens,fieldM_iron);
00284    
00285    medId++;
00286    matId++;
00287    
00288 
00289    // ****************** WATER *******************
00290    MinosMaterial::StdMaterial_t matWater=MinosMaterial::eWater;
00291    TGeoMixture* mat11 = new TGeoMixture("WATER",2,
00292                                    MinosMaterial::Density(matWater) ); // 1.0
00293    mat11 -> SetUniqueID(matId);
00294    mat11 -> DefineElement(0,1.00794,1,0.1118984);
00295    mat11 -> DefineElement(1,15.9994,8,0.8881017);
00296    medId = CreateMediumsforMaterial("WATER",medId,matId,isNotSens,fieldM);
00297 
00298    medId++;
00299    matId++;
00300    
00301 
00302    // *************** GREENSTONE *****************
00303    // GREENSTONE to be used for far detector rock. 
00304    // Composition from K. Ruddick's compilation of data from regional rock 
00305    // samples (c. 1970).  Private communication of 6/21/2007.
00306    // New density of 2.760 introduced on 9/30/2009.  This more precise
00307    // value is from M. Strait's MINOS DocDB #5022.
00308    // Replaces density of 2.8 from K. Ruddick's PDK Note 435, 1990.
00309    TGeoMixture *mat12 = new TGeoMixture("GREENSTONE",13, 2.760); 
00310    mat12->SetUniqueID(matId);
00311    mat12->DefineElement( 0,15.999, 8,0.4586);  // oxygen,    45.86%
00312    mat12->DefineElement( 1,28.086,14,0.2370);  // silicon,   23.70%
00313    mat12->DefineElement( 2,55.845,26,0.0856);  // iron,       8.56%
00314    mat12->DefineElement( 3,26.982,13,0.0797);  // aluminum,   7.97%
00315    mat12->DefineElement( 4,40.078,20,0.0641);  // calcium,    6.41%
00316    mat12->DefineElement( 5,24.305,12,0.0395);  // magnesium,  3.95%
00317    mat12->DefineElement( 6,22.990,11,0.0190);  // sodium,     1.90%
00318    mat12->DefineElement( 7,47.867,22,0.0064);  // titanium,   0.64%
00319    mat12->DefineElement( 8,39.098,19,0.0034);  // potassium,  0.34%
00320    mat12->DefineElement( 9, 1.008, 1,0.0026);  // hydrogen,   0.26%
00321    mat12->DefineElement(10,12.011, 6,0.0020);  // carbon,     0.20%
00322    mat12->DefineElement(11,54.938,25,0.0015);  // manganese,  0.15%
00323    mat12->DefineElement(12,30.974,15,0.0006);  // phosphorus, 0.06%
00324    medId = CreateMediumsforMaterial("GREENSTONE",medId,matId,isNotSens,fieldM);
00325 
00326    medId++;
00327    matId++;
00328    
00329 
00330    // *************** NEARBEDROCK *****************
00331    // NEARBEDROCK to be used for near detector bedrock. 
00332    // Composition (CaMg(CO_3)_2 with 8.0% moisture content by weight) 
00333    // and density from S. Kasahara's NuMI note 2777 and references 
00334    // therein.
00335 
00336    TGeoMixture *mat13 = new TGeoMixture("NEARBEDROCK",5, 2.43); 
00337    mat13->SetUniqueID(matId);
00338    mat13->DefineElement( 0,15.999, 8,0.550);  // oxygen,    55.0%
00339    mat13->DefineElement( 1,40.078,20,0.200);  // calcium,   20.0%
00340    mat13->DefineElement( 2,24.305,12,0.121);  // magnesium, 12.1%
00341    mat13->DefineElement( 3,12.011, 6,0.120);  // carbon,    12.0%
00342    mat13->DefineElement( 4, 1.008, 1,0.009);  // hydrogen,   0.9%
00343    medId = CreateMediumsforMaterial("NEARBEDROCK",medId,matId,isNotSens,
00344                                      fieldM);
00345 
00346    medId++;
00347    matId++;
00348    
00349 
00350    // ************** LEAD ***************
00351    // Lead is used as medium of SOLO (SOudan LOw Background Counting) 
00352    // and "Reeve's Castle" detectors located upstream of the main detector 
00353    // in far detector cavern. Atomic mass and density of lead are taken 
00354    // from 2007 PDG.
00355    TGeoMaterial *mat14 = new TGeoMaterial("LEAD",207.2,82,11.35);
00356    mat14->SetUniqueID(matId);
00357    medId = CreateMediumsforMaterial("LEAD",medId, matId,isNotSens,fieldM);
00358 
00359    medId++;
00360    matId++;
00361 
00362    // At the end of the build of all media, build a duplicate set
00363    // for use with the GeoSwimmer.  These media will have their
00364    // physics process flags and cut thresholds adjusted to suit
00365    // continuous energy loss when using GeoSwimmer with a concrete
00366    // VMC.
00367    Int_t nmed = medId - 1;
00368    TList* medList = gGeoManager -> GetListOfMedia();
00369    for ( int imed = 0; imed < nmed; imed++ ) {
00370      GeoMedium* medium = dynamic_cast<GeoMedium*>(medList->At(imed));
00371      std::string medName = medium->GetName();
00372      Int_t npos = medName.find("_B");
00373      std::string medBaseName = medName.substr(0,npos);
00374      std::string medSwimName 
00375        = medBaseName+"_SWIM"+medName.substr(npos,medName.length()-npos);
00376 
00377      // TGeoMedium copy ctor is protected or we wouldn't have to do this
00378      TGeoMaterial* material = medium -> GetMaterial();
00379      Double_t param[20];
00380      for ( int ipar = 0; ipar < 20; ipar++ ) {
00381        param[ipar] = medium -> GetParam(ipar);
00382      }
00383      param[6] = 0.01; // epsil (Boundary crossing precision (cm)) always 0.01. 
00384      GeoMedium* swimMed= new GeoMedium(medSwimName.c_str(),medId++,material,
00385                                        param);
00386      // Override default physic process and cut flag values
00387      for ( EProcess ip = EProcess(0); ip < kNProcess; ip = EProcess(ip+1) ) {
00388        if ( ip == kLOSS ) swimMed -> SetProcess(kLOSS,4);
00389        else swimMed -> SetProcess(ip,0);
00390      }
00391      for ( ECut ic = ECut(0); ic < kNCut; ic = ECut(ic+1) ) {
00392        if ( ic == kTOFMAX ) continue;
00393        if ( ic == kCUTMUO ) swimMed -> SetCut(kCUTMUO,0.0001);
00394        else swimMed -> SetCut(ic,1000.);
00395      }
00396      ConfigureMedium(swimMed); // configure swimmed if specified by user
00397    }
00398    
00399    
00400    return;
00401 
00402 }
00403 
00404 //_____________________________________________________________________________
00405 void GeoMediumMap::Print(Option_t* /* option */) const {
00406   // Print medium name, swim method associated with each detector
00407   // component
00408 
00409   MediumMapConstItr citr;
00410   cout << setiosflags(ios::left)
00411        << setw(15) << "DetComponent"
00412        << setw(30) << "Medium"
00413        << setw(15) << "FldStrength"
00414        << setw(15) << "SwimMethod"
00415        << endl;
00416   cout << setiosflags(ios::left) 
00417        << setfill('-') << setw(12) << "-" 
00418        << setfill(' ') << setw(3) << " "
00419        << setfill('-') << setw(6) << "-" 
00420        << setfill(' ') << setw(24) << " "
00421        << setfill('-') << setw(11) << "-" 
00422        << setfill(' ') << setw(4) << " "
00423        << setfill('-') << setw(10) << "-" 
00424        << setfill(' ') << setw(5) << " "
00425        << endl;
00426   
00427   for (citr = fMediumMap.begin(); citr != fMediumMap.end(); citr++ ) {
00428     Geo::EFldStrength fldstrength = GetFldStrength(citr->first);
00429     Geo::ESwimMethod swimmethod = GetSwimMethod(fldstrength);
00430     
00431     cout << setiosflags(ios::left) << setfill(' ')
00432          << setw(15) << Geo::AsString(citr->first)
00433          << setw(30) << citr->second->GetName()
00434          << setw(15) << Geo::AsString(fldstrength)
00435          << setw(15) << Geo::AsString(swimmethod) 
00436          << endl;
00437   }
00438   
00439 }
00440 
00441 //_____________________________________________________________________________
00442 int GeoMediumMap::CreateMediumsforMaterial(std::string medbasename,
00443                      int medId, int matId, int isVol, float fieldM) {
00444    // Create 3 Mediums for given Material w/Id matId:
00445    //   1)medbasename_B0 =>No field
00446    //   2)medbasename_B1 =>Swim with Runge-Kutta method (high gradient field)
00447    //   3)medbasename_B2 =>Swim with Stepwise-Helix method (low gradient field)
00448    // Mediums will be assigned Id's incrementally beginning with input medId.
00449    // Argument isVol and fieldM are used to define the medium sensitivity
00450    // and max BField respectively.
00451    // Return value is the last medId assigned to a created medium.
00452 
00453    // TGeoMedium is defined with arguments:
00454    //            medName,medId,matId,
00455    //            isVol  (Sensitive volume flag, set 0 if not sensitive)
00456    //            iField (User defined magnetic field, set
00457    //                    0, if no magnetic field
00458    //                    1, tracking performed with GRKUTA 
00459    //                       (use for strongly inhomogeneous field)
00460    //                    2, tracking performed with GHELIX
00461    //                       (use for inhomogeneous field)
00462    //                    3, tracking performed with GHELX3
00463    //                       (use for uniform magnetic field along z axis)
00464    //            fieldM (Max field value (in kGauss)),
00465    //            tMaxFd (Max angle due to field deflection in one step (deg)),
00466    //            steMax (Max step allowed (cm)),
00467    //            deeMax (Max fractional energy loss in a step),
00468    //            epsil  (Boundary crossing precision (cm)),
00469    //            stMin  (Minimum step due to energy loss, multiple scattering,
00470    //                    Cerenkov or B-fields (cm))
00471 
00472    // The following 4 values will be automatically computed by Geant3 when 
00473    // Auto configuration parameter set to 1 (default).
00474    // Even if the user sets Auto parameter to 0, negative values will still be
00475    // computed automatically.
00476    // See CONS210 section of Geant3 manual. 
00477    Float_t tMaxFd = -1.;   
00478    Float_t steMax = -1.;
00479    Float_t deeMax = -1.;
00480    Float_t  stMin = -1.;
00481 
00482    Float_t epsil       = +0.0001;
00483    Float_t field_epsil = +0.01; // used for all magnetized mediums
00484 
00485     // Non-magnetized medium 
00486    std::string medname;
00487    GeoMedium* med = 0;
00488    
00489    medname  = medbasename + "_B" +UtilString::ToString<int>(Geo::kNoFieldSwim);
00490    med = new GeoMedium(medname.c_str(),medId,matId,isVol, Geo::kNoFieldSwim,
00491                                   fieldM,tMaxFd,steMax,deeMax,epsil,stMin);
00492    ConfigureMedium(med);
00493 
00494    medId++;
00495    
00496    // Magnetized aluminium medium treated as strongly inhomogeneous by swimmer
00497    medname = medbasename + "_B" + UtilString::ToString<int>(Geo::kRungeKutta);
00498    med = new GeoMedium(medname.c_str(),medId,matId,isVol,Geo::kRungeKutta,
00499                        fieldM,tMaxFd,steMax,deeMax,field_epsil,stMin);
00500    ConfigureMedium(med);
00501    
00502    medId++;
00503    
00504    // Magnetized aluminium medium treated as inhomogeneous by swimmer
00505    medname = medbasename + "_B" + UtilString::ToString<int>(Geo::kStepHelix);
00506    med = new GeoMedium(medname.c_str(),medId,matId,isVol,Geo::kStepHelix,
00507                  fieldM,tMaxFd,steMax,deeMax,field_epsil,stMin);
00508    ConfigureMedium(med);
00509    
00510    return medId;
00511    
00512 }
00513 
00514 //_____________________________________________________________________________
00515 void GeoMediumMap::ConfigureMedium(GeoMedium* med) {
00516   // Private method to check UgliLoanPool registry for user specified medium
00517   // definition values that override default.
00518   
00519    std::string medName[2];
00520   
00521    medName[1] = med->GetName(); // specific name
00522    medName[0] = medName[1].substr(0,medName[1].find("_B")); // base name
00523 
00524    // Check UgliLoanPool registry for user specified medium definition values
00525    // that override default.
00526    for ( int i = 0; i < 2; i++ ) {
00527      // base name first, specific name second so that specific name config
00528      // will override that of base name if both present and overlapping.
00529      Registry medConfig 
00530             = UgliLoanPool::Instance()->GetMediumRegistry(medName[i].c_str());
00531    
00532      Registry::RegistryKey keyItr(&medConfig);
00533      while ( const char* newKey = keyItr() ) {
00534        Double_t tmpd;
00535        medConfig.Get(newKey,tmpd);
00536        UtilMCFlag::EProcess process = UtilMCFlag::GetProcess(newKey);
00537        if ( process != UtilMCFlag::kUnknownProcess ) 
00538                                       med -> SetProcess(process,(int)tmpd);
00539        else {
00540          UtilMCFlag::ECut cut = UtilMCFlag::GetCut(newKey);
00541          if ( cut != UtilMCFlag::kUnknownCut ) med -> SetCut(cut,tmpd);
00542          else {
00543            UtilMCFlag::ETracking trk = UtilMCFlag::GetTracking(newKey);
00544            if ( trk != UtilMCFlag::kUnknownTracking ) 
00545              med -> SetTracking(trk,tmpd);
00546            else {
00547              MSG("Geo",Msg::kWarning) << "Unknown medium config key " 
00548                                       << newKey << " for medium " 
00549                                       << med->GetName() << ". Ignored!" 
00550                                       << endl;
00551            }
00552          }
00553        }
00554      }  
00555    }
00556    
00557 }
00558 
00559 //_____________________________________________________________________________
00560 void GeoMediumMap::BuildSwimMethodMap() {
00561   // Private method to build swim method map for given validity context
00562 
00563   // Fill map with default values
00564   fSwimMethodMap[Geo::kNoField]   = Geo::kNoFieldSwim;
00565   fSwimMethodMap[Geo::kLowField]  = Geo::kNoFieldSwim;
00566   fSwimMethodMap[Geo::kHighField] = Geo::kRungeKutta;
00567   
00568   // Override with UgliLoanPool configuration as appropriate
00569   const Registry& r = UgliLoanPool::Instance()->GetConfig();
00570   Int_t tmpi;
00571   if ( r.Get("SwimMethodLowField",tmpi) ) {
00572     Geo::ESwimMethod swimmethod = (Geo::ESwimMethod)tmpi;
00573     if ( swimmethod != fSwimMethodMap[Geo::kLowField] ) {
00574       MSG("Geo",Msg::kInfo) << "SwimMethod for components w/field strength "
00575                             << Geo::AsString(Geo::kLowField) << " set to " 
00576                             << Geo::AsString(swimmethod) 
00577                             << " overriding default " 
00578                             << Geo::AsString(fSwimMethodMap[Geo::kLowField]) 
00579                             << "." << endl;
00580       fSwimMethodMap[Geo::kLowField]  = swimmethod;
00581     }
00582   }
00583   
00584   if ( r.Get("SwimMethodHighField",tmpi) ) {
00585     Geo::ESwimMethod swimmethod = (Geo::ESwimMethod)tmpi;
00586     if ( swimmethod != fSwimMethodMap[Geo::kHighField] ) {
00587       MSG("Geo",Msg::kInfo) << "SwimMethod for components w/field strength "
00588                             << Geo::AsString(Geo::kHighField) << " set to " 
00589                             << Geo::AsString(swimmethod) 
00590                             << " overriding default " 
00591                             << Geo::AsString(fSwimMethodMap[Geo::kHighField]) 
00592                             << "." << endl;
00593       fSwimMethodMap[Geo::kHighField] = swimmethod;
00594     }
00595   }
00596   
00597 }
00598 
00599 //_____________________________________________________________________________
00600 bool GeoMediumMap::BuildMediumMap() {
00601   // Private method to build medium map for given validity context
00602 
00603   bool isOkay = true;
00604   Detector::Detector_t dettype = fVldContext.GetDetector();
00605 
00606   // Set up swim method 
00607   BuildSwimMethodMap();
00608   
00609   switch( dettype ) {
00610 
00611   case Detector::kFar:
00612 
00613     BuildFarMediumMap();
00614     break;
00615    
00616   case Detector::kNear:
00617 
00618     BuildNearMediumMap();
00619     break;
00620 
00621   case Detector::kCalDet:
00622 
00623     BuildCalMediumMap();
00624     break;
00625    
00626   default:
00627 
00628     BuildDefMediumMap();
00629     break;
00630 
00631   } // end of switch
00632 
00633   return isOkay;
00634   
00635 }
00636 
00637 //_____________________________________________________________________________
00638 void GeoMediumMap::BuildFarMediumMap() {
00639   // Private method to build medium map for the far detector.
00640 
00641   // The following regions are never likely to have a field implemented
00642   SetMedium(Geo::kMars,       "GREENSTONE",      Geo::kNoField); // rock 
00643   SetMedium(Geo::kLinr,       "GREENSTONE",      Geo::kNoField); // liner
00644   SetMedium(Geo::kHall,       "AIR",             Geo::kNoField); // hall
00645   SetMedium(Geo::kRCCoil,     "FARCOIL",         Geo::kNoField); // return coil
00646   SetMedium(Geo::kRCAir,      "AIR",             Geo::kNoField); 
00647   SetMedium(Geo::kVSStpScint, "PSTYRENE SCINT",  Geo::kNoField); // veto shield
00648   SetMedium(Geo::kVSStpTiO2,  "COEX TIO2 PSTYRENE", Geo::kNoField);
00649   SetMedium(Geo::kVSMdlAlSkin,"ALUM",               Geo::kNoField);
00650   SetMedium(Geo::kVSMdlAir,   "AIR",                Geo::kNoField);
00651   SetMedium(Geo::kVSPlnScint, "AIR",                Geo::kNoField); 
00652   SetMedium(Geo::kSoLo,       "LEAD",               Geo::kNoField);
00653   
00654   // The following regions in the main detector but outside the steel
00655   // and coil region will have no field or low field implemented
00656   SetMedium(Geo::kStpScint,   "PSTYRENE SCINT",    Geo::kLowField);
00657   SetMedium(Geo::kStpTiO2,    "COEX TIO2 PSTYRENE",Geo::kLowField);
00658   SetMedium(Geo::kMdlAlSkin,  "ALUM",              Geo::kLowField);
00659   SetMedium(Geo::kMdlAir,     "AIR",               Geo::kLowField);
00660   SetMedium(Geo::kPlnScint,   "AIR",               Geo::kLowField);
00661   SetMedium(Geo::kPlnAirGap,  "AIR",               Geo::kLowField);
00662 
00663   // The following regions in the steel and in the coil region have a
00664   // strong field implemented
00665   SetMedium(Geo::kPlnSteel,   "FEPL",         Geo::kHighField); // steel pln 
00666   SetMedium(Geo::kCRGapCoil,  "FARCOIL",      Geo::kHighField); // air gap CR
00667   SetMedium(Geo::kCRGapThroat,"AIR",          Geo::kHighField); 
00668   SetMedium(Geo::kCRGapFlange,"IRON",         Geo::kHighField);
00669   SetMedium(Geo::kCRGapBypass,"AIR",          Geo::kHighField);
00670   SetMedium(Geo::kCRStlCoil,  "FARCOIL",      Geo::kHighField); // steel pln CR
00671   SetMedium(Geo::kCRStlThroat,"AIR",          Geo::kHighField);
00672   SetMedium(Geo::kCRStlNeck,  "IRON",         Geo::kHighField);
00673   SetMedium(Geo::kCRStlHole,  "AIR",          Geo::kHighField);
00674   SetMedium(Geo::kCRStlDetail,"FEPL",         Geo::kHighField);
00675   SetMedium(Geo::kCRSctCoil,  "FARCOIL",      Geo::kHighField); // scint pln CR
00676   SetMedium(Geo::kCRSctThroat,"AIR",          Geo::kHighField);
00677   SetMedium(Geo::kCRSctFlange,"IRON",         Geo::kHighField);
00678   SetMedium(Geo::kCRSctBypass,"AIR",          Geo::kHighField);
00679   
00680 }
00681 
00682 //_____________________________________________________________________________
00683 void GeoMediumMap::BuildNearMediumMap() {
00684   // Private method to build medium map for the near detector.
00685 
00686   // The following regions are never likely to have a field implemented
00687   SetMedium(Geo::kMars,       "NEARBEDROCK",  Geo::kNoField); // rock 
00688   SetMedium(Geo::kLinr,       "NEARBEDROCK",  Geo::kNoField); // liner
00689   SetMedium(Geo::kHall,       "AIR",          Geo::kNoField); // hall
00690   SetMedium(Geo::kRCCoil,     "ALUM",         Geo::kNoField); // return coil
00691   SetMedium(Geo::kRCAir,      "AIR",          Geo::kNoField); 
00692   SetMedium(Geo::kRCShaft,    "IRON",         Geo::kNoField);
00693 
00694   // The following regions in the main detector but outside the steel
00695   // and coil region will have no field or low field implemented
00696   SetMedium(Geo::kStpScint,   "PSTYRENE SCINT",    Geo::kLowField);
00697   SetMedium(Geo::kStpTiO2,    "COEX TIO2 PSTYRENE",Geo::kLowField);
00698   SetMedium(Geo::kMdlAlSkin,  "ALUM",              Geo::kLowField);
00699   SetMedium(Geo::kMdlAir,     "AIR",               Geo::kLowField);
00700   SetMedium(Geo::kPlnScint,   "AIR",               Geo::kLowField);
00701   SetMedium(Geo::kPlnAirGap,  "AIR",               Geo::kLowField);
00702 
00703   // The following regions in the steel and in the coil region have a
00704   // strong field implemented
00705   SetMedium(Geo::kPlnSteel,   "FEPL",        Geo::kHighField); // steel pln
00706   SetMedium(Geo::kCRGapWater, "WATER",       Geo::kHighField); // air gap CR
00707   SetMedium(Geo::kCRGapCoil,  "ALUM",        Geo::kHighField); 
00708   SetMedium(Geo::kCRGapThroat,"AIR",         Geo::kHighField); 
00709   SetMedium(Geo::kCRGapFlange,"IRON",        Geo::kHighField);
00710   SetMedium(Geo::kCRGapBypass,"AIR",         Geo::kHighField);
00711   SetMedium(Geo::kCRStlWater, "WATER",       Geo::kHighField); // steel pln CR
00712   SetMedium(Geo::kCRStlCoil,  "ALUM",        Geo::kHighField); 
00713   SetMedium(Geo::kCRStlThroat,"AIR",         Geo::kHighField);
00714   SetMedium(Geo::kCRStlNeck,  "IRON",        Geo::kHighField);
00715   SetMedium(Geo::kCRStlHole,  "AIR",         Geo::kHighField);
00716   SetMedium(Geo::kCRStlDetail,"FEPL",        Geo::kHighField);
00717   SetMedium(Geo::kCRSctWater, "WATER",       Geo::kHighField); // scint pln CR
00718   SetMedium(Geo::kCRSctCoil,  "ALUM",        Geo::kHighField); 
00719   SetMedium(Geo::kCRSctThroat,"AIR",         Geo::kHighField);
00720   SetMedium(Geo::kCRSctFlange,"IRON",        Geo::kHighField);
00721   SetMedium(Geo::kCRSctBypass,"AIR",         Geo::kHighField);
00722   
00723 }
00724 
00725 //_____________________________________________________________________________
00726 void GeoMediumMap::BuildCalMediumMap() {
00727   // Private method to build medium map for the calibration detector.
00728 
00729   // No magnetic field is present in the calibration detector
00730   SetMedium(Geo::kMars,       "ROCK",               Geo::kNoField); // rock 
00731   SetMedium(Geo::kLinr,       "CONCRETE",           Geo::kNoField); // liner
00732   SetMedium(Geo::kHall,       "AIR",                Geo::kNoField); // hall
00733   SetMedium(Geo::kVSStpScint, "PSTYRENE SCINT",  Geo::kNoField); // veto shield
00734   SetMedium(Geo::kVSStpTiO2,  "COEX TIO2 PSTYRENE", Geo::kNoField);
00735   SetMedium(Geo::kVSMdlAlSkin,"ALUM",               Geo::kNoField);
00736   SetMedium(Geo::kVSMdlAir,   "AIR",                Geo::kNoField);
00737   SetMedium(Geo::kVSPlnScint, "AIR",                Geo::kNoField); 
00738   SetMedium(Geo::kStpScint,   "PSTYRENE SCINT",     Geo::kNoField);
00739   SetMedium(Geo::kStpTiO2,    "COEX TIO2 PSTYRENE", Geo::kNoField);
00740   SetMedium(Geo::kMdlAlSkin,  "ALUM",               Geo::kNoField);
00741   SetMedium(Geo::kMdlAir,     "AIR",                Geo::kNoField);
00742   SetMedium(Geo::kPlnScint,   "AIR",                Geo::kNoField);
00743   SetMedium(Geo::kPlnAirGap,  "AIR",                Geo::kNoField);
00744   SetMedium(Geo::kPlnSteel,   "FEPL",             Geo::kNoField); // steel pln 
00745   
00746 }
00747 
00748 //_____________________________________________________________________________
00749 void GeoMediumMap::BuildDefMediumMap() {
00750   // Private method to build medium map for a default detector.
00751 
00752   // No magnetic field is present in the default detector which is
00753   // an empty hall
00754   SetMedium(Geo::kMars,       "ROCK",          Geo::kNoField); // rock 
00755   SetMedium(Geo::kLinr,       "CONCRETE",      Geo::kNoField); // liner
00756   SetMedium(Geo::kHall,       "AIR",           Geo::kNoField); // hall
00757   
00758 }
00759 
00760 //_____________________________________________________________________________
00761 void GeoMediumMap::SetMedium(Geo::EDetComponent detcomp, 
00762                              std::string medbasename,
00763                              Geo::EFldStrength fldstrength) {
00764   // Private method to place medium of requested type in fMediumMap
00765   // with key detcomp.  
00766 
00767   Geo::ESwimMethod swimmethod = GetSwimMethod(fldstrength);
00768   
00769   std::string usermedname = UgliLoanPool::Instance()->GetMedium(detcomp);  
00770   if (!usermedname.empty()) {
00771     MSG("Geo",Msg::kInfo) << "Medium for component "
00772                           << Geo::AsString(detcomp) << " set to " 
00773                           << usermedname.c_str() << " overriding default " 
00774                           << medbasename.c_str() << "." << endl;
00775     medbasename = usermedname;
00776   }
00777   
00778   std::string medname = medbasename + "_B" 
00779                       + UtilString::ToString<int>(swimmethod);
00780   GeoMedium* medium 
00781           = dynamic_cast<GeoMedium*>(gGeoManager->GetMedium(medname.c_str()));
00782   if ( !medium ) {
00783     MSG("Geo",Msg::kError) << "GeoMediumMap::SetMedium failed to find "
00784                            << "medium of name " << medname.c_str() 
00785                            << "\nfor component " << Geo::AsString(detcomp)
00786                            << ". Abort." << endl;
00787     abort();
00788   }
00789 
00790   fMediumMap[detcomp] = medium;
00791   fFldStrengthMap[detcomp] = fldstrength;
00792   
00793 }