BeamlineMagneticField_service.cc

Go to the documentation of this file.

////////////////////////////////////////////////////////////////////////////
/// \file    BeamlineMagneticField.h
/// \brief   Utility to provide the magnetic field in the beamline for a
///          given run number/run conditions.
///          NB/ Currently this will use a lookup defined in a fhicl file,
///              it should and will be replaced by an interface to a
///              database when the field is stored in such a way.
/// \author  Mike Wallbank (University of Cincinnati) <wallbank@fnal.gov>
/// \date    August 2019
////////////////////////////////////////////////////////////////////////////

#include "BeamlineMagneticField.h"

// -----------------------------------------------------------------------
beamlineutil::BeamlineMagneticField::BeamlineMagneticField(fhicl::ParameterSet const& pset) {
  fhicl::ParameterSet magnet_currents = pset.get<fhicl::ParameterSet>("MagnetCurrents");
  std::vector<std::pair<unsigned int, float> > current_table
    = magnet_currents.get<std::vector<std::pair<unsigned int, float> > >("table");
  for (std::vector<std::pair<unsigned int, float > >::const_iterator current_table_it = current_table.begin();
       current_table_it != current_table.end(); ++current_table_it)
    fMagnetInfo[current_table_it->first] = current_table_it->second;
  fCurrentToField = pset.get<float>("CurrentToField");
}

// -----------------------------------------------------------------------
beamlineutil::BeamlineMagneticField::BeamlineMagneticField(fhicl::ParameterSet const& pset, art::ActivityRegistry&)
  : BeamlineMagneticField(pset) {
}

// -----------------------------------------------------------------------
float beamlineutil::BeamlineMagneticField::MagnetCurrent(unsigned int run) {

  if (run < fMagnetInfo.begin()->first)
    throw cet::exception("BeamlineMagneticField")
      << "Requested run " << run << " has no magnet current information available.";

  unsigned int run_with_current = run;
  while (!fMagnetInfo.count(run_with_current))
    run_with_current -= 1;

  return fMagnetInfo[run_with_current];

}

// -----------------------------------------------------------------------
float beamlineutil::BeamlineMagneticField::MagneticField(unsigned int run) {

  float magnet_current = this->MagnetCurrent(run);
  float b_field = this->FieldFromCurrent(magnet_current);

  return b_field;

}

// -----------------------------------------------------------------------
float beamlineutil::BeamlineMagneticField::FieldFromCurrent(float magnet_current) {

  // eqn from Mark (technote doesn't exist yet - add docdb info here once it does)
  // eqn gives you the aboslute value of the field - use the sign of the current to determine sign of field
  float field_magnitude = -0.0294 + abs(magnet_current)/985.3 - pow(abs(magnet_current)/3451.2,2);
  // returns value of field with sign of current
  float b_field = copysign(field_magnitude, magnet_current);

  return b_field;

}

DEFINE_ART_SERVICE(beamlineutil::BeamlineMagneticField)