Geometry_service.cc

Go to the documentation of this file.

////////////////////////////////////////////////////////////////////////
/// \file    Geometry.cxx
/// \brief
/// \author  messier@indiana.edu
////////////////////////////////////////////////////////////////////////

#include <cassert>
#include <iostream>
#include <algorithm>

// Framework includes
#include "cetlib/search_path.h"
#include "messagefacility/MessageLogger/MessageLogger.h"
#include "art/Framework/Services/Registry/ServiceHandle.h"
#include "art/Framework/Services/Optional/RandomNumberGenerator.h"
#include "CMap/service/DetectorService.h"
#include "RunHistory/service/RunHistoryService.h"
#include "ChannelInfo/BadChanList.h"
#include "GeometryObjects/OfflineChan.h"
#include "GeometryObjects/CellGeo.h"

// NOvA includes
#include "SummaryData/RunData.h"
#include "Geometry/Geometry.h"
#include "GeometryObjects/Geo.h"

namespace geo
{
  //......................................................................
  Geometry::Geometry(const Parameters& params, art::ActivityRegistry& reg)
    : GeometryBase(params())
    , fLoadedGDMLInfo("(none)")
  {
    // This line is not useless as it might first appear. The call to
    // setDetectorID in the GeometryBase constructor does not call the override
    // in this class and so communicate with DetectorService. Now that the
    // derived class is more fully constructed call it again so it will.
    setDetectorID(fDetId);

    // tell the activity registry to check in after runs have been loaded
    // into the event processing so that we can update configurations if necessary
    reg.sPreBeginRun.watch(this, &Geometry::preBeginRun);
  }

  //......................................................................
  Geometry::~Geometry()
  {
    // Clear up temporary files at shutdown, if necessary
    if(fGDMLFile.find("TEMP_GEO_GDML") != std::string::npos)
      RemoveTmpFile(fGDMLFile);
  }

  //......................................................................
  // 5.15.12 BJR: use the gdml file for both the fGDMLFile and fROOTFile
  // variables as ROOT v5.30.06 is once again able to read in gdml files
  // during batch operation, in this case think of fROOTFile meaning the
  // file used to make the ROOT TGeoManager.  I don't want to remove
  // the separate variables in case ROOT breaks again
  void Geometry::preBeginRun(art::Run const& run)
  {
    // check the current run number, do nothing if it is the same

    if(run.run() == fRunNumber) return;
    fRunNumber = run.run();

     // check here to see if we need to load a new geometry.
      // get the detector id from the run object
    std::vector< art::Handle<sumdata::RunData> > rdcol;
    run.getManyByType(rdcol);
    if(!rdcol.empty() && !fParams.ForceUseFCLOnly()){
      const std::string newGDMLInfo = rdcol[0]->DetGDMLInfo();

      if(newGDMLInfo == fLoadedGDMLInfo){
        std::cout << "Geometry for run " << fRunNumber
                  << " unchanged from previous run." << std::endl;
        // So there's nothing left to do
        return;
      }

      //if there is a gdml string in the file, use that
      if(!newGDMLInfo.empty()){

        //if geometry for last run came from a temp file, remove it
        if (fGDMLFile.find("TEMP_GEO_GDML") != std::string::npos){
          this->RemoveTmpFile(fGDMLFile);
        }

        //gdml info is nonzero, attempt to load geometry from here
        if (this->MakeTmpFile(newGDMLInfo) == 0){
          this->LoadGeometryFile(fGDMLFile, rdcol[0]->DetId());
          fLoadedGDMLInfo = newGDMLInfo;
          std::cout<<"Successfully loaded geometry from in file "<<fGDMLFile<<std::endl;
          std::cout<<"The geometry stored in the file was originally named: "<<rdcol[0]->DetGeometryBaseName()<<std::endl;
          fROOTFile = fGDMLFile; //hook left for root file
          return;
        }
        else {
          std::cout<<"WARNING: Could not load geometry information "
              <<"from within the file. The temporary file was not "
              <<"successfully created, something is wrong.  Will "
              <<"attempt to load from the gdml file directly."
              <<std::endl;
        }
      }

      // get strings with the relative paths for the detector geometry files
      std::string relpathgdml("Geometry/gdml/");
      relpathgdml += rdcol[0]->DetGeometryBaseName();
      relpathgdml += ".gdml";

      if(fGDMLFile.find(relpathgdml) != std::string::npos) return;

      // constructor decides if initialized value is a path or an environment variable
      fGDMLFile = FindGDMLFile(relpathgdml);

      mf::LogWarning("LoadNewGeometry") << "loading new geometry files\n"
                                        << fGDMLFile << "\n";
  
      this->LoadGeometryFile(fGDMLFile);
      fROOTFile = fGDMLFile; //hook left for root file  
      return;
    }
    else{
      //if we made it this far without a load, check to see if the original
      //gdml file specified in the fcl file is still stored in fGDMLFile
      //if this is not the case, make it so and load the geometry
      if (fGDMLFile.find(fGDMLFromFCL) == std::string::npos){
        fGDMLFile = fGDMLFromFCL;
        this->LoadGeometryFile(fGDMLFile);
      }
      mf::LogWarning("LoadNewGeometry") << "cannot find sumdata::RunData object to grab detid\n"
                                        << "or user has determined to only use the value "
                                        << "from the fcl file\n"
                                        << "using currently configured geometry:\n"
                                        << "\t" << fGDMLFile << "\n";
      fROOTFile = fGDMLFile; //hook left for root file
      return;
    }
    
  }// end of preBeginRun
  
  //......................................................................
  bool Geometry::calculateMassesLong(const unsigned int number_of_points) const
  {
    art::ServiceHandle<art::RandomNumberGenerator> rng;
    CLHEP::HepRandomEngine& engine = rng->getEngine();
    return GeometryBase::calculateMassesLong(number_of_points, engine);
  }


  //----------------------------------------------------------------------------
  TVector3 Geometry::getDirectionInBeamCoordinates(TVector3 const& startPos,
                                                   TVector3 const& dir) const
  {

    /// x0 in Detector coordinates. All dimensions are in cm
    const TVector3 x0_det_cm = startPos;

    /// Small displacement in the Prong's direction
    const TVector3 dr        = dir * 0.001;

    /// x1 in Detector coordinates. All dimensions are in cm
    const TVector3 x1_det_cm = x0_det_cm + dr;

    /// x0, x1 in Detector coordinates. All dimensions are CLHEP
    double x0_det[3] = {x0_det_cm.X()*CLHEP::cm
                       ,x0_det_cm.Y()*CLHEP::cm
                       ,x0_det_cm.Z()*CLHEP::cm
                       };
    double x1_det[3] = {x1_det_cm.X()*CLHEP::cm
                       ,x1_det_cm.Y()*CLHEP::cm
                       ,x1_det_cm.Z()*CLHEP::cm
                       };

    /// x0 and x1 in beam coordinates
    double x0_beam[3] = {0.0, 0.0, 0.0};
    double x1_beam[3] = {0.0, 0.0, 0.0};

    /// Calculate beam coordinates for x0 and x1
    const geo::CoordinateTransformation* coordinate_transformation = getCoordinateTransformation();
    coordinate_transformation->transformDetectorToBeamCoordinates(x0_det, x0_beam);
    coordinate_transformation->transformDetectorToBeamCoordinates(x1_det, x1_beam);

    /// Transform into a vector class
    const TVector3 x0_beam_vec(x0_beam);

    /// Direction in beam coordinates. Not normalized
    TVector3 dir_beam(x1_beam);
    dir_beam -= x0_beam_vec;

    /// return Unit vector
    return dir_beam.Unit();
  }

  //......................................................................
  void Geometry::setDetectorID(novadaq::cnv::DetId detid)
  {
    // Intercept all calls to setDetectorID. First set the parameter in
    // GeometryBase as usual.
    GeometryBase::setDetectorID(detid);

    // Make sure to also inform the detector service
    art::ServiceHandle<ds::DetectorService> ds;
    ds->SetDetector(fDetId);
  }

  DEFINE_ART_SERVICE(Geometry)
} // end namespace geo
////////////////////////////////////////////////////////////////////////