UtilityApps/src/print_materials.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 //  Simple program to print all the materials in a detector on
0012 //  a straight line between two given points
0013 //
0014 //  Author     : F.Gaede, DESY
0015 //
0016 //==========================================================================
0017
0018 // Framework include files
0019 #include "DD4hep/Detector.h"
0020 #include "DD4hep/DD4hepUnits.h"
0021 #include "DDRec/MaterialManager.h"
0022
0023 using namespace std ;
0024 using namespace dd4hep::detail;
0025 using namespace dd4hep::rec;
0026 using namespace dd4hep ;
0027
0028
0029 //=============================================================================
0030 #include "main.h"
0031
0032 int main_wrapper(int argc, char** argv ){
0033
0034   if( argc != 8 ) {
0035     std::cout << " usage: print_materials compact.xml x0 y0 z0 x1 y1 z1 " << std::endl
0036               << "        -> prints the materials on a straight line between the two given points ( unit is cm) "
0037               << std::endl ;
0038     exit(1) ;
0039   }
0040
0041   std::string inFile =  argv[1] ;
0042
0043   std::stringstream sstr ;
0044   sstr << argv[2] << " " << argv[3] << " " << argv[4] << " " << argv[5] << " " << argv[6] << " " << argv[7] ;
0045
0046   double x0,y0,z0,x1,y1,z1 ;
0047   sstr >> x0 ;
0048   sstr >> y0 ;
0049   sstr >> z0 ;
0050   sstr >> x1 ;
0051   sstr >> y1 ;
0052   sstr >> z1 ;
0053
0054   Detector& description = Detector::getInstance();
0055
0056   description.fromCompact( inFile );
0057
0058   Vector3D p0( x0, y0, z0 ) ;
0059   Vector3D p1( x1, y1, z1 ) ;
0060
0061   MaterialManager matMgr( description.world().volume() ) ;
0062
0063   const MaterialVec& materials = matMgr.materialsBetween( p0 , p1  ) ;
0064
0065   std::cout  << std::endl  << " #######  materials between the two  points : " << p0 << "*cm  and " << p1 << "*cm :  "  << std::endl ;
0066
0067   double sum_x0 = 0 ;
0068   double sum_lambda = 0 ;
0069   double path_length = 0 ;
0070   for( unsigned i=0,n=materials.size();i<n;++i){
0071
0072     Material mat =  materials[i].first  ;
0073     double length = materials[i].second  ;
0074
0075     double nx0 = length / mat.radLength()  ;
0076     sum_x0 += nx0 ;
0077
0078     double nLambda = length / mat.intLength()  ;
0079     sum_lambda += nLambda ;
0080
0081     path_length += length ;
0082
0083     std::cout <<  "      " << mat << " thickness: " <<  length << " path_length:" << path_length  << " integrated_X0: " <<  sum_x0 << " integrated_lambda: " <<  sum_lambda << std::endl ;
0084   }
0085
0086   std::cout << "############################################################################### "  << std::endl  << std::endl  ;
0087
0088
0089   const MaterialData& avMat = matMgr.createAveragedMaterial( materials ) ;
0090
0091   std::cout << "     averaged Material : " << " Z: " << avMat.Z() << " A: " << avMat.A() << " densitiy: " << avMat.density()
0092             << " radiationLength: " <<  avMat.radiationLength()
0093             << " interactionLength: " << avMat.interactionLength()  << std::endl << std::endl  ;
0094
0095
0096   std::cout << "     Total length : "  << path_length / dd4hep::mm << " mm "  << std::endl  ;
0097
0098   std::cout << "     Integrated radiation lengths : "     << path_length / avMat.radiationLength()  << " X0 " << std::endl  ;
0099
0100   std::cout << "     Integrated interaction lengths : "  << path_length / avMat.interactionLength() << " lambda "   << std::endl  ;
0101
0102   std::cout << "############################################################################### "  << std::endl  << std::endl  ;
0103
0104
0105   return 0;
0106 }
0107
0108 //=============================================================================