UtilityApps/src/dumpBfield.cpp

0001 //==========================================================================
0002 //  AIDA Detector description implementation
0003 //--------------------------------------------------------------------------
0004 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
0005 // All rights reserved.
0006 //
0007 // For the licensing terms see $DD4hepINSTALL/LICENSE.
0008 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
0009 //
0010 //==========================================================================
0011 //
0012 //  Simple program to dump the B-field of the world volume in a cartesian grid
0013 //
0014 //  Author     : F.Gaede, DESY
0015 //  Date       : 29 Mar 2017
0016 //==========================================================================
0017
0018 // Framework include files
0019 #include "DD4hep/Detector.h"
0020 #include "DD4hep/DD4hepUnits.h"
0021
0022 using namespace std ;
0023 using namespace dd4hep ;
0024 using namespace dd4hep::detail;
0025
0026 //=============================================================================
0027
0028
0029 static int invoke_dump_B_field(int argc, char** argv ){
0030
0031   if( argc != 8 ) {
0032     std::cout << " usage: dumpBfield compact.xml xmin[:xmax] ymin[:ymax] zmin[:zmax] dx dy dz [in cm]" << std::endl
0033           << "    will dump the B-field in volume (xmin:xmax, ymin:ymax, zmin:zmax) with steps (dx,dy,dz). All values are in cm."
0034               << "    If a single value is given for a range, symmetric boundaries are used"
0035           << std::endl ;
0036
0037     exit(1) ;
0038   }
0039
0040   std::string inFile =  argv[1] ;
0041
0042   std::stringstream sstr ;
0043   sstr << argv[2] << " " << argv[3] << " " << argv[4] << " " << argv[5] << " " << argv[6] << " " << argv[7] ;
0044
0045   std::string RangeX , RangeY , RangeZ;
0046   float dx , dy, dz ;
0047
0048   sstr >> RangeX >> RangeY >> RangeZ >> dx >> dy >> dz;
0049
0050   size_t colon_posX = RangeX.find(':');
0051   size_t colon_posY = RangeY.find(':');
0052   size_t colon_posZ = RangeZ.find(':');
0053
0054   float minX=0, maxX=0, minY=0, maxY=0, minZ=0, maxZ=0;
0055
0056   if( colon_posX == std::string::npos ) {
0057     std::cout << "X Interval not specified as xmin:xmax" << std::endl
0058       << "  setting xmin = -xmax " << std::endl;
0059     maxX = std::stof( RangeX );
0060     minX = -maxX;
0061   }
0062   else {
0063     minX = std::stof( RangeX.substr(0, colon_posX) );
0064     maxX = std::stof( RangeX.substr(colon_posX+1) );
0065   }
0066
0067   if( colon_posY == std::string::npos ) {
0068     std::cout << "Y Interval not specified as ymin:ymax" << std::endl
0069       << "  setting ymin = -ymax " << std::endl;
0070     maxY = std::stof( RangeY );
0071     minY = -maxY;
0072   }
0073   else {
0074     minY = std::stof( RangeY.substr(0, colon_posY) );
0075     maxY = std::stof( RangeY.substr(colon_posY+1) );
0076   }
0077
0078 if( colon_posZ == std::string::npos ) {
0079     std::cout << "Z Interval not specified as zmin:zmax" << std::endl
0080       << "  setting zmin = -zmax " << std::endl;
0081     maxZ = std::stof( RangeZ );
0082     minZ = -maxZ;
0083   }
0084   else {
0085     minZ = std::stof( RangeZ.substr(0, colon_posZ) );
0086     maxZ = std::stof( RangeZ.substr(colon_posZ+1) );
0087   }
0088
0089   minX *= dd4hep::cm;
0090   maxX *= dd4hep::cm;
0091   minY *= dd4hep::cm;
0092   maxY *= dd4hep::cm;
0093   minZ *= dd4hep::cm;
0094   maxZ *= dd4hep::cm;
0095   dx *= dd4hep::cm;
0096   dy *= dd4hep::cm;
0097   dz *= dd4hep::cm;
0098
0099   Detector& description = Detector::getInstance();
0100   description.fromCompact( inFile );
0101
0102   printf("#######################################################################################################\n");
0103   printf("      x[cm]            y[cm]            z[cm]          Bx[Tesla]        By[Tesla]        Bz[Tesla]     \n");
0104
0105   for( float x = minX ; x <= maxX ; x += dx ){
0106     for( float y = minY ; y <= maxY ; y += dy ){
0107       for( float z = minZ ; z <= maxZ ; z += dz ){
0108
0109     double posV[3] = { x, y, z }  ;
0110     double bfieldV[3] ;
0111     description.field().magneticField( posV  , bfieldV  ) ;
0112
0113     printf(" %+15.8e  %+15.8e  %+15.8e  %+15.8e  %+15.8e  %+15.8e  \n",
0114          posV[0]/dd4hep::cm, posV[1]/dd4hep::cm,  posV[2]/dd4hep::cm,
0115          bfieldV[0]/dd4hep::tesla , bfieldV[1]/dd4hep::tesla, bfieldV[2]/dd4hep::tesla ) ;
0116
0117       }
0118     }
0119   }
0120   printf("#######################################################################################################\n");
0121   return 0;
0122 }
0123
0124
0125
0126 int main(int argc, char** argv ){
0127   try {
0128     return invoke_dump_B_field(argc,argv);
0129   }
0130   catch(const std::exception& e)  {
0131     std::cout << "Got uncaught exception: " << e.what() << std::endl;
0132   }
0133   catch (...)  {
0134     std::cout << "Got UNKNOWN uncaught exception." << std::endl;
0135   }
0136   return EINVAL;
0137 }
0138
0139 //=============================================================================