#include <G4UIdirectory.hh>
#include <G4UIcmdWithADoubleAndUnit.hh>
#include <G4UIcmdWithABool.hh>
#include <G4UIcmdWithAString.hh>
#include <G4EmCalculator.hh>

#include <RAT/PhysicsListMessenger.hh>
#include <RAT/PhysicsList.hh>
#include <RAT/Log.hh>

#include <CLHEP/Units/SystemOfUnits.h>
#include <sstream>
using namespace CLHEP;
using namespace std;

RAT::PhysicsListMessenger::PhysicsListMessenger( PhysicsList* physicsList )
  :fPhysicsList( physicsList )
{
  fPhysDir = new G4UIdirectory("/rat/physics_list/");
  fPhysDir->SetGuidance("Physics list commands");

  fOmitAllCmd = new G4UIcmdWithABool("/rat/physics_list/OmitAll",this);
  fOmitAllCmd->SetGuidance("Omits all physics. NOT for simulation ");
  fOmitAllCmd->SetParameterName("OmitAll",false);

  fOmitOpticalBoundaryEffects = new G4UIcmdWithABool("/rat/physics_list/Optical/OmitBoundaryEffects",this);
  fOmitOpticalBoundaryEffects->SetGuidance("Omits optical boundary effects between different optical media i.e. reflection/refraction ");
  fOmitOpticalBoundaryEffects->SetParameterName("OmitBoundaryEffects",false);

  fOmitOpticalAttenuationCmd = new G4UIcmdWithABool("/rat/physics_list/Optical/OmitAttenuation",this);
  fOmitOpticalAttenuationCmd->SetGuidance("Omits optical absorption processes");
  fOmitOpticalAttenuationCmd->SetParameterName("OmitAttenuation",false);

  fOmitOpticalAbsorptionCmd = new G4UIcmdWithABool("/rat/physics_list/Optical/OmitAbsorption",this);
  fOmitOpticalAbsorptionCmd->SetGuidance("Omits optical absorption processes");
  fOmitOpticalAbsorptionCmd->SetParameterName("OmitAbsorption",false);

  fOmitOpticalRayleighCmd = new G4UIcmdWithABool("/rat/physics_list/Optical/OmitRayleigh",this);
  fOmitOpticalRayleighCmd->SetGuidance("Omits optical absorption processes");
  fOmitOpticalRayleighCmd->SetParameterName("OmitRayleigh",false);

  fOmitOpticalMieCmd = new G4UIcmdWithABool("/rat/physics_list/Optical/OmitMie",this);
  fOmitOpticalMieCmd->SetGuidance("Omits optical absorption processes");
  fOmitOpticalMieCmd->SetParameterName("OmitMie",false);

  fOmitCerenkovCmd = new G4UIcmdWithABool("/rat/physics_list/OmitCerenkov",this);
  fOmitCerenkovCmd->SetGuidance("Omits all Cerenkov processes");
  fOmitCerenkovCmd->SetParameterName("OmitCerenkov",false);

  fOmitHadronListCmd = new G4UIcmdWithABool("/rat/physics_list/OmitHadronicPhysicsList",this);
  fOmitHadronListCmd->SetGuidance("Omits the hadronic physics list (only)");
  fOmitHadronListCmd->SetParameterName("OmitHadronicPhysicsList",false);

  fOmitHadronicCmd = new G4UIcmdWithABool("/rat/physics_list/OmitHadronicProcesses",this);
  fOmitHadronicCmd->SetGuidance("Omita all Hadronic Processes (including EM lists).");
  fOmitHadronicCmd->SetParameterName("OmitHadronic",false);

  fOmitMuonicCmd = new G4UIcmdWithABool("/rat/physics_list/OmitMuonicProcesses",this);
  fOmitMuonicCmd->SetGuidance("Omit Muonic Processes.");
  fOmitMuonicCmd->SetParameterName("OmitMuonic",false);

  fTimeCutCmd = new G4UIcmdWithADoubleAndUnit("/rat/physics_list/SetMaxGlobalTimeCut",this);
  fTimeCutCmd->SetGuidance("Set maximum global time of a step, cut.");
  fTimeCutCmd->SetParameterName("TimeCut",false);
  fTimeCutCmd->SetUnitCategory("Time");
  fTimeCutCmd->SetRange("TimeCut>0.0");
  fTimeCutCmd->AvailableForStates(G4State_PreInit,G4State_Idle);

  fGammaCutCmd = new G4UIcmdWithADoubleAndUnit("/rat/physics_list/setGammaCut",this);
  fGammaCutCmd->SetGuidance("Set global gamma cut.");
  fGammaCutCmd->SetParameterName("GammaCut",false);
  fGammaCutCmd->SetUnitCategory("Length");
  fGammaCutCmd->SetRange("GammaCut>0.0");
  fGammaCutCmd->AvailableForStates(G4State_PreInit,G4State_Idle);

  fElectronCutCmd = new G4UIcmdWithADoubleAndUnit("/rat/physics_list/setElectronCut",this);
  fElectronCutCmd->SetGuidance("Set global electron cut.");
  fElectronCutCmd->SetParameterName("ElectronCut",false);
  fElectronCutCmd->SetUnitCategory("Length");
  fElectronCutCmd->SetRange("ElectronCut>0.0");
  fElectronCutCmd->AvailableForStates(G4State_PreInit,G4State_Idle);

  fMuonCutCmd = new G4UIcmdWithADoubleAndUnit("/rat/physics_list/setMuonCut",this);
  fMuonCutCmd->SetGuidance("Set global electron cut.");
  fMuonCutCmd->SetParameterName("MuonCut",false);
  fMuonCutCmd->SetUnitCategory("Length");
  fMuonCutCmd->SetRange("MuonCut>0.0");
  fMuonCutCmd->AvailableForStates(G4State_PreInit,G4State_Idle);

  fAllCutCmd = new G4UIcmdWithADoubleAndUnit("/rat/physics_list/setAllCuts",this);
  fAllCutCmd->SetGuidance("Set all global cuts.");
  fAllCutCmd->SetParameterName("cut",false);
  fAllCutCmd->SetUnitCategory("Length");
  fAllCutCmd->SetRange("cut>0.0");
  fAllCutCmd->AvailableForStates(G4State_PreInit,G4State_Idle);

  fEMListCmd = new G4UIcmdWithAString("/rat/physics_list/SelectEMPhysicsList",this);
  fEMListCmd->SetGuidance("Choose the EM Physics List, options: V1, V2, V3, V4, V5, standard, penelope, livermore.");
  fEMListCmd->SetParameterName("EMList",false);
  fEMListCmd->AvailableForStates(G4State_PreInit);

  fHadronicListCmd = new G4UIcmdWithAString("/rat/physics_list/SelectHadronicPhysicsList",this);
  fHadronicListCmd->SetGuidance("Choose the Hadronic Physics List, options: V1, V2, V3, V4");
  fHadronicListCmd->SetParameterName("HadronicList",false);
  fHadronicListCmd->AvailableForStates(G4State_PreInit);

  fCerenkovProcessCmd = new G4UIcmdWithAString("/rat/physics_list/SelectCerenkovProcess",this);
  fCerenkovProcessCmd->SetGuidance("Choose the Cerenkov Process, options: Default, SNOMAN.");
  fCerenkovProcessCmd->SetParameterName("CerenkovProcess",false);
  fCerenkovProcessCmd->AvailableForStates(G4State_PreInit);

  fOpRayleighProcessCmd = new G4UIcmdWithAString("/rat/physics_list/SelectOpRayleighProcess",this);
  fOpRayleighProcessCmd->SetGuidance("Choose the Optical Rayleigh Scattering Process, options: Default, SNOMAN.");
  fOpRayleighProcessCmd->SetParameterName("OpRayleighProcess",false);
  fOpRayleighProcessCmd->AvailableForStates(G4State_PreInit);

  fPrintEletronicDEDXCmd = new G4UIcommand( "/rat/physics_list/Print/ElectronicDEDX", this );
  fPrintEletronicDEDXCmd->SetGuidance( "Print the dE/dx for the specified particle in the specified material for electronic processes." );
  fPrintEletronicDEDXCmd->SetParameter( new G4UIparameter( "particle", 's', false ) );
  fPrintEletronicDEDXCmd->SetParameter( new G4UIparameter( "material", 's', false ) );

  fPrintTotalDEDXCmd = new G4UIcommand( "/rat/physics_list/Print/TotalDEDX", this );
  fPrintTotalDEDXCmd->SetGuidance( "Print the dE/dx for the specified particle in the specified material for all processes." );
  fPrintTotalDEDXCmd->SetParameter( new G4UIparameter( "particle", 's', false ) );
  fPrintTotalDEDXCmd->SetParameter( new G4UIparameter( "material", 's', false ) );

  fPrintMeanPathCmd = new G4UIcommand( "/rat/physics_list/Print/GammaMeanFreePath", this );
  fPrintMeanPathCmd->SetGuidance( "Print the gamma mean free path for the specified material." );
  fPrintMeanPathCmd->SetParameter( new G4UIparameter( "material", 's', false ) );
}

RAT::PhysicsListMessenger::~PhysicsListMessenger()
{
  delete fOmitAllCmd;
  delete fOmitOpticalBoundaryEffects;
  delete fOmitOpticalAttenuationCmd;
  delete fOmitOpticalAbsorptionCmd;
  delete fOmitOpticalRayleighCmd;
  delete fOmitOpticalMieCmd;
  delete fOmitCerenkovCmd;
  delete fOmitHadronListCmd;
  delete fOmitHadronicCmd;
  delete fOmitMuonicCmd;
  delete fTimeCutCmd;
  delete fGammaCutCmd;
  delete fElectronCutCmd;
  delete fMuonCutCmd;
  delete fAllCutCmd;
  delete fEMListCmd;
  delete fHadronicListCmd;
  delete fCerenkovProcessCmd;
  delete fOpRayleighProcessCmd;
  delete fPrintEletronicDEDXCmd;
  delete fPrintTotalDEDXCmd;
  delete fPrintMeanPathCmd;
  delete fPhysDir;
}

void RAT::PhysicsListMessenger::SetNewValue( G4UIcommand* command,
                                             G4String newValue )
{
  if( command == fOmitAllCmd )
    fPhysicsList->OmitProcesses( "all", fOmitAllCmd->GetNewBoolValue( newValue ) );
  else if( command == fTimeCutCmd )
    fPhysicsList->SetCut( "Time", fTimeCutCmd->GetNewDoubleValue( newValue ) );
  else if( command == fGammaCutCmd )
    fPhysicsList->SetCut( "Gamma", fGammaCutCmd->GetNewDoubleValue( newValue ) );
  else if( command == fElectronCutCmd )
    fPhysicsList->SetCut( "Electron", fElectronCutCmd->GetNewDoubleValue( newValue ) );
  else if( command == fMuonCutCmd )
    fPhysicsList->SetCut( "Muon", fMuonCutCmd->GetNewDoubleValue( newValue ) );
  else if( command == fAllCutCmd )
    fPhysicsList->SetCut( "All", fAllCutCmd->GetNewDoubleValue( newValue ) );
  else if( command == fEMListCmd )
    {
      if( newValue( 0, 3 ) == G4String( "RAT" ) )
        newValue = newValue( 3 );
      fPhysicsList->SelectEMPhysicsList( newValue );
    }
  else if( command == fHadronicListCmd )
    fPhysicsList->SelectHadronicPhysicsList( newValue );
  else if( command == fCerenkovProcessCmd )
    fPhysicsList->SelectCerenkovProcess( newValue );
  else if( command == fOpRayleighProcessCmd )
    fPhysicsList->SelectOpRayleighProcess( newValue );
  else if( command == fOmitCerenkovCmd )
    fPhysicsList->OmitProcesses("cerenkov",fOmitCerenkovCmd->GetNewBoolValue( newValue ));
  else if( command == fOmitHadronListCmd )
    fPhysicsList->OmitProcesses("hadronicListOnly",fOmitHadronListCmd->GetNewBoolValue( newValue ));
  else if( command == fOmitHadronicCmd )
    fPhysicsList->OmitProcesses("hadronic", fOmitHadronicCmd->GetNewBoolValue( newValue ));
  else if( command == fOmitMuonicCmd )
    fPhysicsList->OmitProcesses("muonic", fOmitMuonicCmd->GetNewBoolValue( newValue ));
  else if( command == fOmitOpticalAttenuationCmd )
    fPhysicsList->OmitProcesses("attenuation", fOmitOpticalAttenuationCmd->GetNewBoolValue( newValue ));
  else if ( command == fOmitOpticalAbsorptionCmd )
    fPhysicsList->OmitProcesses("absorption", fOmitOpticalAbsorptionCmd->GetNewBoolValue( newValue ));
  else if ( command == fOmitOpticalRayleighCmd )
    fPhysicsList->OmitProcesses("rayleigh", fOmitOpticalRayleighCmd->GetNewBoolValue( newValue ));
  else if ( command == fOmitOpticalMieCmd )
    fPhysicsList->OmitProcesses("mie", fOmitOpticalMieCmd->GetNewBoolValue( newValue ));
  else if ( command == fOmitOpticalBoundaryEffects )
    fPhysicsList->OmitProcesses("boundary", fOmitOpticalBoundaryEffects->GetNewBoolValue( newValue ));
  else if( command == fPrintEletronicDEDXCmd )
    {
      istringstream args(newValue);
      string particle, material;
      args >> particle >> material;
      G4EmCalculator emCalculator;
      const double energy[] = { 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0 };
      warn << "PhysicsListMessenger::SetNewValue: Printing electronic dE/dx for " << particle << " in " << material << ".\n{Energy [MeV], dEdx [mm/MeV]}" << newline;
      for( int iEnergy = 0; iEnergy < 19; iEnergy++ )
        warn << "{" << energy[iEnergy] << ", " << emCalculator.ComputeElectronicDEDX( energy[iEnergy], particle, material ) << "}, ";
      warn << newline;
    }
  else if( command == fPrintTotalDEDXCmd )
    {
      istringstream args(newValue);
      string particle, material;
      args >> particle >> material;
      G4EmCalculator emCalculator;
      const double energy[] = { 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0 };
      warn << "PhysicsListMessenger::SetNewValue: Printing total dE/dx for " << particle << " in " << material << ".\n{Energy [MeV], dEdx [mm/MeV]}" << newline;
      for( int iEnergy = 0; iEnergy < 19; iEnergy++ )
        warn << "{" << energy[iEnergy] << ", " << emCalculator.ComputeTotalDEDX( energy[iEnergy], particle, material ) << "}, ";
      warn << newline;
    }
  else if( command == fPrintMeanPathCmd )
    {
      const string material = newValue;
      G4EmCalculator emCalculator;
      const double energy[] = { 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0 };
      warn << "PhysicsListMessenger::SetNewValue: Printing gamma mean free paths in " << material <<
        "\n{Energy [MeV], photoelectron [mm], conversion [mm], compton [mm], rayleigh [mm]}" << newline;
      for( int iEnergy = 0; iEnergy < 19; iEnergy++ )
        {
          warn << "{" << energy[iEnergy] << ", ";
          warn << emCalculator.ComputeMeanFreePath( energy[iEnergy], "gamma", "phot", material ) << ", ";
          warn << emCalculator.ComputeMeanFreePath( energy[iEnergy], "gamma", "conv", material ) << ", ";
          warn << emCalculator.ComputeMeanFreePath( energy[iEnergy], "gamma", "compt", material ) << ", ";
          warn << emCalculator.ComputeMeanFreePath( energy[iEnergy], "gamma", "Rayl", material ) << "}, ";
        }
      warn << newline;
    }
  else
  {
    RAT::Log::Die("PhysicsListMessenger::SetNewValue: ERROR > Invalid PhysicsListMessenger \"PhysicsList\" command ");
  }
}