MagneticField/src/FieldMapReadBack.cc

0001
0002 #include "FieldMapReadBack.h"
0003
0004 #include <phfield/PHField.h>
0005 #include <phfield/PHField3DCartesian.h>
0006 #include <phfield/PHFieldUtility.h>
0007
0008 #include <fun4all/Fun4AllReturnCodes.h>
0009
0010 #include <phool/PHCompositeNode.h>
0011
0012 #include <Geant4/G4SystemOfUnits.hh>
0013
0014 #include <iostream>
0015
0016 //____________________________________________________________________________..
0017 int FieldMapReadBack::InitRun(PHCompositeNode *topNode)
0018 {
0019   fieldmap = PHFieldUtility::GetFieldMapNode(nullptr, topNode);
0020   if (fieldmap)
0021   {
0022     std::cout << "Found or created fieldmap" << std::endl;
0023   }
0024   else
0025   {
0026     std::cout << "Fieldmap not found or created" << std::endl;
0027   }
0028   return Fun4AllReturnCodes::EVENT_OK;
0029 }
0030
0031 //____________________________________________________________________________..
0032 int FieldMapReadBack::process_event(PHCompositeNode *topNode)
0033 {
0034   double bfield[3];
0035   // the output needs to be converted into the unit you want - here it is tesla
0036   fieldmap->GetFieldValue(Point, bfield);
0037   std::cout << "bx: " << bfield[0] / tesla
0038             << " by: " << bfield[1] / tesla
0039             << " bz: " << bfield[2] / tesla
0040             << std::endl;
0041   return Fun4AllReturnCodes::EVENT_OK;
0042 }
0043
0044 void FieldMapReadBack::Load3dCartMap(const std::string &fname, const float magfield_rescale)
0045 {
0046   fieldmap = new PHField3DCartesian(fname, magfield_rescale);
0047 }
0048
0049 void FieldMapReadBack::PrintField(const double x, const double y, const double z, const double t)
0050 {
0051   SetFieldPoint(x, y, z, t);
0052
0053   double Bf[3];
0054   fieldmap->GetFieldValue(Point, Bf);
0055   std::cout << "Point: " << x << "/" << y << "/" << z << " cm" << std::endl;
0056   std::cout << "BField: " << Bf[0] / tesla << "/" << Bf[1] / tesla << "/" << Bf[2] / tesla << std::endl;
0057   return;
0058 }
0059
0060 void FieldMapReadBack::SetFieldPoint(const double x, const double y, const double z, const double t)
0061 {
0062   Point[0] = x * cm;
0063   Point[1] = y * cm;
0064   Point[2] = z * cm;
0065   Point[3] = t * s;
0066 }
0067
0068 void FieldMapReadBack::Verbosity(const int i)
0069 {
0070   if (fieldmap)
0071   {
0072     fieldmap->Verbosity(i);
0073   }
0074   Fun4AllBase::Verbosity(i);
0075   return;
0076 }