DataUtil/EnergyCorrections.cxx

Go to the documentation of this file.

//
// EnergyCorrections.h
//
// Corrections to energy based quantities
//
// Created:  M. Kordosky -- December, 2005
//
// $Author: xbhuang $
//
// $Revision: 1.25 $
//
// $Name:  $
//
// $Id: EnergyCorrections.cxx,v 1.25 2009/12/02 15:07:14 xbhuang Exp $
//
// $Log: EnergyCorrections.cxx,v $
// Revision 1.25  2009/12/02 15:07:14  xbhuang
// Include the fudge factor (1.0000) for the Dogwood1 calibration.
//
// Revision 1.24  2009/08/05 16:09:40  xbhuang
//
// Modified in order to add corrections for the Dogwood0 reconstruction.
//
// Revision 1.23  2009/07/22 17:48:33  xbhuang
//
// The far detector shower energy tuning function for Dogwood has been updated.
//
// Revision 1.22  2009/07/04 01:53:52  rhatcher
// Remove tabs (replaced w/ spaces per coding conventions).  Did a bit of
// other coding convention cleanup, but there is still a lot of ugliness
// about the formatting that makes the code very hard to read.
//
// Revision 1.21  2009/07/03 14:12:39  xbhuang
//
// Defined 4 functions to do CC shower energy correction for Dogwood reconstruction.
//
// Revision 1.20  2007/12/21 15:36:08  rjn
// Added switch for the ability to go back to the May version of the Masaki conversion
//
// Revision 1.19  2007/12/21 11:12:32  rjn
// Added a WeightedShowerEnergyConversionDogwood function
//
// Revision 1.18  2007/12/19 16:39:48  rjn
// Added Masaki's correction from DocDB 3895_v4
//
// Revision 1.17  2007/11/16 11:43:21  rjn
// Added ShowerEnergyConversionDogwood which is the function to be called by CandShowerHandle to apply ourbest knowledge of the conversion from MEU to GeV for the dogwood reconstruction, such that to first order the GeV number in the tree is correct (until we know more at least).
//
// Revision 1.16  2007/11/11 06:41:26  rhatcher
// Part of the ongoing purge of "DetectorType" in favor of "Detector"
// (the former is a synonym of the latter ...)
//
// Revision 1.15  2007/10/09 14:25:52  hartnell
//
// One of the useful new messages was missing a "new line" at the end of it.
//
// Revision 1.14  2007/09/28 17:37:10  rjn
// Adding the latest version of EnergyCorrections to the repository. I'm sure this is going to upset somebody, but then again it always does.
//
// Revision 1.13  2007/06/22 12:38:12  rjn
// Hack for Release Type problem. The auto generated release types do not pick up the reconstruction sub versions at present
//
// Revision 1.12  2007/06/04 18:28:23  rjn
// The preliminary version of the Calibration Group final MEU scale fudge factors are now available.
//
// Things to note for those users who are only interested in calibrating shower energy scale for the CC analysis, you can use the FullyCorrectShowerEnergy function, to apply both Calibration Group corrections and Masaki's shower energy tuning. Note that it is very important that you get the release type correct, for the analysis you should be using:
// Data -- Cedar_phy
// MC -- Cedar R1.24.1
// These are the only data sets for which calibration numbers are provided.
//
//
// If you are interested in making data or MC plots in fully calibrated MEU numbers you can use the CalibrationGroupEnergyCorrections function. This will convert MEU's from the sntp into fully calibrated MEU's. NB: It is very inmportant to be careful not to apply the corrections twice by calling both CalibrationGroupEnergyCorrections and FullyCorrectShowerEnergy. If you wish to do this investigate the WhichCorrection_t options and particular the kNoCalGroup option.
//
// Revision 1.11  2007/05/21 09:55:48  rjn
// Just a tidying up sort of commit. Moved the Calibration Group corrections from the obliquely named CorrectionsForMasaki to the rather easier to understand EnergyCorrections::CalibrationGroupEnergyCorrections. NB: This function is called automatically by FullyCorrectShowerEnergy, so do not pass already corrected values to FullyCorrectShowerEnergy.
//
// Revision 1.10  2007/05/19 12:18:39  rjn
// Added to the FullyCorrect family of EnergyCorrections. Now have corrections for shower energy and track momentum/energy from range. At this point users are encouraged to use these. There will still be a few under the hood changes but at this point I think there is a workable interface.
//
// Revision 1.9  2007/05/18 16:24:54  rjn
// Added FullyCorrectShowerEnergy function that should replaces the previously used shower energy corrections. Should work with birch and cedar. It also applies the appropriate calibration group corrections. Everyone should use this who used any of the previous functions.
//
// Revision 1.8  2007/03/16 19:48:41  kordosky
// Energy corrections for the FD!
//
// Revision 1.7  2007/03/13 16:36:15  kordosky
// Energy Corrections for Cedar reconstruction
// ND ONLY!!!
// See 2785-v1 for the correction and 2767 for a description of the work (done by Masaki).
//
// Revision 1.6  2007/02/01 20:06:50  ishi
// Semicolons at the end of namespace{} were removed.
//
// Revision 1.5  2007/02/01 17:47:33  kordosky
// New version of EnergyCorrections. EnergyCorrections is now a namespace and calling code should invoke "using namespace EnergyCorrections". This is not expected
// to cause problems since in the previous version all the correction functions were at global scope. Indeed, code including EnergyCorrections.h have been modified. Will handle cedar and birch by using static variables fVersion and fSubVersion, modified via SetCorrectionVersion(). Users need to do this, perhaps in their job macro. So far there is only one correction of each type and fSubVersion should be held at the default 0. Also, note the utility VersionFromFilename().
//
// Revision 1.4  2006/02/15 15:57:29  kordosky
// add correction from for momentum via curvature as advocated in 1430-v2. Only modifies momenta of positively charged tracks.
//
// Revision 1.3  2006/02/13 00:20:02  kordosky
// correction of +1.8% to FD data shower energy
//
// Revision 1.2  2006/02/12 22:09:58  kordosky
// Update momentum correction to account for new density and thickness measurements.
//
// Revision 1.1  2005/12/16 03:38:51  kordosky
// place implementation in a cxx file
//
// Revision 1.4  2005/12/16 02:52:43  kordosky
// Forgot the include guards.  Thansk Niki.
//
// Revision 1.3  2005/12/15 22:44:31  kordosky
// retweaked A Culling formulae in reaction to his presentation this morning
//
// Revision 1.2  2005/12/15 22:26:19  kordosky
// E=sqrt(p*p+m*m)  != p   (almost negligible)
//
// Revision 1.1  2005/12/15 21:57:09  kordosky
// moved from Mad
//
// Revision 1.2  2005/12/14 20:35:59  kordosky
// Updated energy corrections.
//
// Revision 1.1  2005/12/13 19:41:35  kordosky
// Standalone routines to correct the momentum via range to to groom et al range tables. Also, correct data to 7.755 g/cc which is our current best estimate of the detector density.
//

#include "DataUtil/EnergyCorrections.h"
//#include "DataUtil/CalibratedShowerEnergy.h"
#include "Calibrator/Calibrator.h"
#include "MessageService/MsgService.h"
#include <cmath>
#include <iostream>
#include <sstream>

CVSID( " $Id: EnergyCorrections.cxx,v 1.25 2009/12/02 15:07:14 xbhuang Exp $");



// Shower Energy Correction
// ==========================
//
// All users are recommended to use the FullyCorrectShowerEnergy function to
// apply corrections to shower energy. This function can be used for all
// shower types and performs Calibration Group + other corrections. NB:
// one now needs access to the offline database to call this function.
//

// New era of time dependent energy corrections
float EnergyCorrections::FullyCorrectShowerEnergy(float E,// Input calorimetric shower energy in GeV
                                                  const CandShowerHandle::ShowerType_t& st, //Weighted not supported
                                                  VldContext vc, //Context of this specific event
                                                  ReleaseType::Release_t release, // Release version of the batch processing
                                                  EnergyCorrections::WhichCorrection_t whichCor)
{

  //First step is to apply Calibration Group Corrections
  float eCor=EnergyCorrections::CalibrationGroupEnergyCorrections(E,vc,release,whichCor);

  // Get the reconstruction version.
  ReleaseType::Release_t recoVers = ReleaseType::GetRecoInfo(release);

  //Now need to apply Masaki's correction
  if (ReleaseType::IsBirch(release)) {
    int mode=1;
    if (whichCor==EnergyCorrections::kVersion2)
      mode=2;
    if (vc.GetDetector()==Detector::kNear)
      return EnergyCorrections::CorrectShowerEnergyNear_Birch(eCor,st,mode);
     else if (vc.GetDetector()==Detector::kFar)
       return EnergyCorrections::CorrectShowerEnergyFar_Birch(eCor,st,mode);

  }
  else if (ReleaseType::IsCedar(release)) {
    //Now need to check version
    if (recoVers >= ReleaseType::kR1_24_3) {
      if (vc.GetDetector()==Detector::kNear)
    return EnergyCorrections::ShowerEnergyCorrectionNearCedarPhyBhcurve(eCor,st,whichCor);
      else if (vc.GetDetector()==Detector::kFar)
    return EnergyCorrections::ShowerEnergyCorrectionFarCedarPhyBhcurve(eCor,st,whichCor);

    }
    else {
      if (vc.GetDetector()==Detector::kNear)
    return EnergyCorrections::ShowerEnergyCorrectionNearCedar(eCor,st,whichCor);
      else if (vc.GetDetector()==Detector::kFar)
    return EnergyCorrections::ShowerEnergyCorrectionFarCedar(eCor,st,whichCor);
    }

  }
  else if (ReleaseType::IsDogwood(release)) {
    if (vc.GetDetector()==Detector::kNear)
      return EnergyCorrections::ShowerEnergyCorrectionNearDogwood(eCor,st,whichCor);
    else if (vc.GetDetector()==Detector::kFar)
      return EnergyCorrections::ShowerEnergyCorrectionFarDogwood(eCor,st,whichCor);
  }

  return E;
}



//
// All users are recommended to use the FullyCorrect family of track functions
// apply corrections to track momentum/energy.
//
//The following routines should be used to correct track energy
//and momentum.
float EnergyCorrections::FullyCorrectMomentumFromRange(float p, //Momentum from range
                                                       VldContext vc, //Context of this specific event
                                                       ReleaseType::Release_t release, //Release verion of the batch processing
                                                       EnergyCorrections::WhichCorrection_t whichCor)
{
   float pcor=p;
   if (ReleaseType::IsBirch(release)) {
      pcor=EnergyCorrections::MomentumRangeCorrectionBirch(p,vc,whichCor);
   }
   else if (ReleaseType::IsCedar(release)) {
      pcor=EnergyCorrections::MomentumRangeCorrectionCedar(p,vc,whichCor);
   }
   return pcor;

}

float EnergyCorrections::FullyCorrectEnergyFromRange(float E, //Energy from range
                                                     VldContext vc, //Context of this specific event
                                                     ReleaseType::Release_t release, //Release verion of the batch processing
                                                     EnergyCorrections::WhichCorrection_t whichCor)
{
  if (ReleaseType::IsBirch(release)) {
    const float m=0.1057;// muon mass
    float p = sqrt(E*E -m*m);
    float pcor = EnergyCorrections::FullyCorrectMomentumFromRange(p,vc,release,whichCor);
    return sqrt(pcor*pcor +m*m);
  }
  else if (ReleaseType::IsCedar(release)) {
    return EnergyCorrections::EnergyRangeCorrectionCedar(E,vc,whichCor);
  }
  return E;
}

float EnergyCorrections::FullyCorrectSignedMomentumFromCurvature(float p, //Momentum from curvature
                                                                 VldContext vc, //Context of this specific event
                                                                 ReleaseType::Release_t release, //Release verion of the batch processing
                                                                 EnergyCorrections::WhichCorrection_t whichCor)
{
   float pcor=p;
   if (ReleaseType::IsBirch(release)) {
      pcor=EnergyCorrections::SignedMomentumCurvatureCorrectionBirch(p,vc,whichCor);
   }
   else if (ReleaseType::IsCedar(release)) {
      pcor=EnergyCorrections::SignedMomentumCurvatureCorrectionCedar(p,vc,whichCor);
   }
   return pcor;

}






//The Calibration Group Fudge Factors
namespace EnergyCorrections {

  static const float cgffBirchDataFD=1.018;

  //These Cedar_Phy numbers and CedarDaikon numbers are taken from
  //DocDb 3322_v1
  static const float cgffCedarPhyDataFD=0.9988;
  static const float cgffCedarPhyDataND=1.0039;
  static const float cgffCedarR1_24_1MCFD=(0.001737/0.001786);
  static const float cgffCedarR1_24_1MCND=(0.001759/0.001792);

  //These are the old Cedar_Phy numbers and CedarDaikon numbers that were
  //used in the Summer 2007 CC analysis for the summer conferences
  static const float cgffCedarPhyDataFDOld=0.996;
  static const float cgffCedarPhyDataNDOld=1.003;
  static const float cgffCedarR1_24_1MCFDOld=0.986;
  static const float cgffCedarR1_24_1MCNDOld=0.994;

  //These are currently just place holders for the fudge factors that may
  // be necessary for CedarPhyDaikon -- with the new GevPerMip
  static const float cgffCedarR1_24_2MCFD=cgffCedarR1_24_1MCFD*(0.001786/0.001737);
  static const float cgffCedarR1_24_2MCND=cgffCedarR1_24_1MCND*(0.001792/0.001759);

  static const float cgffCedarPhyDaikonND=1;
  static const float cgffCedarPhyDaikonFD=1;

  //Dogwood1 calibration
  static const float cgffDogwood1DataFD=1.0;
  static const float cgffDogwood1DataND=1.0;
  static const float cgffDogwood1DaikonFD=1.0;
  static const float cgffDogwood1DaikonND=1.0;
}

float EnergyCorrections::CalibrationGroupEnergyCorrections(float E, //Energy in MEU, GeV or SigMap
                                                           VldContext vc, //Context of this specific event
                                                           ReleaseType::Release_t release, //Release verion of the batch processing
                                                           EnergyCorrections::WhichCorrection_t whichCor)
{
  //Nb: This function is called by FullyCorrectShowerEnergy, do not use it for CC Shower's

  float retval = E; // Default return value.

  MAXMSG("DataUtil",Msg::kInfo,1)
    << "CalibrationGroupEnergyCorrections:: Using Release Type: "
    << ReleaseType::AsString(release) << "\tusing correction version: "
    << whichCor << "\n";

  static Calibrator* customCalibrator = 0;
  if (customCalibrator ==0 ) {
    //Initialize. Disable warning that we're about to do something 'tricksy'.
    // This tricksyness is just to make a custom instance of the calibrator, so we don't
    // interfere with any other custom settings that the user (i.e. Jeff H) is doing
    // This is actually pretty safe, if an expert like me (Nathaniel) does it.
    MsgService::Instance()->GetStream("Calib")->SetLogLevel(Msg::kFatal);
    customCalibrator = Calibrator::CreateCustomCalibrator();

    // Turn off the units we aren't using, just to make them go faster.
    // "Speed holes!" -Homer Simpson

    customCalibrator->Set("Thermometer=SimpleCalScheme "
                          "PeCalibrator=SimpleCalScheme "
                          "LinCalibrator=SimpleCalScheme "
                          "StripCalibrator=SimpleCalScheme "
                          "AttenCalibrator=SimpleCalScheme "
                          "MIPCalibrator=SimpleCalScheme "
                          "TimeCalibrator=SimpleCalScheme ");
    MsgService::Instance()->GetStream("Calib")->SetLogLevel(Msg::kWarning);
  }



  // Get the reconstruction version.
  ReleaseType::Release_t recoVers = ReleaseType::GetRecoInfo(release);

  // Fudge #0:
  // For PRL-era data, we fudged the FD data by 1.8%
  if (    ReleaseType::IsBirch(release)
     && vc.GetSimFlag()==SimFlag::kData
     && vc.GetDetector()==Detector::kFar )
    {
      MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Birch Far Detector Factor ( "
    << EnergyCorrections::cgffBirchDataFD << ")\n";
      retval = retval*EnergyCorrections::cgffBirchDataFD;
    }



  if (ReleaseType::IsCedar(release) && whichCor!=EnergyCorrections::kNoCalGroup) {
    //
    // Fudge #1
    // Attempt to remove the MC bug, where the drift decalibration was applied with the wrong sign.
    // This is a time-dependent correction which should fix the problem for cedar-daikon MC R1.24.1
    // This problem was fixed in R1.24.2 and R1.24.calB, but NOT R1.24.calA.
    //
    if (   ReleaseType::IsCedar(release)
      && ReleaseType::IsDaikon(release)
      && (recoVers < ReleaseType::kR1_24_2 || recoVers == ReleaseType::kR1_24_Cal )
      && vc.GetSimFlag()==SimFlag::kMC)
      {
    // Removes MC bug where drift correction was applied twice.
    // Fix it by de-applying the drift twice.
    MAXMSG("DataUtil",Msg::kInfo,1)
      << "EnergyCorrections -- Applying Infamous MC Drift Bug Correction\n";

    customCalibrator->Reset(vc);
    float drift = customCalibrator->GetDriftCorrected(1.0,PlexStripEndId());
    retval = retval / (drift*drift);
      }


    // Correction #2
    // These are the preliminary final corrections for the June 2007 cc box opening.
    // The numbers are taken from DocDB 3139 by Jeff Hartnell and Tingjun Yang
    // These are only applied if using the default (or version3) correction
    if (!(whichCor==EnergyCorrections::kVersion1 || whichCor==EnergyCorrections::kVersion2)) {

      if (vc.GetSimFlag()==SimFlag::kMC) {
    //For the Cedar MC we have the following possible datasets:
    // i) CedarDaikon -- used in Summer analysis 2007
    // ii) CedarPhyDaikon -- reprocessed by batch group over summer 2007

    if (ReleaseType::IsCedar(release) && recoVers<=ReleaseType::kR1_24_1) {
      //Now here we have two options either use the latest greatest numbers
      // or we use the old numbers from the 2007 summer analysis
      if (whichCor==kVersion3) {
        //Use the old numbers
        if (vc.GetDetector()==Detector::kFar) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying R_1_24_1 (Summer 2007) Far MC Correction Factor"
        << " (" << EnergyCorrections::cgffCedarR1_24_1MCFDOld << ")\n";
          retval*=EnergyCorrections::cgffCedarR1_24_1MCFDOld;
        }
        else if (vc.GetDetector()==Detector::kNear) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying R_1_24_1 (Summer 2007) Near MC Correction Factor"
        << " (" << EnergyCorrections::cgffCedarR1_24_1MCNDOld << ")\n";
          retval*=EnergyCorrections::cgffCedarR1_24_1MCNDOld;
        }
      }
      else {
        //Use the new numbers
        if (vc.GetDetector()==Detector::kFar) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying R_1_24_1 (Fall 2007) Far MC Correction Factor"
        << " (" << EnergyCorrections::cgffCedarR1_24_1MCFD << ")\n";
          retval*=EnergyCorrections::cgffCedarR1_24_1MCFD;
        }
        else if (vc.GetDetector()==Detector::kNear) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying R_1_24_1 (Fall 2007) Near MC Correction Factor"
        << " (" << EnergyCorrections::cgffCedarR1_24_1MCND << ")\n";
          retval*=EnergyCorrections::cgffCedarR1_24_1MCND;
        }

      }
    }
    else if (ReleaseType::IsCedar(release) && recoVers>=ReleaseType::kR1_24_2)
      {
        //Changed because we updated the GevPerMIP number in the MC
        // reprocessing. Bloody annoying.
        if (vc.GetDetector()==Detector::kFar) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying R_1_24_2 (Fall 2007) Far MC Correction Factor"
        << " (" << EnergyCorrections::cgffCedarR1_24_2MCFD << ")\n";
          retval*=EnergyCorrections::cgffCedarR1_24_2MCFD;
        }
        else if (vc.GetDetector()==Detector::kNear) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying R_1_24_2 (Fall 2007) Near MC Correction Factor"
        << " (" << EnergyCorrections::cgffCedarR1_24_2MCND << ")\n";
          retval*=EnergyCorrections::cgffCedarR1_24_2MCND;
        }
      }
      }
      else if (vc.GetSimFlag()==SimFlag::kData) {
    //Data corrections only available for Cedar_Phy
    //Once again we are going to have the option to use the old or
    // the new numbers
    if (whichCor==kVersion3) {
      //Use the old numbers from the Summer 2007 analysis
      if (recoVers>=ReleaseType::kCedarPhy) {
        if (vc.GetDetector()==Detector::kFar) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying CedarPhy (Summer 2007) Far Correction Factor"
        << " (" << EnergyCorrections::cgffCedarPhyDataFDOld << ")\n";
          retval*=EnergyCorrections::cgffCedarPhyDataFDOld;
        }
        else if (vc.GetDetector()==Detector::kNear) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying CedarPhy (Summer 2007) Near Correction Factor"
        << " (" << EnergyCorrections::cgffCedarPhyDataNDOld << ")\n";
          retval*=EnergyCorrections::cgffCedarPhyDataNDOld;
        }
      }
    }
    else {
      //Use the latest greatest numbers
      if (recoVers>=ReleaseType::kCedarPhy) {
        if (vc.GetDetector()==Detector::kFar) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying CedarPhy (Fall 2007) Far Correction Factor"
        << " (" << EnergyCorrections::cgffCedarPhyDataFD << ")\n";
          retval*=EnergyCorrections::cgffCedarPhyDataFD;
        }
        else if (vc.GetDetector()==Detector::kNear) {
          MAXMSG("DataUtil",Msg::kInfo,1)
        << "EnergyCorrections -- Applying CedarPhy (Fall 2007) Near Correction Factor"
        << " (" << EnergyCorrections::cgffCedarPhyDataND << ")\n";
          retval*=EnergyCorrections::cgffCedarPhyDataND;
        }
      }
    }
      }
    }
  }

  return retval;

}



// Dogwood Shower Energy Conversion
// =================================
//
//
//New function to be called directly from CandShowerHandle to set the GeV
// value in the reco tree to be as close as possible to correct for the
// forthcoming Dogwood release
float EnergyCorrections::ShowerEnergyConversionDogwood(float E, VldContext vc)
{
  //Two things to do:
  //i) apply Calibration Group tweaks
  //ii) apply Maskai's Reco->Truth formula
  float eCor=E;
  if (vc.GetDetector()==Detector::kFar) {
    //Cal group corrections
    if (vc.GetSimFlag()==SimFlag::kData) {
      eCor*=EnergyCorrections::cgffDogwood1DataFD;
    }
    else if (vc.GetSimFlag()==SimFlag::kMC) {
      eCor*=EnergyCorrections::cgffDogwood1DaikonFD;
    }
    //Reco-Truth Conversion
    eCor=EnergyCorrections::ShowerEnergyCorrectionFarDogwood(eCor,CandShowerHandle::kCC,EnergyCorrections::kDefault);
  }
  else if (vc.GetDetector()==Detector::kNear) {
    //Cal group corrections
    if (vc.GetSimFlag()==SimFlag::kData) {
      eCor*=EnergyCorrections::cgffDogwood1DataND;
    }
    else if (vc.GetSimFlag()==SimFlag::kMC) {
      eCor*=EnergyCorrections::cgffDogwood1DaikonND;
    }
    //Reco-Truth Conversion
    eCor=EnergyCorrections::ShowerEnergyCorrectionNearDogwood(eCor,CandShowerHandle::kCC,EnergyCorrections::kDefault);
  }

  return eCor;

}



//New function to be called directly from CandShowerHandle to set the GeV
// value in the reco tree to be as close as possible to correct for the
// forthcoming Dogwood release
float EnergyCorrections::WeightedShowerEnergyConversionDogwood(float E, VldContext vc)
{
  //Two things to do:
  //i) apply Calibration Group tweaks
  //ii) apply Maskai's Reco->Truth formula
  float eCor=E;
  if (vc.GetDetector()==Detector::kFar) {
    //Cal group corrections
    if (vc.GetSimFlag()==SimFlag::kData) {
      eCor*=EnergyCorrections::cgffDogwood1DataFD;
    }
    else if (vc.GetSimFlag()==SimFlag::kMC) {
      eCor*=EnergyCorrections::cgffDogwood1DaikonFD;
    }
    //Reco-Truth Conversion
    eCor=EnergyCorrections::ShowerEnergyCorrectionFarDogwood(eCor,CandShowerHandle::kWtCC,EnergyCorrections::kDefault);
  }
  else if (vc.GetDetector()==Detector::kNear) {
    //Cal group corrections
    if (vc.GetSimFlag()==SimFlag::kData) {
      eCor*=EnergyCorrections::cgffDogwood1DataND;
    }
    else if (vc.GetSimFlag()==SimFlag::kMC) {
      eCor*=EnergyCorrections::cgffDogwood1DaikonND;
    }
    //Reco-Truth Conversion
    eCor=EnergyCorrections::ShowerEnergyCorrectionNearDogwood(eCor,CandShowerHandle::kWtCC,EnergyCorrections::kDefault);
  }

  return eCor;

}


//**************** Legacy *************************************************************************************
//*************************************************************************************************************
//EnergyCorrections::CorrectionVersion_t EnergyCorrections::fVersion = EnergyCorrections::kUnknown;
//Short_t EnergyCorrections::fSubVersion = 0;

namespace EnergyCorrections {
  static CorrectionVersion_t fVersion = EnergyCorrections::kUnknown;
  static Short_t fSubVersion = 0;
}

void EnergyCorrections::SetCorrectionVersion(const CorrectionVersion_t& ver, 
                                             Short_t subver)
{
  fVersion=ver;
  fSubVersion=subver;
}

std::string EnergyCorrections::GetCorrectionAsString() 
{
  std::string s;
  switch(fVersion){
  case kBirch:
    s+="BIRCH"; break;
  case kCedar:
    s+="CEDAR"; break;
  case kUnknown:
  default:
    s+="???";
    break;
  }
  std::ostringstream os; os<<"-v"<<fSubVersion<<std::ends;
  s+=os.str();

  return s;
}

void EnergyCorrections::WarnUnknownVersion(const char* caller_routine)
{
  static Short_t nwarn=0;
  if (nwarn<=9) {
    std::cerr<<"Energy Corrections: In "<<caller_routine
         <<"Energy Corrections: Warning, unknown correction version "
         <<GetCorrectionAsString()
         <<"Energy Corrections: Defaulting to Birch era corrections.\n"
         <<"Energy Corrections: Please Call SetCorrectionVersion() in the future.\n"
         <<std::endl;
    nwarn++;
  }
  if (nwarn==9) {
    std::cerr<<"Energy Corrections: last message of this type..."<<std::endl;
  }

}

EnergyCorrections::CorrectionVersion_t EnergyCorrections::VersionFromFilename(const char* name)
{
  CorrectionVersion_t ver = kUnknown;
  std::string s=name;

  if (s.find("R1_18")!=std::string::npos) {
    ver=EnergyCorrections::kBirch;
  }
  else if (s.find("cedar")!=std::string::npos) {
    ver=EnergyCorrections::kCedar;
  }
  return ver;
}



float EnergyCorrections::CorrectMomentumFromRange(float p, bool isdata,
                                                  Detector::Detector_t det)
{
  float pcor=p;
    switch(fVersion){
    case kCedar:
      pcor=CorrectMomentumFromRange_Cedar(p,isdata,det);
      break;
    case kBirch:
      pcor=CorrectMomentumFromRange_Birch(p,isdata,det);
      break;
    case kUnknown:
    default:
      WarnUnknownVersion("CorrectMomentumFromRange()");
      pcor=CorrectMomentumFromRange_Birch(p,isdata,det);
      break;
    }
    return pcor;
}

float EnergyCorrections::CorrectSignedMomentumFromCurvature(float p, 
                                                            bool isdata , 
                                                            Detector::Detector_t det )
{

  float pcor=p;
    switch(fVersion){
    case kCedar:
      pcor=CorrectSignedMomentumFromCurvature_Cedar(p,isdata,det);
      break;
    case kBirch:
      pcor=CorrectSignedMomentumFromCurvature_Birch(p,isdata,det);
      break;
    case kUnknown:
    default:
      WarnUnknownVersion("CorrectSignedMomentumFromCurvature()");
      pcor=CorrectSignedMomentumFromCurvature_Birch(p,isdata,det);
      break;
    }
    return pcor;
}


float EnergyCorrections::CorrectEnergyFromRange(float E, bool isdata, Detector::Detector_t det)
{
   const float m=0.1057;// mon mass
   float p = sqrt(E*E -m*m);
   float pcor = CorrectMomentumFromRange(p,isdata,det);
   return sqrt(pcor*pcor +m*m);
}


float EnergyCorrections::CorrectShowerEnergyNear(float E, const CandShowerHandle::ShowerType_t& st, int mode, bool isdata )
{
   //People should not be using this function

  float ecor=E;
  switch(fVersion){
  case kCedar:
    ecor=CorrectShowerEnergyNear_Cedar(E,st,mode,isdata);
    break;
  case kBirch:
    ecor=CorrectShowerEnergyNear_Birch(E,st,mode,isdata);
    break;
  case kUnknown:
  default:
    WarnUnknownVersion("CorrectShowerEnergyNear()");
    ecor=CorrectShowerEnergyNear_Birch(E,st,mode,isdata);
    break;
  }
  return ecor;
}



float EnergyCorrections::CorrectShowerEnergyFar(float E, const CandShowerHandle::ShowerType_t& st, int mode, bool isdata)
{
   //People should not be using this function

  float ecor=E;
  switch(fVersion){
  case kCedar:
    ecor=CorrectShowerEnergyFar_Cedar(E,st,mode,isdata);
    break;
  case kBirch:
     if (isdata) {
    // a correction for the FD MIP scale
    // for the beam data, one measured MIP
    // actually corresponds to 1.018 MIPs
    // so we must correct shower energy up
    const float mip_scale_correction=1.018;
    E*=mip_scale_correction;
     }
    ecor=CorrectShowerEnergyFar_Birch(E,st,mode,isdata);
    break;
  case kUnknown:
  default:
    WarnUnknownVersion("CorrectShowerEnergyFar()");
    if (isdata) {
       // a correction for the FD MIP scale
       // for the beam data, one measured MIP
       // actually corresponds to 1.018 MIPs
       // so we must correct shower energy up
       const float mip_scale_correction=1.018;
    E*=mip_scale_correction;
    }
    ecor=CorrectShowerEnergyFar_Birch(E,st,mode,isdata);
    break;
  }
  return ecor;

}

float EnergyCorrections::CorrectShowerEnergy(float E, 
                                             const Detector::Detector_t& det,
                                             const CandShowerHandle::ShowerType_t& st,
                                             int mode, bool isdata)
{
   //People should not be using this function
  float ecor=E;
   if (det==Detector::kNear) {
      ecor = CorrectShowerEnergyNear(E,st,mode,isdata);
   }
   else if (det==Detector::kFar) {
      ecor = CorrectShowerEnergyFar(E,st,mode,isdata);
   }

   return ecor;
}


float EnergyCorrections::CorrectMomentumFromRange_Birch(float p, bool isdata, 
                                                        Detector::Detector_t det)
{
  static const float c[4]={1.01334,0.05563,-0.05346,0.01205};

  // correction for difference in data mc steel density
  if (isdata){
    // inital correction, pre-Oxford 2006
    //static const float dcor=7.755/7.87;// data/mc density
    float dcor=1;
    if (det==Detector::kNear) dcor=(7.85*2.563)/(7.87*2.54);
    else if (det==Detector::kFar) dcor=(7.85*2.558)/(7.87*2.54);

    p*=dcor;
  }
  //
  float pcor=p/(c[0] + c[1]*log(p) + c[2]*pow(log(p), 2) + c[3]*pow(log(p),3));
  return pcor;
}

float EnergyCorrections::CorrectSignedMomentumFromCurvature_Birch(float p, bool /* isdata */, Detector::Detector_t /* det */)
{
  // input is the signed momentum (e.g. p/q)
  // isdata and det are not used... but maybe in the future
  float pcor=p;
  if (pcor!=0) {
    // correction advertised in 1430-v2, J. Musser
    // note: for 1/p < 0   C=1, so correction only matters for mu+
    float C = (1.01+0.1*fabs(1/p))/(1.01-0.1*(1/p));
    pcor*=(1.0/C);
  }
  return pcor;
}

float EnergyCorrections::CorrectShowerEnergyNear_Birch(float E, const CandShowerHandle::ShowerType_t& st, int mode, bool /* isdata */)
{

   //   std::cout << "CorrectShowerEnergyNear_Birch: " << E << "\t" << mode <<  std::endl;
   float ecor=E;
   if (st==CandShowerHandle::kCC){
      if (mode==1){
     // Niki Correction
     ecor=E/1.18;
      }
      else if (mode==2){
     // Andy Correction
     ecor=((E/1.05)*(1.-0.35*exp(-0.18*E/1.06)));
      }
   }
   else if (st==CandShowerHandle::kWtCC){
      if (mode==1){
     // Niki Correction
     ecor=((E)*(1.+0.50*exp(-1.00*E)));
      }
      else if (mode==2){
     // Andy Correction
     ecor=E/1.03;
      }
   }
   return ecor;
}

float EnergyCorrections::CorrectShowerEnergyFar_Birch(float E, const CandShowerHandle::ShowerType_t& st, int mode, bool /*isdata*/)
{
   //   std::cout << "CorrectShowerEnergyFar_Birch: " << E << std::endl;
   float ecor=E;
   if (st==CandShowerHandle::kCC){
     if (mode==1){
       // Niki Correction
       ecor=((E)*(1.-0.12*exp(-0.12*E)));
     }
     else if (mode==2){
       // Andy Correction
       //     ecor=(E)*(1.-0.2*exp(-0.2*E));
       ecor=E*(0.99-0.035*E*exp(-0.25*E));
     }
   }
   else if (st==CandShowerHandle::kWtCC){
      if (mode==1){
     // Niki Correction
     ecor=((E)*(1.+0.18*exp(-0.35*E)));
      }
      else if (mode==2){
     // Andy Correction
     E=ecor;
      }
   }
   return ecor;
}



float EnergyCorrections::CorrectMomentumFromRange_Cedar(float p, bool /* isdata */, Detector::Detector_t /* det */)
{
  return p;
}

float EnergyCorrections::CorrectSignedMomentumFromCurvature_Cedar(float p, bool /* isdata */, Detector::Detector_t /* det */)
{
  return p;
}

float EnergyCorrections::CorrectShowerEnergyNear_Cedar(float E, const CandShowerHandle::ShowerType_t&  st , int /* mode */, bool /* isdata */)
{
   float ecor=E;
   if (st==CandShowerHandle::kCC){
     ecor = ecor*(0.921+0.231*exp(-ecor*1.63));
   }
   else if (st==CandShowerHandle::kWtCC){
     ecor = ecor*(0.924+0.235*exp(-ecor*1.63));
   }
   return ecor;
}

float EnergyCorrections::CorrectShowerEnergyFar_Cedar(float E, const CandShowerHandle::ShowerType_t&  st , int /* mode */, bool /* isdata */)
{
   float ecor=E;
   if (st==CandShowerHandle::kCC){
     ecor = ecor*(0.950+0.277*exp(-ecor*1.64));
   }
   else if (st==CandShowerHandle::kWtCC){
     ecor = ecor*(0.957+0.271*exp(-ecor*1.64));
   }
   return ecor;
}

//for Dogwood reconstruction
float EnergyCorrections::ShowerEnergyCorrectionNearDogwood(float energy,
                                   const CandShowerHandle::ShowerType_t& st,
                                   EnergyCorrections::WhichCorrection_t whichCor)
{
  switch(whichCor) {
  case EnergyCorrections::kVersion4:
    return EnergyCorrections::MasakiNearJune30_2009(energy,st);
  case EnergyCorrections::kDefault:
  default:
    return EnergyCorrections::MasakiNearJune30_2009(energy,st);
  }
  return energy;
}

float EnergyCorrections::ShowerEnergyCorrectionFarDogwood(float energy,
                                const CandShowerHandle::ShowerType_t& st,
                                EnergyCorrections::WhichCorrection_t whichCor)
{
  switch(whichCor) {
  case EnergyCorrections::kVersion4:
    return EnergyCorrections::MasakiFarJune30_2009(energy,st);
  case EnergyCorrections::kDefault:
  default:
    return EnergyCorrections::MasakiFarJune30_2009(energy,st);
  }
  return energy;
}

float EnergyCorrections::ShowerEnergyCorrectionNearDogwood0(float energy,
                                                           const CandShowerHandle::ShowerType_t& st,
                                                           EnergyCorrections::WhichCorrection_t whichCor)
{
  switch(whichCor) {
  case EnergyCorrections::kVersion4:
    return EnergyCorrections::MasakiNear_forDogwood0(energy,st);
  case EnergyCorrections::kDefault:
  default:
    return EnergyCorrections::MasakiNear_forDogwood0(energy,st);
  }
  return energy;
}

float EnergyCorrections::ShowerEnergyCorrectionFarDogwood0(float energy,
                                                          const CandShowerHandle::ShowerType_t& st,
                                                          EnergyCorrections::WhichCorrection_t whichCor)
{
  switch(whichCor) {
  case EnergyCorrections::kVersion4:
    return EnergyCorrections::MasakiFar_forDogwood0(energy,st);
  case EnergyCorrections::kDefault:
  default:
    return EnergyCorrections::MasakiFar_forDogwood0(energy,st);
  }
  return energy;
}


//


// for Cedar_phy_bhcurv reconstruction
float EnergyCorrections::ShowerEnergyCorrectionNearCedarPhyBhcurve(float energy,
                             const CandShowerHandle::ShowerType_t& st,
                             EnergyCorrections::WhichCorrection_t whichCor)
{
   switch(whichCor) {
   case EnergyCorrections::kVersion4:
      return EnergyCorrections::MasakiNearMay17thScaled(energy,st);
   case EnergyCorrections::kDefault:
   default:
      return EnergyCorrections::MasakiNearDec15th(energy,st);
   }
   return energy;
}

float EnergyCorrections::ShowerEnergyCorrectionFarCedarPhyBhcurve(float energy,
                             const CandShowerHandle::ShowerType_t& st,
                             EnergyCorrections::WhichCorrection_t whichCor)
{
   switch(whichCor) {
   case EnergyCorrections::kVersion4:
      return EnergyCorrections::MasakiFarMay17thScaled(energy,st);
   case EnergyCorrections::kDefault:
   default:
      return EnergyCorrections::MasakiFarDec15th(energy,st);
   }
   return energy;
}
//

// for cedar reconstruction
float EnergyCorrections::ShowerEnergyCorrectionNearCedar(float energy,
                             const CandShowerHandle::ShowerType_t& st,
                             EnergyCorrections::WhichCorrection_t whichCor)
{
   switch(whichCor) {
   case EnergyCorrections::kVersion2:
      return EnergyCorrections::CorrectShowerEnergyNear_Cedar(energy,st);
   case EnergyCorrections::kVersion1:
      return EnergyCorrections::MasakiNearMay17th(energy,st);
   case EnergyCorrections::kVersion3:
   case EnergyCorrections::kDefault:
   default:
      return EnergyCorrections::MasakiNearMay17thScaled(energy,st);
   }
   return energy;
}

float EnergyCorrections::ShowerEnergyCorrectionFarCedar(float energy,
                             const CandShowerHandle::ShowerType_t& st,
                             EnergyCorrections::WhichCorrection_t whichCor)
{
   switch(whichCor) {
   case EnergyCorrections::kVersion2:
      return EnergyCorrections::CorrectShowerEnergyFar_Cedar(energy,st);
   case EnergyCorrections::kVersion1:
      return EnergyCorrections::MasakiFarMay17th(energy,st);
   case EnergyCorrections::kVersion3:
   case EnergyCorrections::kDefault:
   default:
      return EnergyCorrections::MasakiFarMay17thScaled(energy,st);
   }
   return energy;
}
//


float EnergyCorrections::MasakiNearJune30_2009(float energy,
                                               const CandShowerHandle::ShowerType_t& st)
{
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Near Detector from DocDB XXXX_v4\n";
  float recoE=energy;
  float le = log10(fmin(fmax(energy,0.3),20));
  float we=0;
  if (st==CandShowerHandle::kCC){
    recoE = energy*(0.96922+0.175773*le -0.0684406*(2.*pow(le,2)-1)+0.00940122*(4.*pow(le,3)-3.*le));
  }
  //Warning: weight shower energy is not used after Cedar_phy_bhcurv
  else if (st==CandShowerHandle::kWtCC) {
    we = fmin(fmax(energy,0.3),20);
    recoE = energy*(0.96922+0.175773*le -0.0684406*(2.*pow(le,2)-1)+0.00940122*(4.*pow(le,3)-3.*le));
  }
  return recoE;
}


float EnergyCorrections::MasakiFarJune30_2009(float energy,
                           const CandShowerHandle::ShowerType_t& st)
{
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Far Detector from DocDB XXXX_v4\n";
  float recoE=energy;
  float le = log10(fmin(fmax(energy,0.3),20));
  float we=0;
  if (st==CandShowerHandle::kCC){
    recoE = energy*(1.01397 + 0.0646697*le -0.0258817*(2.*pow(le,2)-1) + 0.00117583*(4.*pow(le,3)-3.*le));
  }
  //Warning: weight shower energy is not used after Cedar_phy_bhcurv
  else if (st==CandShowerHandle::kWtCC) {
    we = fmin(fmax(energy,0.3),20);
    recoE= energy*(1.01397 + 0.0646697*le -0.0258817*(2.*pow(le,2)-1) + 0.00117583*(4.*pow(le,3)-3.*le));
  }
  return recoE;
}

// for Dogwood0 below
float EnergyCorrections::MasakiNear_forDogwood0(float energy,
                                               const CandShowerHandle::ShowerType_t& st)
{
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Near Detector Dogwood0\n";
  float recoE=energy;
  float le = log10(fmin(fmax(energy,0.3),20));
  float we=0;
  if (st==CandShowerHandle::kCC){
    recoE = energy*(0.978739 + 0.156093*le + -0.0608388*(2.*pow(le,2)-1) + 0.00818386*(4.*pow(le,3)-3.*le));
  }
  //Warning: weight shower energy is not used after Cedar_phy_bhcurv                                                                                         
  else if (st==CandShowerHandle::kWtCC) {
    we = fmin(fmax(energy,0.3),20);
    recoE = energy*(0.978739 + 0.156093*le + -0.0608388*(2.*pow(le,2)-1) + 0.00818386*(4.*pow(le,3)-3.*le));
  }
  return recoE;
}

float EnergyCorrections::MasakiFar_forDogwood0(float energy,
                                              const CandShowerHandle::ShowerType_t& st)
{
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Far Detector Dogwood0\n";
  float recoE=energy;
  float le = log10(fmin(fmax(energy,0.3),20));
  float we=0;
  if (st==CandShowerHandle::kCC){
    recoE = energy*(1.02473 + 0.0429276*le + -0.016319*(2.*pow(le,2)-1) + -0.000380156*(4.*pow(le,3)-3.*le));
  }
  //Warning: weight shower energy is not used after Cedar_phy_bhcurv                                                                                         
  else if (st==CandShowerHandle::kWtCC) {
    we = fmin(fmax(energy,0.3),20);
    recoE= energy*(1.02473 + 0.0429276*le + -0.016319*(2.*pow(le,2)-1) + -0.000380156*(4.*pow(le,3)-3.*le));
  }
  return recoE;
}
//



float EnergyCorrections::MasakiNearDec15th(float energy,
                       const CandShowerHandle::ShowerType_t& st)
{
   //From DocDB 3077_v3
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Near Detector from DocDB 3895_v4\n";
  float recoE=energy;
  float le = log10(fmin(fmax(energy,0.3),20));
  float we=0;
  if (st==CandShowerHandle::kCC){
    recoE = energy*(1.10973-0.248714*le +0.116769*(2.*pow(le,2)-1)-0.0200268*(4.*pow(le,3)-3.*le));
  }
  else if (st==CandShowerHandle::kWtCC) {
    we = fmin(fmax(energy,0.3),20);
    recoE= energy*(0.999461-0.00334628*we+0.0000563316*pow(we,2)+0.35232*TMath::Exp(-we*1.76594));
  }
  return recoE;
}

float EnergyCorrections::MasakiFarDec15th(float energy,
                       const CandShowerHandle::ShowerType_t& st)
{
   //From DocDB 3077_v3
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Far Detector from DocDB 3895_v4\n";
  float recoE=energy;
  float le = log10(fmin(fmax(energy,0.3),20));
  float we=0;
  if (st==CandShowerHandle::kCC){
    recoE = energy*( 1.15566-0.286445*le+ 0.122705*(2.*pow(le,2)-1)-0.0183855*(4.*pow(le,3)-3.*le));
  }
  else if (st==CandShowerHandle::kWtCC) {
    we = fmin(fmax(energy,0.3),20);
    recoE= energy*(0.971346+0.00314063*we-0.000135242*pow(we,2)+0.626512*TMath::Exp(-we*3.26053));
  }
  return recoE;
}



float EnergyCorrections::MasakiNearDec15thScaled(float energy,
                       const CandShowerHandle::ShowerType_t& st)
{
   //From DocDB 3077_v3
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Near Detector from DocDB 3895_v4\n";
  float tempE=energy/EnergyCorrections::cgffCedarPhyDaikonND;
  float recoE=tempE;
  float le = log10(fmin(fmax(tempE,0.3),20));
  float we=0;
  if (st==CandShowerHandle::kCC){
    recoE = tempE*(1.10973-0.248714*le +0.116769*(2.*pow(le,2)-1)-0.0200268*(4.*pow(le,3)-3.*le));
  }
  else if (st==CandShowerHandle::kWtCC) {
    we = fmin(fmax(energy,0.3),20);
    recoE= energy*(0.999461-0.00334628*we+0.0000563316*pow(we,2)+0.35232*TMath::Exp(-we*1.76594));
  }
  return recoE;
}

float EnergyCorrections::MasakiFarDec15thScaled(float energy,
                       const CandShowerHandle::ShowerType_t& st)
{
   //From DocDB 3077_v3
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Far Detector from DocDB 3895_v4\n";
  float tempE=energy/EnergyCorrections::cgffCedarPhyDaikonFD;
  float recoE=tempE;
  float le = log10(fmin(fmax(tempE,0.3),20));
  float we=0;
  if (st==CandShowerHandle::kCC){
    recoE = tempE*( 1.15566-0.286445*le+ 0.122705*(2.*pow(le,2)-1)-0.0183855*(4.*pow(le,3)-3.*le));
  }
  else if (st==CandShowerHandle::kWtCC) {
    we = fmin(fmax(energy,0.3),20);
    recoE= energy*(0.971346+0.00314063*we-0.000135242*pow(we,2)+0.626512*TMath::Exp(-we*3.26053));
  }
  return recoE;
}


float EnergyCorrections::MasakiNearMay17th(float energy,
                       const CandShowerHandle::ShowerType_t& st)
{
   //From DocDB 3077_v3
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Near Detector from DocDB 3077_v3\n";
  float recoE=energy;
   float le = log10(fmax(energy,0.2));
   if (st==CandShowerHandle::kCC){
      recoE = energy*(1.078984-0.249843*le+0.134518*(2.*pow(le,2)-1)-0.025613*(4.*pow(le,3)-3.*le));
   }
   else if (st==CandShowerHandle::kWtCC) {
      recoE= energy*(1.070553-0.207148*le+0.0943124*(2.*pow(le,2)-1)-0.0153231*(4.*pow(le,3)-3.*le));
   }
   return recoE;
}

float EnergyCorrections::MasakiFarMay17th(float energy,
                       const CandShowerHandle::ShowerType_t& st)
{
   //From DocDB 3077_v3
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Far Detector from DocDB 3077_v3\n";
   float recoE=energy;
   //   std::cout << "Old:\t" << recoE << std::endl;
   float le = log10(fmax(energy,0.2));
   if (st==CandShowerHandle::kCC){
      recoE = energy*(1.113584-0.299139*le+0.145169*(2.*pow(le,2)-1)-0.0232853*(4.*pow(le,3)-3.*le));
   }
   else if (st==CandShowerHandle::kWtCC){
      recoE= energy*(1.052872-0.19185*le+0.102449*(2.*pow(le,2)-1)-0.0182317*(4.*pow(le,3)-3.*le));
   }
   return recoE;
}



float EnergyCorrections::MasakiNearMay17thScaled(float energy,
                           const CandShowerHandle::ShowerType_t& st)
{
   //From DocDB 3077_v3
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Near Detector from DocDB 3077_v3 -- Scaled by Appropriate CG factor for uncalibrated data he used\n";
   float tempE=energy/EnergyCorrections::cgffCedarR1_24_1MCND;
   float recoE=tempE;
   float le = log10(fmax(tempE,0.2));
   if (st==CandShowerHandle::kCC){
      recoE = tempE*(1.078984-0.249843*le+0.134518*(2.*pow(le,2)-1)-0.025613*(4.*pow(le,3)-3.*le));
   }
   else if (st==CandShowerHandle::kWtCC) {
      recoE= tempE*(1.070553-0.207148*le+0.0943124*(2.*pow(le,2)-1)-0.0153231*(4.*pow(le,3)-3.*le));
   }
   return recoE;
}


float EnergyCorrections::MasakiFarMay17thScaled(float energy,
                        const CandShowerHandle::ShowerType_t& st)
{
   //From DocDB 3077_v3
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Masaki (Reco to True) Correction for Far Detector from DocDB 3077_v3 -- Scaled by Appropriate CG factor for uncalibrated data he used\n";
   float tempE=energy/EnergyCorrections::cgffCedarR1_24_1MCFD;
   float recoE=tempE;
   //   std::cout <<  "New:\t" <<recoE << std::endl;
   float le = log10(fmax(tempE,0.2));
   if (st==CandShowerHandle::kCC){
      recoE = tempE*(1.113584-0.299139*le+0.145169*(2.*pow(le,2)-1)-0.0232853*(4.*pow(le,3)-3.*le));
   }
   else if (st==CandShowerHandle::kWtCC){
      recoE= tempE*(1.052872-0.19185*le+0.102449*(2.*pow(le,2)-1)-0.0182317*(4.*pow(le,3)-3.*le));
   }
   return recoE;
}

float EnergyCorrections::MomentumRangeCorrectionBirch(float p, VldContext vc, EnergyCorrections::WhichCorrection_t /*whichCor*/)
{
  static const float c[4]={1.01334,0.05563,-0.05346,0.01205};

  // correction for difference in data mc steel density
  if (vc.GetSimFlag()==SimFlag::kData){
    // inital correction, pre-Oxford 2006
    //static const float dcor=7.755/7.87;// data/mc density
    float dcor=1;
    if (vc.GetDetector()==Detector::kNear) dcor=(7.85*2.563)/(7.87*2.54);
    else if (vc.GetDetector()==Detector::kFar) dcor=(7.85*2.558)/(7.87*2.54);

    p*=dcor;
  }
  //
  float pcor=p/(c[0] + c[1]*log(p) + c[2]*pow(log(p), 2) + c[3]*pow(log(p),3));
  return pcor;
}

float EnergyCorrections::SignedMomentumCurvatureCorrectionBirch(float p, VldContext /*vc*/, WhichCorrection_t /*whichCor*/){
  // input is the signed momentum (e.g. p/q)
  // isdata and det are not used... but maybe in the future
  float pcor=p;
  if (pcor!=0) {
    // correction advertised in 1430-v2, J. Musser
    // note: for 1/p < 0   C=1, so correction only matters for mu+
    float C = (1.01+0.1*fabs(1/p))/(1.01-0.1*(1/p));
    pcor*=(1.0/C);
  }
  return pcor;
}


float EnergyCorrections::MomentumRangeCorrectionCedar(float p, VldContext vc, WhichCorrection_t whichCor)
{
  //return p;
    const float m=0.1057;// muon mass
    float E = sqrt(p*p+m*m);
    float eCor = EnergyCorrections::EnergyRangeCorrectionCedar(E,vc,whichCor);
    return sqrt(eCor*eCor-m*m);
}

float EnergyCorrections::EnergyRangeCorrectionCedar(float E, VldContext vc, WhichCorrection_t /*whichCor*/)
{
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Applying Energy from Range Correction for Cedar (1.018*E)-0.009 for ND Data and (1.010*E)-0.009 for everything else.\n";
  float eCor=E;
  if (vc.GetSimFlag()==SimFlag::kData && vc.GetDetector()==Detector::kNear) {
    eCor=(1.018*E)-0.009;
  }
  else {
    eCor=(1.010*E)-0.009;
  }

  return eCor;
}

float EnergyCorrections::SignedMomentumCurvatureCorrectionCedar(float p, VldContext /*vc*/, WhichCorrection_t /*whichCor*/)
{
  MAXMSG("DataUtil",Msg::kInfo,1)
    << "EnergyCorrections -- Not applying momentum from curvature correction for Cedar\n";
  return p;
}

---