PassiveProtonBeamLine.cc

Go to the documentation of this file.

//
// ********************************************************************
// * License and Disclaimer                                           *
// *                                                                  *
// * The  Geant4 software  is  copyright of the Copyright Holders  of *
// * the Geant4 Collaboration.  It is provided  under  the terms  and *
// * conditions of the Geant4 Software License,  included in the file *
// * LICENSE and available at  http://cern.ch/geant4/license .  These *
// * include a list of copyright holders.                             *
// *                                                                  *
// * Neither the authors of this software system, nor their employing *
// * institutes,nor the agencies providing financial support for this *
// * work  make  any representation or  warranty, express or implied, *
// * regarding  this  software system or assume any liability for its *
// * use.  Please see the license in the file  LICENSE  and URL above *
// * for the full disclaimer and the limitation of liability.         *
// *                                                                  *
// * This  code  implementation is the result of  the  scientific and *
// * technical work of the GEANT4 collaboration.                      *
// * By using,  copying,  modifying or  distributing the software (or *
// * any work based  on the software)  you  agree  to acknowledge its *
// * use  in  resulting  scientific  publications,  and indicate your *
// * acceptance of all terms of the Geant4 Software license.          *
// ********************************************************************
//
// Hadrontherapy (both basic and full version) are supported by the Italian INFN
// Institute in the framework of the MC-INFN Group
//
#include "globals.hh"
#include "G4SystemOfUnits.hh"
#include "G4Box.hh"
#include "G4Tubs.hh"
#include "G4VisAttributes.hh"
#include "G4Colour.hh"
#include "G4RunManager.hh"
#include "G4LogicalVolume.hh"
#include "G4PVPlacement.hh"
#include "G4RotationMatrix.hh"
#include "G4NistManager.hh"
#include "G4NistElementBuilder.hh"
#include "HadrontherapyDetectorConstruction.hh"
#include "HadrontherapyModulator.hh"
#include "PassiveProtonBeamLine.hh"
#include "PassiveProtonBeamLineMessenger.hh"
//#include "FaradayCup.hh"

//G4bool PassiveProtonBeamLine::doCalculation = false;
PassiveProtonBeamLine::PassiveProtonBeamLine():
modulator(0), physicalTreatmentRoom(0),hadrontherapyDetectorConstruction(0),
physiBeamLineSupport(0), physiBeamLineCover(0), physiBeamLineCover2(0),
firstScatteringFoil(0), physiFirstScatteringFoil(0), physiKaptonWindow(0),
solidStopper(0), physiStopper(0), secondScatteringFoil(0), physiSecondScatteringFoil(0),
physiFirstCollimator(0), solidRangeShifterBox(0), logicRangeShifterBox(0),
physiRangeShifterBox(0), physiSecondCollimator(0), physiFirstCollimatorModulatorBox(0),
physiHoleFirstCollimatorModulatorBox(0), physiSecondCollimatorModulatorBox(0),
physiHoleSecondCollimatorModulatorBox(0), physiMOPIMotherVolume(0),
physiFirstMonitorLayer1(0), physiFirstMonitorLayer2(0), physiFirstMonitorLayer3(0),
physiFirstMonitorLayer4(0), physiSecondMonitorLayer1(0), physiSecondMonitorLayer2(0),
physiSecondMonitorLayer3(0), physiSecondMonitorLayer4(0), physiNozzleSupport(0), //physiHoleNozzleSupport(0),
physiBrassTube(0), solidFinalCollimator(0), physiFinalCollimator(0)
{
    // Messenger to change parameters of the passiveProtonBeamLine geometry
    passiveMessenger = new PassiveProtonBeamLineMessenger(this);
    
    //***************************** PW ***************************************
    static G4String ROGeometryName = "DetectorROGeometry";
    RO = new HadrontherapyDetectorROGeometry(ROGeometryName);
    
    G4cout << "Going to register Parallel world...";
    RegisterParallelWorld(RO);
    G4cout << "... done" << G4endl;
    //***************************** PW ***************************************
    /* // For the Faraday Cup activation
    if (name)
    {
        doCalculation = true;
        }*/
}
PassiveProtonBeamLine::~PassiveProtonBeamLine()
{
    delete passiveMessenger;
    delete hadrontherapyDetectorConstruction;
    /* if (!pFaradayCup)
    {
        delete pFaradayCup;
        }*/
}

G4VPhysicalVolume* PassiveProtonBeamLine::Construct()
{
    // Sets default geometry and materials
    SetDefaultDimensions();
    
    // Construct the whole Passive Beam Line
    ConstructPassiveProtonBeamLine();
    
    //***************************** PW ***************************************
    if (!hadrontherapyDetectorConstruction)
        
        //***************************** PW ***************************************
        
        // HadrontherapyDetectorConstruction builds ONLY the phantom and the detector with its associated ROGeometry
        hadrontherapyDetectorConstruction = new HadrontherapyDetectorConstruction(physicalTreatmentRoom);
    
    
    //***************************** PW ***************************************
    
    hadrontherapyDetectorConstruction->InitializeDetectorROGeometry(RO,hadrontherapyDetectorConstruction->GetDetectorToWorldPosition());
    
    //***************************** PW ***************************************
    /*  if (doCalculation)
    {
        pFaradayCup = new FaradayCup(physicalTreatmentRoom);
        G4cout << "Faraday Cup detector has been activated" << G4endl;
        }*/
    return physicalTreatmentRoom;
}

// In the following method the DEFAULTS used in the geometry of
// passive beam line are provided
// HERE THE USER CAN CHANGE THE GEOMETRY CHARACTERISTICS OF BEAM
// LINE ELEMENTS, ALTERNATIVELY HE/SHE CAN USE THE MACRO FILE (IF A
// MESSENGER IS PROVIDED)
//
// DEFAULT MATERIAL ARE ALSO PROVIDED
// and COLOURS ARE ALSO DEFINED
// ----------------------------------------------------------
void PassiveProtonBeamLine::SetDefaultDimensions()
{
    // Set of coulors that can be used
    white = new G4VisAttributes( G4Colour());
    white -> SetVisibility(true);
    white -> SetForceSolid(true);
        
    blue = new G4VisAttributes(G4Colour(0. ,0. ,1.));
    blue -> SetVisibility(true);
    blue -> SetForceSolid(true);
        
    gray = new G4VisAttributes( G4Colour(0.5, 0.5, 0.5 ));
    gray-> SetVisibility(true);
    gray-> SetForceSolid(true);
        
    red = new G4VisAttributes(G4Colour(1. ,0. ,0.));
    red-> SetVisibility(true);
    red-> SetForceSolid(true);
        
    yellow = new G4VisAttributes(G4Colour(1., 1., 0. ));
    yellow-> SetVisibility(true);
    yellow-> SetForceSolid(true);
        
    green = new G4VisAttributes( G4Colour(25/255. , 255/255. ,  25/255. ));
    green -> SetVisibility(true);
    green -> SetForceSolid(true);
        
    darkGreen = new G4VisAttributes( G4Colour(0/255. , 100/255. ,  0/255. ));
    darkGreen -> SetVisibility(true);
    darkGreen -> SetForceSolid(true);
    
    darkOrange3 = new G4VisAttributes( G4Colour(205/255. , 102/255. ,  000/255. ));
    darkOrange3 -> SetVisibility(true);
    darkOrange3 -> SetForceSolid(true);
        
    skyBlue = new G4VisAttributes( G4Colour(135/255. , 206/255. ,  235/255. ));
    skyBlue -> SetVisibility(true);
    skyBlue -> SetForceSolid(true);
    
    
    // VACUUM PIPE: first track of the beam line is inside vacuum;
    // The PIPE contains the FIRST SCATTERING FOIL and the KAPTON WINDOW
    G4double defaultVacuumZoneXSize = 100.0 *mm;
    vacuumZoneXSize = defaultVacuumZoneXSize;
    
    G4double defaultVacuumZoneYSize = 52.5 *mm;
    vacuumZoneYSize = defaultVacuumZoneYSize;
    
    G4double defaultVacuumZoneZSize = 52.5 *mm;
    vacuumZoneZSize = defaultVacuumZoneZSize;
    
    G4double defaultVacuumZoneXPosition = -3010.0 *mm;
    vacuumZoneXPosition = defaultVacuumZoneXPosition;
    
    // FIRST SCATTERING FOIL: a thin foil performing a first scattering
    // of the original beam
    G4double defaultFirstScatteringFoilXSize = 0.0075 *mm;
    firstScatteringFoilXSize = defaultFirstScatteringFoilXSize;
    
    G4double defaultFirstScatteringFoilYSize = 52.5   *mm;
    firstScatteringFoilYSize = defaultFirstScatteringFoilYSize;
    
    G4double defaultFirstScatteringFoilZSize = 52.5   *mm;
    firstScatteringFoilZSize = defaultFirstScatteringFoilZSize;
    
    G4double defaultFirstScatteringFoilXPosition = 0.0 *mm;
    firstScatteringFoilXPosition = defaultFirstScatteringFoilXPosition;
    
    // KAPTON WINDOW: it prmits the passage of the beam from vacuum to air
    G4double defaultKaptonWindowXSize = 0.010*mm;
    kaptonWindowXSize = defaultKaptonWindowXSize;
    
    G4double defaultKaptonWindowYSize = 5.25*cm;
    kaptonWindowYSize = defaultKaptonWindowYSize;
    
    G4double defaultKaptonWindowZSize = 5.25*cm;
    kaptonWindowZSize = defaultKaptonWindowZSize;
    
    G4double defaultKaptonWindowXPosition = 100.0*mm - defaultKaptonWindowXSize;
    kaptonWindowXPosition = defaultKaptonWindowXPosition;
    
    // STOPPER: is a small cylinder able to stop the central component
    // of the beam (having a gaussian shape). It is connected to the SECON SCATTERING FOIL
    // and represent the second element of the scattering system
    G4double defaultInnerRadiusStopper = 0.*cm;
    innerRadiusStopper = defaultInnerRadiusStopper;
    
    G4double defaultHeightStopper = 3.5*mm;
    heightStopper = defaultHeightStopper;
    
    G4double defaultStartAngleStopper = 0.*deg;
    startAngleStopper = defaultStartAngleStopper;
    
    G4double defaultSpanningAngleStopper = 360.*deg;
    spanningAngleStopper = defaultSpanningAngleStopper;
    
    G4double defaultStopperXPosition = -2705.0 *mm;
    stopperXPosition = defaultStopperXPosition;
    
    G4double defaultStopperYPosition = 0.*m;
    stopperYPosition = defaultStopperYPosition;
    
    G4double defaultStopperZPosition = 0.*m;
    stopperZPosition = defaultStopperZPosition;
    
    G4double defaultOuterRadiusStopper = 2 *mm;
    outerRadiusStopper = defaultOuterRadiusStopper;
    
    // SECOND SCATTERING FOIL: it is another thin foil and provides the
    // final diffusion of the beam. It represents the third element of the scattering
    // system;
    G4double defaultSecondScatteringFoilXSize = 0.0125 *mm;
    secondScatteringFoilXSize = defaultSecondScatteringFoilXSize;
    
    G4double defaultSecondScatteringFoilYSize = 52.5   *mm;
    secondScatteringFoilYSize = defaultSecondScatteringFoilYSize;
    
    G4double defaultSecondScatteringFoilZSize = 52.5   *mm;
    secondScatteringFoilZSize = defaultSecondScatteringFoilZSize;
    
    G4double defaultSecondScatteringFoilXPosition = defaultStopperXPosition + defaultHeightStopper + defaultSecondScatteringFoilXSize;
    secondScatteringFoilXPosition = defaultSecondScatteringFoilXPosition;
    
    G4double defaultSecondScatteringFoilYPosition =  0 *mm;
    secondScatteringFoilYPosition = defaultSecondScatteringFoilYPosition;
    
    G4double defaultSecondScatteringFoilZPosition =  0 *mm;
    secondScatteringFoilZPosition = defaultSecondScatteringFoilZPosition;
    
    // RANGE SHIFTER: is a slab of PMMA acting as energy degreader of
    // primary beam
    
    //Default material of the range shifter
    
    G4double defaultRangeShifterXSize = 5. *mm;
    rangeShifterXSize = defaultRangeShifterXSize;
    
    G4double defaultRangeShifterYSize = 176. *mm;
    rangeShifterYSize = defaultRangeShifterYSize;
    
    G4double defaultRangeShifterZSize = 176. *mm;
    rangeShifterZSize = defaultRangeShifterZSize;
    
    G4double defaultRangeShifterXPosition = -2393.0 *mm;
    rangeShifterXPosition = defaultRangeShifterXPosition;
    
    G4double defaultRangeShifterYPosition = 0. *mm;
    rangeShifterYPosition = defaultRangeShifterYPosition;
    
    G4double defaultRangeShifterZPosition = 0. *mm;
    rangeShifterZPosition = defaultRangeShifterZPosition;
    
    // MOPI DETECTOR: two orthogonal microstrip gas detectors developed
    // by the INFN Section of Turin in collaboration with some
    // of the author of this example. It permits the
    // on-line check of the beam simmetry via the signal
    // integration of the collected charge for each strip.
    
    // Mother volume of MOPI
    
    G4double defaultMOPIMotherVolumeXSize = 12127.0 *um;
    MOPIMotherVolumeXSize = defaultMOPIMotherVolumeXSize;
    
    G4double defaultMOPIMotherVolumeYSize = 40.0 *cm;
    MOPIMotherVolumeYSize = defaultMOPIMotherVolumeYSize;
    
    G4double defaultMOPIMotherVolumeZSize = 40.0 *cm;
    MOPIMotherVolumeZSize = defaultMOPIMotherVolumeZSize;
    
    G4double defaultMOPIMotherVolumeXPosition = -1000.0 *mm;
    MOPIMotherVolumeXPosition = defaultMOPIMotherVolumeXPosition;
    
    G4double defaultMOPIMotherVolumeYPosition = 0.0 *mm;
    MOPIMotherVolumeYPosition = defaultMOPIMotherVolumeYPosition;
    
    G4double defaultMOPIMotherVolumeZPosition = 0.0 *mm;
    MOPIMotherVolumeZPosition = defaultMOPIMotherVolumeZPosition;
    
    // First Kapton Layer of MOPI
    G4double defaultMOPIFirstKaptonLayerXSize = 35 *um;
    MOPIFirstKaptonLayerXSize = defaultMOPIFirstKaptonLayerXSize;
    
    G4double defaultMOPIFirstKaptonLayerYSize = 30 *cm;
    MOPIFirstKaptonLayerYSize = defaultMOPIFirstKaptonLayerYSize;
    
    G4double defaultMOPIFirstKaptonLayerZSize = 30 *cm;
    MOPIFirstKaptonLayerZSize = defaultMOPIFirstKaptonLayerZSize;
    
    G4double defaultMOPIFirstKaptonLayerXPosition = -(MOPIMotherVolumeXSize/2 - (MOPIFirstKaptonLayerXSize/2));
    MOPIFirstKaptonLayerXPosition = defaultMOPIFirstKaptonLayerXPosition;
    
    G4double defaultMOPIFirstKaptonLayerYPosition = 0.0 *mm;
    MOPIFirstKaptonLayerYPosition = defaultMOPIFirstKaptonLayerYPosition;
    
    G4double defaultMOPIFirstKaptonLayerZPosition = 0.0 *mm;
    MOPIFirstKaptonLayerZPosition = defaultMOPIFirstKaptonLayerZPosition;
    
    //First Aluminum  Layer of MOPI
    G4double defaultMOPIFirstAluminumLayerXSize = 15 *um;
    MOPIFirstAluminumLayerXSize = defaultMOPIFirstAluminumLayerXSize;
    
    G4double defaultMOPIFirstAluminumLayerYSize = 30 *cm;
    MOPIFirstAluminumLayerYSize = defaultMOPIFirstAluminumLayerYSize;
    
    G4double defaultMOPIFirstAluminumLayerZSize = 30 *cm;
    MOPIFirstAluminumLayerZSize = defaultMOPIFirstAluminumLayerZSize;
    
    G4double defaultMOPIFirstAluminumLayerXPosition =
    MOPIFirstKaptonLayerXPosition + MOPIFirstKaptonLayerXSize/2 + MOPIFirstAluminumLayerXSize/2;
    MOPIFirstAluminumLayerXPosition = defaultMOPIFirstAluminumLayerXPosition;
    
    G4double defaultMOPIFirstAluminumLayerYPosition = 0.0 *mm;
    MOPIFirstAluminumLayerYPosition = defaultMOPIFirstAluminumLayerYPosition;
    
    G4double defaultMOPIFirstAluminumLayerZPosition = 0.0 *mm;
    MOPIFirstAluminumLayerZPosition = defaultMOPIFirstAluminumLayerZPosition;
    
    // First Air gap of MOPI
    G4double defaultMOPIFirstAirGapXSize = 6000 *um;
    MOPIFirstAirGapXSize = defaultMOPIFirstAirGapXSize;
    
    G4double defaultMOPIFirstAirGapYSize = 30 *cm;
    MOPIFirstAirGapYSize = defaultMOPIFirstAirGapYSize;
    
    G4double defaultMOPIFirstAirGapZSize = 30 *cm;
    MOPIFirstAirGapZSize = defaultMOPIFirstAirGapZSize;
    
    G4double defaultMOPIFirstAirGapXPosition =
    MOPIFirstAluminumLayerXPosition + MOPIFirstAluminumLayerXSize/2 + MOPIFirstAirGapXSize/2;
    MOPIFirstAirGapXPosition = defaultMOPIFirstAirGapXPosition;
    
    G4double defaultMOPIFirstAirGapYPosition = 0.0 *mm;
    MOPIFirstAirGapYPosition = defaultMOPIFirstAirGapYPosition;
    
    G4double defaultMOPIFirstAirGapZPosition = 0.0 *mm;
    MOPIFirstAirGapZPosition = defaultMOPIFirstAirGapZPosition;
    
    // Cathode of MOPI
    G4double defaultMOPICathodeXSize = 25.0 *um;
    MOPICathodeXSize = defaultMOPICathodeXSize;
    
    G4double defaultMOPICathodeYSize = 30.0 *cm;
    MOPICathodeYSize = defaultMOPICathodeYSize;
    
    G4double defaultMOPICathodeZSize = 30.0 *cm;
    MOPICathodeZSize = defaultMOPICathodeZSize;
    
    G4double defaultMOPICathodeXPosition =
    MOPIFirstAirGapXPosition + MOPIFirstAirGapXSize/2 + MOPICathodeXSize/2;
    MOPICathodeXPosition = defaultMOPICathodeXPosition;
    
    G4double defaultMOPICathodeYPosition = 0.0 *mm;
    MOPICathodeYPosition = defaultMOPICathodeYPosition;
    
    G4double defaultMOPICathodeZPosition = 0.0 *mm;
    MOPICathodeZPosition = defaultMOPICathodeZPosition;
    
    // Second Air gap of MOPI
    G4double defaultMOPISecondAirGapXSize = 6000 *um;
    MOPISecondAirGapXSize = defaultMOPISecondAirGapXSize;
    
    G4double defaultMOPISecondAirGapYSize = 30 *cm;
    MOPISecondAirGapYSize = defaultMOPISecondAirGapYSize;
    
    G4double defaultMOPISecondAirGapZSize = 30 *cm;
    MOPISecondAirGapZSize = defaultMOPISecondAirGapZSize;
    
    G4double defaultMOPISecondAirGapXPosition =
    MOPICathodeXPosition + MOPICathodeXSize/2 + MOPISecondAirGapXSize/2;
    MOPISecondAirGapXPosition = defaultMOPISecondAirGapXPosition;
    
    G4double defaultMOPISecondAirGapYPosition = 0.0 *mm;
    MOPISecondAirGapYPosition = defaultMOPISecondAirGapYPosition;
    
    G4double defaultMOPISecondAirGapZPosition = 0.0 *mm;
    MOPISecondAirGapZPosition = defaultMOPISecondAirGapZPosition;
    
    //Second Aluminum  Layer of MOPI
    G4double defaultMOPISecondAluminumLayerXSize = 15 *um;
    MOPISecondAluminumLayerXSize = defaultMOPISecondAluminumLayerXSize;
    
    G4double defaultMOPISecondAluminumLayerYSize = 30 *cm;
    MOPISecondAluminumLayerYSize = defaultMOPISecondAluminumLayerYSize;
    
    G4double defaultMOPISecondAluminumLayerZSize = 30 *cm;
    MOPISecondAluminumLayerZSize = defaultMOPISecondAluminumLayerZSize;
    
    G4double defaultMOPISecondAluminumLayerXPosition =
    MOPISecondAirGapXPosition + MOPISecondAirGapXSize/2 + MOPISecondAluminumLayerXSize/2;
    MOPISecondAluminumLayerXPosition = defaultMOPISecondAluminumLayerXPosition;
    
    G4double defaultMOPISecondAluminumLayerYPosition = 0.0 *mm;
    MOPISecondAluminumLayerYPosition = defaultMOPISecondAluminumLayerYPosition;
    
    G4double defaultMOPISecondAluminumLayerZPosition = 0.0 *mm;
    MOPISecondAluminumLayerZPosition = defaultMOPISecondAluminumLayerZPosition;
    
    // Second Kapton Layer of MOPI
    G4double defaultMOPISecondKaptonLayerXSize = 35 *um;
    MOPISecondKaptonLayerXSize = defaultMOPISecondKaptonLayerXSize;
    
    G4double defaultMOPISecondKaptonLayerYSize = 30 *cm;
    MOPISecondKaptonLayerYSize = defaultMOPISecondKaptonLayerYSize;
    
    G4double defaultMOPISecondKaptonLayerZSize = 30 *cm;
    MOPISecondKaptonLayerZSize = defaultMOPISecondKaptonLayerZSize;
    
    G4double defaultMOPISecondKaptonLayerXPosition =
    MOPISecondAluminumLayerXPosition + MOPISecondAluminumLayerXSize/2 + MOPISecondKaptonLayerXSize/2;
    MOPISecondKaptonLayerXPosition = defaultMOPISecondKaptonLayerXPosition;
    
    G4double defaultMOPISecondKaptonLayerYPosition = 0.0 *mm;
    MOPISecondKaptonLayerYPosition = defaultMOPISecondKaptonLayerYPosition;
    
    G4double defaultMOPISecondKaptonLayerZPosition = 0.0 *mm;
    MOPISecondKaptonLayerZPosition = defaultMOPISecondKaptonLayerZPosition;
    
    
    // FINAL COLLIMATOR: is the collimator giving the final transversal shape
    // of the beam
    G4double defaultinnerRadiusFinalCollimator = 7.5 *mm;
    innerRadiusFinalCollimator = defaultinnerRadiusFinalCollimator;
    
    // DEFAULT DEFINITION OF THE MATERIALS
    // All elements and compound definition follows the NIST database
    
    // ELEMENTS
    G4bool isotopes = false;
    G4Material* aluminumNist = G4NistManager::Instance()->FindOrBuildMaterial("G4_Al", isotopes);
    G4Material* tantalumNist = G4NistManager::Instance()->FindOrBuildMaterial("G4_Ta", isotopes);
    G4Material* copperNistAsMaterial = G4NistManager::Instance()->FindOrBuildMaterial("G4_Cu", isotopes);
    G4Element* zincNist = G4NistManager::Instance()->FindOrBuildElement("Zn");
    G4Element* copperNist = G4NistManager::Instance()->FindOrBuildElement("Cu");
    
    // COMPOUND
    G4Material* airNist =  G4NistManager::Instance()->FindOrBuildMaterial("G4_AIR", isotopes);
    G4Material* kaptonNist = G4NistManager::Instance()->FindOrBuildMaterial("G4_KAPTON", isotopes);
    G4Material* galacticNist = G4NistManager::Instance()->FindOrBuildMaterial("G4_Galactic", isotopes);
    G4Material* PMMANist = G4NistManager::Instance()->FindOrBuildMaterial("G4_PLEXIGLASS", isotopes);
    G4Material* mylarNist = G4NistManager::Instance()->FindOrBuildMaterial("G4_MYLAR", isotopes);
    
    G4double d; // Density
    G4int nComponents;// Number of components
    G4double fractionmass; // Fraction in mass of an element in a material
    
    d = 8.40*g/cm3;
    nComponents = 2;
    G4Material* brass = new G4Material("Brass", d, nComponents);
    brass -> AddElement(zincNist, fractionmass = 30 *perCent);
    brass -> AddElement(copperNist, fractionmass = 70 *perCent);
    
    
    //***************************** PW ***************************************
    
    // DetectorROGeometry Material
    new G4Material("dummyMat", 1., 1.*g/mole, 1.*g/cm3);
    
    //***************************** PW ***************************************
    
    
    
    // MATERIAL ASSIGNMENT
    // Range shifter
    rangeShifterMaterial = airNist;
    
    // Support of the beam line
    beamLineSupportMaterial = aluminumNist;
    
    // Vacuum pipe
    vacuumZoneMaterial = galacticNist;
    
    // Material of the fisrt scattering foil
    firstScatteringFoilMaterial = tantalumNist;
    
    // Material of kapton window
    kaptonWindowMaterial = kaptonNist;
    
    // Material of the stopper
    stopperMaterial = brass;
    
    // Material of the second scattering foil
    secondScatteringFoilMaterial = tantalumNist;
    
    // Materials of the collimators
    firstCollimatorMaterial = PMMANist;
    holeFirstCollimatorMaterial = airNist;
    
    // Box containing the modulator wheel
    modulatorBoxMaterial = aluminumNist;
    holeModulatorBoxMaterial = airNist;
    
    // Materials of the monitor chamber
    layer1MonitorChamberMaterial = kaptonNist;
    layer2MonitorChamberMaterial = copperNistAsMaterial;
    layer3MonitorChamberMaterial = airNist;
    layer4MonitorChamberMaterial = copperNistAsMaterial;
    
    // Mother volume of the MOPI detector
    MOPIMotherVolumeMaterial = airNist;
    MOPIFirstKaptonLayerMaterial = kaptonNist;
    MOPIFirstAluminumLayerMaterial = aluminumNist;
    MOPIFirstAirGapMaterial = airNist;
    MOPICathodeMaterial = mylarNist;
    MOPISecondAirGapMaterial = airNist;
    MOPISecondAluminumLayerMaterial = aluminumNist;
    MOPISecondKaptonLayerMaterial = kaptonNist;
    
    // material of the final nozzle
    nozzleSupportMaterial = PMMANist;
    brassTubeMaterial = brassTube2Material = brassTube3Material = brass;
    holeNozzleSupportMaterial = airNist;
    
    // Material of the final collimator
    finalCollimatorMaterial = brass;
}

void PassiveProtonBeamLine::ConstructPassiveProtonBeamLine()
{
    // -----------------------------
    // Treatment room - World volume
    //------------------------------
    // Treatment room sizes
    const G4double worldX = 400.0 *cm;
    const G4double worldY = 400.0 *cm;
    const G4double worldZ = 400.0 *cm;
    G4bool isotopes = false;
    
    G4Material* airNist =  G4NistManager::Instance()->FindOrBuildMaterial("G4_AIR", isotopes);
    G4Box* treatmentRoom = new G4Box("TreatmentRoom",worldX,worldY,worldZ);
    G4LogicalVolume* logicTreatmentRoom = new G4LogicalVolume(treatmentRoom,
                                                              airNist,
                                                              "logicTreatmentRoom",
                                                              0,0,0);
    physicalTreatmentRoom = new G4PVPlacement(0,
                                              G4ThreeVector(),
                                              "physicalTreatmentRoom",
                                              logicTreatmentRoom,
                                              0,false,0);
    
    
    // The treatment room is invisible in the Visualisation
    logicTreatmentRoom -> SetVisAttributes (G4VisAttributes::Invisible);
    
    // Components of the Passive Proton Beam Line
    HadrontherapyBeamLineSupport();
    HadrontherapyBeamScatteringFoils();
    HadrontherapyRangeShifter();
    HadrontherapyBeamCollimators();
    HadrontherapyBeamMonitoring();
    HadrontherapyMOPIDetector();
    HadrontherapyBeamNozzle();
    HadrontherapyBeamFinalCollimator();
    
    // The following lines construc a typical modulator wheel inside the Passive Beam line.
    // Please remember to set the nodulator material (default is air, i.e. no modulator!)
    // in the HadrontherapyModulator.cc file
    modulator = new HadrontherapyModulator();
    modulator -> BuildModulator(physicalTreatmentRoom);
}

void PassiveProtonBeamLine::HadrontherapyBeamLineSupport()
{
    // ------------------//
    // BEAM LINE SUPPORT //
    //-------------------//
    const G4double beamLineSupportXSize = 1.5*m;
    const G4double beamLineSupportYSize = 20.*mm;
    const G4double beamLineSupportZSize = 600.*mm;
    
    const G4double beamLineSupportXPosition = -1745.09 *mm;
    const G4double beamLineSupportYPosition = -230. *mm;
    const G4double beamLineSupportZPosition = 0.*mm;
    
    G4Box* beamLineSupport = new G4Box("BeamLineSupport",
                                       beamLineSupportXSize,
                                       beamLineSupportYSize,
                                       beamLineSupportZSize);
    
    G4LogicalVolume* logicBeamLineSupport = new G4LogicalVolume(beamLineSupport,
                                                                beamLineSupportMaterial,
                                                                "BeamLineSupport");
    physiBeamLineSupport = new G4PVPlacement(0, G4ThreeVector(beamLineSupportXPosition,
                                                              beamLineSupportYPosition,
                                                              beamLineSupportZPosition),
                                             "BeamLineSupport",
                                             logicBeamLineSupport,
                                             physicalTreatmentRoom, false, 0);
    
    // Visualisation attributes of the beam line support
    
    logicBeamLineSupport -> SetVisAttributes(gray);
    
    //---------------------------------//
    //  Beam line cover 1 (left panel) //
    //---------------------------------//
    const G4double beamLineCoverXSize = 1.5*m;
    const G4double beamLineCoverYSize = 750.*mm;
    const G4double beamLineCoverZSize = 10.*mm;
    
    const G4double beamLineCoverXPosition = -1745.09 *mm;
    const G4double beamLineCoverYPosition = -1000.*mm;
    const G4double beamLineCoverZPosition = 600.*mm;
    
    G4Box* beamLineCover = new G4Box("BeamLineCover",
                                     beamLineCoverXSize,
                                     beamLineCoverYSize,
                                     beamLineCoverZSize);
    
    G4LogicalVolume* logicBeamLineCover = new G4LogicalVolume(beamLineCover,
                                                              beamLineSupportMaterial,
                                                              "BeamLineCover");
    
    physiBeamLineCover = new G4PVPlacement(0, G4ThreeVector(beamLineCoverXPosition,
                                                            beamLineCoverYPosition,
                                                            beamLineCoverZPosition),
                                           "BeamLineCover",
                                           logicBeamLineCover,
                                           physicalTreatmentRoom,
                                           false,
                                           0);
    
    // ---------------------------------//
    //  Beam line cover 2 (rigth panel) //
    // ---------------------------------//
    // It has the same characteristic of beam line cover 1 but set in a different position
    physiBeamLineCover2 = new G4PVPlacement(0, G4ThreeVector(beamLineCoverXPosition,
                                                             beamLineCoverYPosition,
                                                             - beamLineCoverZPosition),
                                            "BeamLineCover2",
                                            logicBeamLineCover,
                                            physicalTreatmentRoom,
                                            false,
                                            0);
    
    logicBeamLineCover -> SetVisAttributes(blue);
}

void PassiveProtonBeamLine::HadrontherapyBeamScatteringFoils()
{
    // ------------//
    // VACUUM PIPE //
    //-------------//
    //
    // First track of the beam line is inside vacuum;
    // The PIPE contains the FIRST SCATTERING FOIL and the KAPTON WINDOW
    G4Box* vacuumZone = new G4Box("VacuumZone", vacuumZoneXSize, vacuumZoneYSize, vacuumZoneZSize);
    G4LogicalVolume* logicVacuumZone = new G4LogicalVolume(vacuumZone, vacuumZoneMaterial, "VacuumZone");
    G4VPhysicalVolume* physiVacuumZone = new G4PVPlacement(0, G4ThreeVector(vacuumZoneXPosition, 0., 0.),
                                                           "VacuumZone", logicVacuumZone, physicalTreatmentRoom, false, 0);
    // --------------------------//
    // THE FIRST SCATTERING FOIL //
    // --------------------------//
    // A thin foil performing a first scattering
    // of the original beam
    firstScatteringFoil = new G4Box("FirstScatteringFoil",
                                    firstScatteringFoilXSize,
                                    firstScatteringFoilYSize,
                                    firstScatteringFoilZSize);
    
    G4LogicalVolume* logicFirstScatteringFoil = new G4LogicalVolume(firstScatteringFoil,
                                                                    firstScatteringFoilMaterial,
                                                                    "FirstScatteringFoil");
    
    physiFirstScatteringFoil = new G4PVPlacement(0, G4ThreeVector(firstScatteringFoilXPosition, 0.,0.),
                                                 "FirstScatteringFoil", logicFirstScatteringFoil, physiVacuumZone,
                                                 false, 0);
    
    logicFirstScatteringFoil -> SetVisAttributes(skyBlue);
    // -------------------//
    // THE KAPTON WINDOWS //
    //--------------------//
    //It prmits the passage of the beam from vacuum to air
    G4Box* solidKaptonWindow = new G4Box("KaptonWindow",
                                         kaptonWindowXSize,
                                         kaptonWindowYSize,
                                         kaptonWindowZSize);
    
    G4LogicalVolume* logicKaptonWindow = new G4LogicalVolume(solidKaptonWindow,
                                                             kaptonWindowMaterial,
                                                             "KaptonWindow");
    
    physiKaptonWindow = new G4PVPlacement(0, G4ThreeVector(kaptonWindowXPosition, 0., 0.),
                                          "KaptonWindow", logicKaptonWindow,
                                          physiVacuumZone, false,       0);
    
    logicKaptonWindow -> SetVisAttributes(darkOrange3);
    
    // ------------//
    // THE STOPPER //
    //-------------//
    // Is a small cylinder able to stop the central component
    // of the beam (having a gaussian shape). It is connected to the SECON SCATTERING FOIL
    // and represent the second element of the scattering system
    G4double phi = 90. *deg;
    // Matrix definition for a 90 deg rotation with respect to Y axis
    G4RotationMatrix rm;
    rm.rotateY(phi);
    
    solidStopper = new G4Tubs("Stopper",
                              innerRadiusStopper,
                              outerRadiusStopper,
                              heightStopper,
                              startAngleStopper,
                              spanningAngleStopper);
    
    logicStopper = new G4LogicalVolume(solidStopper,
                                       stopperMaterial,
                                       "Stopper",
                                       0, 0, 0);
    
    physiStopper = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector(stopperXPosition,
                                                                     stopperYPosition,
                                                                     stopperZPosition)),
                                     "Stopper",
                                     logicStopper,
                                     physicalTreatmentRoom,
                                     false,
                                     0);
    
    logicStopper -> SetVisAttributes(red);
    
    // ---------------------------//
    // THE SECOND SCATTERING FOIL //
    // ---------------------------//
    // It is another thin foil and provides the
    // final diffusion of the beam. It represents the third element of the scattering
    // system;
    
    secondScatteringFoil = new G4Box("SecondScatteringFoil",
                                     secondScatteringFoilXSize,
                                     secondScatteringFoilYSize,
                                     secondScatteringFoilZSize);
    
    G4LogicalVolume* logicSecondScatteringFoil = new G4LogicalVolume(secondScatteringFoil,
                                                                     secondScatteringFoilMaterial,
                                                                     "SecondScatteringFoil");
    
    physiSecondScatteringFoil = new G4PVPlacement(0, G4ThreeVector(secondScatteringFoilXPosition,
                                                                   secondScatteringFoilYPosition,
                                                                   secondScatteringFoilZPosition),
                                                  "SeconScatteringFoil",
                                                  logicSecondScatteringFoil,
                                                  physicalTreatmentRoom,
                                                  false,
                                                  0);
    
    logicSecondScatteringFoil -> SetVisAttributes(skyBlue);
}
void PassiveProtonBeamLine::HadrontherapyRangeShifter()
{
    // ---------------------------- //
    //         THE RANGE SHIFTER    //
    // -----------------------------//
    // It is a slab of PMMA acting as energy degreader of
    // primary beam
    solidRangeShifterBox = new G4Box("RangeShifterBox",
                                     rangeShifterXSize,
                                     rangeShifterYSize,
                                     rangeShifterZSize);
    
    logicRangeShifterBox = new G4LogicalVolume(solidRangeShifterBox,
                                               rangeShifterMaterial,
                                               "RangeShifterBox");
    physiRangeShifterBox = new G4PVPlacement(0,
                                             G4ThreeVector(rangeShifterXPosition, 0., 0.),
                                             "RangeShifterBox",
                                             logicRangeShifterBox,
                                             physicalTreatmentRoom,
                                             false,
                                             0);
    
    
    logicRangeShifterBox -> SetVisAttributes(yellow);
}
void PassiveProtonBeamLine::HadrontherapyBeamCollimators()
{
    // -----------------//
    // FIRST COLLIMATOR //
    // -----------------//
    // It is a slab of PMMA with an hole in its center
    const G4double firstCollimatorXSize = 20.*mm;
    const G4double firstCollimatorYSize = 100.*mm;
    const G4double firstCollimatorZSize = 100.*mm;
    
    const G4double firstCollimatorXPosition = -2673.00*mm;
    const G4double firstCollimatorYPosition = 0.*mm;
    const G4double firstCollimatorZPosition = 0.*mm;
    
    G4Box* solidFirstCollimator = new G4Box("FirstCollimator",
                                            firstCollimatorXSize,
                                            firstCollimatorYSize,
                                            firstCollimatorZSize);
    
    G4LogicalVolume* logicFirstCollimator = new G4LogicalVolume(solidFirstCollimator,
                                                                firstCollimatorMaterial,
                                                                "FirstCollimator");
    
    physiFirstCollimator = new G4PVPlacement(0, G4ThreeVector(firstCollimatorXPosition,
                                                              firstCollimatorYPosition,
                                                              firstCollimatorZPosition),
                                             "FirstCollimator",
                                             logicFirstCollimator,
                                             physicalTreatmentRoom,
                                             false,
                                             0);
    // ----------------------------//
    // Hole of the first collimator//
    //-----------------------------//
    G4double innerRadiusHoleFirstCollimator   = 0.*mm;
    G4double outerRadiusHoleFirstCollimator   = 15.*mm;
    G4double hightHoleFirstCollimator         = 20.*mm;
    G4double startAngleHoleFirstCollimator    = 0.*deg;
    G4double spanningAngleHoleFirstCollimator = 360.*deg;
    
    G4Tubs* solidHoleFirstCollimator = new G4Tubs("HoleFirstCollimator",
                                                  innerRadiusHoleFirstCollimator,
                                                  outerRadiusHoleFirstCollimator,
                                                  hightHoleFirstCollimator,
                                                  startAngleHoleFirstCollimator,
                                                  spanningAngleHoleFirstCollimator);
    
    G4LogicalVolume* logicHoleFirstCollimator = new G4LogicalVolume(solidHoleFirstCollimator,
                                                                    holeFirstCollimatorMaterial,
                                                                    "HoleFirstCollimator",
                                                                    0, 0, 0);
    G4double phi = 90. *deg;
    // Matrix definition for a 90 deg rotation. Also used for other volumes
    G4RotationMatrix rm;
    rm.rotateY(phi);
    
    physiHoleFirstCollimator = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector()),
                                                 "HoleFirstCollimator",
                                                 logicHoleFirstCollimator,
                                                 physiFirstCollimator,
                                                 false,
                                                 0);
    // ------------------//
    // SECOND COLLIMATOR //
    //-------------------//
    // It is a slab of PMMA with an hole in its center
    const G4double secondCollimatorXPosition = -1900.00*mm;
    const G4double secondCollimatorYPosition =  0*mm;
    const G4double secondCollimatorZPosition =  0*mm;
    
    physiSecondCollimator = new G4PVPlacement(0, G4ThreeVector(secondCollimatorXPosition,
                                                               secondCollimatorYPosition,
                                                               secondCollimatorZPosition),
                                              "SecondCollimator",
                                              logicFirstCollimator,
                                              physicalTreatmentRoom,
                                              false,
                                              0);
    
    // ------------------------------//
    // Hole of the second collimator //
    // ------------------------------//
    physiHoleSecondCollimator = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector()),
                                                  "HoleSecondCollimator",
                                                  logicHoleFirstCollimator,
                                                  physiSecondCollimator,
                                                  false,
                                                  0);
    
    // --------------------------------------//
    // FIRST SIDE OF THE MODULATOR BOX      //
    // --------------------------------------//
    // The modulator box is an aluminum box in which
    // the range shifter and the energy modulator are located
    // In this example only the entrance and exit
    // faces of the box are simulated.
    // Each face is an aluminum slab with an hole in its center
    
    const G4double firstCollimatorModulatorXSize = 10.*mm;
    const G4double firstCollimatorModulatorYSize = 200.*mm;
    const G4double firstCollimatorModulatorZSize = 200.*mm;
    
    const G4double firstCollimatorModulatorXPosition = -2523.00*mm;
    const G4double firstCollimatorModulatorYPosition = 0.*mm;
    const G4double firstCollimatorModulatorZPosition = 0.*mm;
    
    G4Box* solidFirstCollimatorModulatorBox = new G4Box("FirstCollimatorModulatorBox",
                                                        firstCollimatorModulatorXSize,
                                                        firstCollimatorModulatorYSize,
                                                        firstCollimatorModulatorZSize);
    
    G4LogicalVolume* logicFirstCollimatorModulatorBox = new G4LogicalVolume(solidFirstCollimatorModulatorBox,
                                                                            modulatorBoxMaterial,
                                                                            "FirstCollimatorModulatorBox");
    
    physiFirstCollimatorModulatorBox = new G4PVPlacement(0, G4ThreeVector(firstCollimatorModulatorXPosition,
                                                                          firstCollimatorModulatorYPosition,
                                                                          firstCollimatorModulatorZPosition),
                                                         "FirstCollimatorModulatorBox",
                                                         logicFirstCollimatorModulatorBox,
                                                         physicalTreatmentRoom, false, 0);
    
    // ----------------------------------------------------//
    //   Hole of the first collimator of the modulator box //
    // ----------------------------------------------------//
    const G4double innerRadiusHoleFirstCollimatorModulatorBox = 0.*mm;
    const G4double outerRadiusHoleFirstCollimatorModulatorBox = 31.*mm;
    const G4double hightHoleFirstCollimatorModulatorBox = 10.*mm;
    const G4double startAngleHoleFirstCollimatorModulatorBox = 0.*deg;
    const G4double spanningAngleHoleFirstCollimatorModulatorBox = 360.*deg;
    
    G4Tubs* solidHoleFirstCollimatorModulatorBox  = new G4Tubs("HoleFirstCollimatorModulatorBox",
                                                               innerRadiusHoleFirstCollimatorModulatorBox,
                                                               outerRadiusHoleFirstCollimatorModulatorBox,
                                                               hightHoleFirstCollimatorModulatorBox ,
                                                               startAngleHoleFirstCollimatorModulatorBox,
                                                               spanningAngleHoleFirstCollimatorModulatorBox);
    
    G4LogicalVolume* logicHoleFirstCollimatorModulatorBox = new G4LogicalVolume(solidHoleFirstCollimatorModulatorBox,
                                                                                holeModulatorBoxMaterial,
                                                                                "HoleFirstCollimatorModulatorBox",
                                                                                0, 0, 0);
    
    physiHoleFirstCollimatorModulatorBox = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector()),
                                                             "HoleFirstCollimatorModulatorBox",
                                                             logicHoleFirstCollimatorModulatorBox,
                                                             physiFirstCollimatorModulatorBox, false, 0);
    
    // --------------------------------------------------//
    //       SECOND SIDE OF THE MODULATOR BOX            //
    // --------------------------------------------------//
    const G4double secondCollimatorModulatorXSize = 10.*mm;
    const G4double secondCollimatorModulatorYSize = 200.*mm;
    const G4double secondCollimatorModulatorZSize = 200.*mm;
    
    const G4double secondCollimatorModulatorXPosition = -1953.00 *mm;
    
    const G4double secondCollimatorModulatorYPosition = 0.*mm;
    const G4double secondCollimatorModulatorZPosition = 0.*mm;
    
    G4Box* solidSecondCollimatorModulatorBox = new G4Box("SecondCollimatorModulatorBox",
                                                         secondCollimatorModulatorXSize,
                                                         secondCollimatorModulatorYSize,
                                                         secondCollimatorModulatorZSize);
    
    G4LogicalVolume* logicSecondCollimatorModulatorBox = new G4LogicalVolume(solidSecondCollimatorModulatorBox,
                                                                             modulatorBoxMaterial,
                                                                             "SecondCollimatorModulatorBox");
    
    physiSecondCollimatorModulatorBox = new G4PVPlacement(0, G4ThreeVector(secondCollimatorModulatorXPosition,
                                                                           secondCollimatorModulatorYPosition,
                                                                           secondCollimatorModulatorZPosition),
                                                          "SecondCollimatorModulatorBox",
                                                          logicSecondCollimatorModulatorBox,
                                                          physicalTreatmentRoom, false, 0);
    
    // ----------------------------------------------//
    //   Hole of the second collimator modulator box //
    // ----------------------------------------------//
    const G4double innerRadiusHoleSecondCollimatorModulatorBox = 0.*mm;
    const G4double outerRadiusHoleSecondCollimatorModulatorBox = 31.*mm;
    const G4double hightHoleSecondCollimatorModulatorBox = 10.*mm;
    const G4double startAngleHoleSecondCollimatorModulatorBox = 0.*deg;
    const G4double spanningAngleHoleSecondCollimatorModulatorBox = 360.*deg;
    
    G4Tubs* solidHoleSecondCollimatorModulatorBox  = new G4Tubs("HoleSecondCollimatorModulatorBox",
                                                                innerRadiusHoleSecondCollimatorModulatorBox,
                                                                outerRadiusHoleSecondCollimatorModulatorBox,
                                                                hightHoleSecondCollimatorModulatorBox ,
                                                                startAngleHoleSecondCollimatorModulatorBox,
                                                                spanningAngleHoleSecondCollimatorModulatorBox);
    
    G4LogicalVolume* logicHoleSecondCollimatorModulatorBox = new G4LogicalVolume(solidHoleSecondCollimatorModulatorBox,
                                                                                 holeModulatorBoxMaterial,
                                                                                 "HoleSecondCollimatorModulatorBox",
                                                                                 0, 0, 0);
    
    physiHoleSecondCollimatorModulatorBox = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector()),
                                                              "HoleSecondCollimatorModulatorBox",
                                                              logicHoleSecondCollimatorModulatorBox,
                                                              physiSecondCollimatorModulatorBox, false, 0);
    
    logicFirstCollimator -> SetVisAttributes(yellow);
    logicFirstCollimatorModulatorBox -> SetVisAttributes(blue);
    logicSecondCollimatorModulatorBox -> SetVisAttributes(blue);
}

void PassiveProtonBeamLine::HadrontherapyBeamMonitoring()
{
    // ----------------------------
    //   THE FIRST MONITOR CHAMBER
    // ----------------------------
    // A monitor chamber is a free-air  ionisation chamber
    // able to measure do proton fluence during the treatment.
    // Here its responce is not simulated in terms of produced
    // charge but only the energy losses are taked into account.
    // Each chamber consist of 9 mm of air in a box
    // that has two layers one of kapton and one
    // of copper
    const G4double monitor1XSize = 4.525022*mm;
    const G4double monitor2XSize = 0.000011*mm;
    const G4double monitor3XSize = 4.5*mm;
    const G4double monitorYSize = 10.*cm;
    const G4double monitorZSize = 10.*cm;
    const G4double monitor1XPosition = -1262.47498 *mm;
    const G4double monitor2XPosition = -4.500011*mm;
    const G4double monitor4XPosition = 4.500011*mm;
    
    G4Box* solidFirstMonitorLayer1 = new G4Box("FirstMonitorLayer1",
                                               monitor1XSize,
                                               monitorYSize,
                                               monitorZSize);
    
    G4LogicalVolume* logicFirstMonitorLayer1 = new G4LogicalVolume(solidFirstMonitorLayer1,
                                                                   layer1MonitorChamberMaterial,
                                                                   "FirstMonitorLayer1");
    
    physiFirstMonitorLayer1 = new G4PVPlacement(0,
                                                G4ThreeVector(monitor1XPosition,0.*cm,0.*cm),
                                                "FirstMonitorLayer1",
                                                logicFirstMonitorLayer1,
                                                physicalTreatmentRoom,
                                                false,
                                                0);
    
    G4Box* solidFirstMonitorLayer2 = new G4Box("FirstMonitorLayer2",
                                               monitor2XSize,
                                               monitorYSize,
                                               monitorZSize);
    
    G4LogicalVolume* logicFirstMonitorLayer2 = new G4LogicalVolume(solidFirstMonitorLayer2,
                                                                   layer2MonitorChamberMaterial,
                                                                   "FirstMonitorLayer2");
    
    physiFirstMonitorLayer2 = new G4PVPlacement(0, G4ThreeVector(monitor2XPosition,0.*cm,0.*cm),
                                                "FirstMonitorLayer2",
                                                logicFirstMonitorLayer2,
                                                physiFirstMonitorLayer1,
                                                false,
                                                0);
    
    G4Box* solidFirstMonitorLayer3 = new G4Box("FirstMonitorLayer3",
                                               monitor3XSize,
                                               monitorYSize,
                                               monitorZSize);
    
    G4LogicalVolume* logicFirstMonitorLayer3 = new G4LogicalVolume(solidFirstMonitorLayer3,
                                                                   layer3MonitorChamberMaterial,
                                                                   "FirstMonitorLayer3");
    
    physiFirstMonitorLayer3 = new G4PVPlacement(0,
                                                G4ThreeVector(0.*mm,0.*cm,0.*cm),
                                                "MonitorLayer3",
                                                logicFirstMonitorLayer3,
                                                physiFirstMonitorLayer1,
                                                false,
                                                0);
    
    G4Box* solidFirstMonitorLayer4 = new G4Box("FirstMonitorLayer4",
                                               monitor2XSize,
                                               monitorYSize,
                                               monitorZSize);
    
    G4LogicalVolume* logicFirstMonitorLayer4 = new G4LogicalVolume(solidFirstMonitorLayer4,
                                                                   layer4MonitorChamberMaterial,
                                                                   "FirstMonitorLayer4");
    
    physiFirstMonitorLayer4 = new G4PVPlacement(0, G4ThreeVector(monitor4XPosition,0.*cm,0.*cm),
                                                "FirstMonitorLayer4",
                                                logicFirstMonitorLayer4,
                                                physiFirstMonitorLayer1, false, 0);
    // ----------------------------//
    // THE SECOND MONITOR CHAMBER  //
    // ----------------------------//
    physiSecondMonitorLayer1 = new G4PVPlacement(0, G4ThreeVector(-1131.42493 *mm,0.*cm,0.*cm),
                                                 "SecondMonitorLayer1", logicFirstMonitorLayer1,physicalTreatmentRoom, false, 0);
    
    physiSecondMonitorLayer2 = new G4PVPlacement(0, G4ThreeVector( monitor2XPosition,0.*cm,0.*cm), "SecondMonitorLayer2",
                                                 logicFirstMonitorLayer2, physiSecondMonitorLayer1, false, 0);
    
    physiSecondMonitorLayer3 = new G4PVPlacement(0, G4ThreeVector(0.*mm,0.*cm,0.*cm), "MonitorLayer3",
                                                 logicFirstMonitorLayer3, physiSecondMonitorLayer1, false, 0);
    
    physiSecondMonitorLayer4 = new G4PVPlacement(0, G4ThreeVector(monitor4XPosition,0.*cm,0.*cm), "SecondMonitorLayer4",
                                                 logicFirstMonitorLayer4, physiSecondMonitorLayer1, false, 0);
    
    logicFirstMonitorLayer3 -> SetVisAttributes(white);
    
}
void PassiveProtonBeamLine::HadrontherapyMOPIDetector()
{
    // --------------------------------//
    //        THE MOPI DETECTOR        //
    // --------------------------------//
    // MOPI DETECTOR: two orthogonal microstrip gas detectors developed
    // by the INFN Section of Turin in collaboration with some
    // of the author of this example. It permits the
    // on-line check of the beam simmetry via the signal
    // integration of the collected charge for each strip.
    //
    // In this example it is simulated as:
    // 1. First anode: 35 mu of kapton + 15 mu of aluminum,
    // 2. First air gap: 6 mm of air,
    // 3. The cathode: 1 mu Al + 25 mu mylar + 1 mu Al
    //    (in common with the two air gap),
    // 4. Second air gap: 6 mm of air,
    // 5  Second anode: 15 mu Al + 35 mu kapton
    // Color used in the graphical output
    
    
    // Mother volume
    solidMOPIMotherVolume = new G4Box("MOPIMotherVolume",
                                      MOPIMotherVolumeXSize/2,
                                      MOPIMotherVolumeYSize/2,
                                      MOPIMotherVolumeYSize/2);
    
    logicMOPIMotherVolume = new G4LogicalVolume(solidMOPIMotherVolume,
                                                MOPIMotherVolumeMaterial,
                                                "MOPIMotherVolume");
    physiMOPIMotherVolume = new G4PVPlacement(0,
                                              G4ThreeVector(MOPIMotherVolumeXPosition,
                                                            MOPIMotherVolumeYPosition,
                                                            MOPIMotherVolumeZPosition),
                                              "MOPIMotherVolume",
                                              logicMOPIMotherVolume,
                                              physicalTreatmentRoom,
                                              false,
                                              0);
    
    // First Kapton layer
    solidMOPIFirstKaptonLayer = new G4Box("MOPIFirstKaptonLayer",
                                          MOPIFirstKaptonLayerXSize/2,
                                          MOPIFirstKaptonLayerYSize/2 ,
                                          MOPIFirstKaptonLayerZSize/2);
    
    logicMOPIFirstKaptonLayer = new G4LogicalVolume(solidMOPIFirstKaptonLayer,
                                                    MOPIFirstKaptonLayerMaterial,
                                                    "MOPIFirstKaptonLayer");
    
    physiMOPIFirstKaptonLayer = new G4PVPlacement(0,
                                                  G4ThreeVector(MOPIFirstKaptonLayerXPosition,
                                                                MOPIFirstKaptonLayerYPosition ,
                                                                MOPIFirstKaptonLayerZPosition),
                                                  "MOPIFirstKaptonLayer",
                                                  logicMOPIFirstKaptonLayer,
                                                  physiMOPIMotherVolume,
                                                  false,
                                                  0);
    
    // First Aluminum layer
    solidMOPIFirstAluminumLayer = new G4Box("MOPIFirstAluminumLayer",
                                            MOPIFirstAluminumLayerXSize/2,
                                            MOPIFirstAluminumLayerYSize/2 ,
                                            MOPIFirstAluminumLayerZSize/2);
    
    logicMOPIFirstAluminumLayer = new G4LogicalVolume(solidMOPIFirstAluminumLayer,
                                                      MOPIFirstAluminumLayerMaterial,
                                                      "MOPIFirstAluminumLayer");
    
    physiMOPIFirstAluminumLayer = new G4PVPlacement(0,
                                                    G4ThreeVector(MOPIFirstAluminumLayerXPosition,
                                                                  MOPIFirstAluminumLayerYPosition ,
                                                                  MOPIFirstAluminumLayerZPosition),
                                                    "MOPIFirstAluminumLayer",
                                                    logicMOPIFirstAluminumLayer, physiMOPIMotherVolume, false, 0);
    
    // First Air GAP
    solidMOPIFirstAirGap = new G4Box("MOPIFirstAirGap",
                                     MOPIFirstAirGapXSize/2,
                                     MOPIFirstAirGapYSize/2,
                                     MOPIFirstAirGapZSize/2);
    
    logicMOPIFirstAirGap = new G4LogicalVolume(solidMOPIFirstAirGap,
                                               MOPIFirstAirGapMaterial,
                                               "MOPIFirstAirgap");
    
    physiMOPIFirstAirGap = new G4PVPlacement(0,
                                             G4ThreeVector(MOPIFirstAirGapXPosition,
                                                           MOPIFirstAirGapYPosition ,
                                                           MOPIFirstAirGapZPosition),
                                             "MOPIFirstAirGap",
                                             logicMOPIFirstAirGap, physiMOPIMotherVolume, false, 0);
    
    
    // The Cathode
    solidMOPICathode = new G4Box("MOPICathode",
                                 MOPICathodeXSize/2,
                                 MOPICathodeYSize/2,
                                 MOPICathodeZSize/2);
    
    logicMOPICathode = new G4LogicalVolume(solidMOPICathode,
                                           MOPICathodeMaterial,
                                           "MOPICathode");
    
    physiMOPICathode = new G4PVPlacement(0,
                                         G4ThreeVector(MOPICathodeXPosition,
                                                       MOPICathodeYPosition ,
                                                       MOPICathodeZPosition),
                                         "MOPICathode",
                                         logicMOPICathode,
                                         physiMOPIMotherVolume, false, 0);
    
    // Second Air GAP
    solidMOPISecondAirGap = new G4Box("MOPISecondAirGap",
                                      MOPISecondAirGapXSize/2,
                                      MOPISecondAirGapYSize/2,
                                      MOPISecondAirGapZSize/2);
    
    logicMOPISecondAirGap = new G4LogicalVolume(solidMOPISecondAirGap,
                                                MOPISecondAirGapMaterial,
                                                "MOPISecondAirgap");
    
    physiMOPISecondAirGap = new G4PVPlacement(0,
                                              G4ThreeVector(MOPISecondAirGapXPosition,
                                                            MOPISecondAirGapYPosition ,
                                                            MOPISecondAirGapZPosition),
                                              "MOPISecondAirGap",
                                              logicMOPISecondAirGap, physiMOPIMotherVolume, false, 0);
    
    // Second Aluminum layer
    solidMOPISecondAluminumLayer = new G4Box("MOPISecondAluminumLayer",
                                             MOPISecondAluminumLayerXSize/2,
                                             MOPISecondAluminumLayerYSize/2 ,
                                             MOPISecondAluminumLayerZSize/2);
    
    logicMOPISecondAluminumLayer = new G4LogicalVolume(solidMOPISecondAluminumLayer,
                                                       MOPISecondAluminumLayerMaterial,
                                                       "MOPISecondAluminumLayer");
    
    physiMOPISecondAluminumLayer = new G4PVPlacement(0,
                                                     G4ThreeVector(MOPISecondAluminumLayerXPosition,
                                                                   MOPISecondAluminumLayerYPosition ,
                                                                   MOPISecondAluminumLayerZPosition),
                                                     "MOPISecondAluminumLayer",
                                                     logicMOPISecondAluminumLayer,
                                                     physiMOPIMotherVolume,
                                                     false,
                                                     0);
    
    // Second Kapton layer
    solidMOPISecondKaptonLayer = new G4Box("MOPISecondKaptonLayer",
                                           MOPISecondKaptonLayerXSize/2,
                                           MOPISecondKaptonLayerYSize/2 ,
                                           MOPISecondKaptonLayerZSize/2);
    
    logicMOPISecondKaptonLayer = new G4LogicalVolume(solidMOPISecondKaptonLayer,
                                                     MOPIFirstKaptonLayerMaterial,
                                                     "MOPISecondKaptonLayer");
    
    physiMOPISecondKaptonLayer = new G4PVPlacement(0,
                                                   G4ThreeVector(MOPISecondKaptonLayerXPosition,
                                                                 MOPISecondKaptonLayerYPosition ,
                                                                 MOPISecondKaptonLayerZPosition),
                                                   "MOPISecondKaptonLayer",
                                                   logicMOPISecondKaptonLayer,
                                                   physiMOPIMotherVolume,
                                                   false,
                                                   0);
    
    logicMOPIFirstAirGap -> SetVisAttributes(darkGreen);
    logicMOPISecondAirGap -> SetVisAttributes(darkGreen);
    
}
void PassiveProtonBeamLine::HadrontherapyBeamNozzle()
{
    // ------------------------------//
    // THE FINAL TUBE AND COLLIMATOR //
    //-------------------------------//
    // The last part of the transport beam line consists of
    // a 59 mm thick PMMA slab (to stop all the diffused radiation), a 370 mm brass tube
    // (to well collimate the proton beam) and a final collimator with 25 mm diameter
    // aperture (that provide the final trasversal shape of the beam)
    
    // -------------------//
    //     PMMA SUPPORT   //
    // -------------------//
    const G4double nozzleSupportXSize = 29.5 *mm;
    const G4double nozzleSupportYSize = 180. *mm;
    const G4double nozzleSupportZSize = 180. *mm;
    
    const G4double nozzleSupportXPosition = -397.50 *mm;
    
    G4double phi = 90. *deg;
    // Matrix definition for a 90 deg rotation. Also used for other volumes
    G4RotationMatrix rm;
    rm.rotateY(phi);
    
    G4Box* solidNozzleSupport = new G4Box("NozzleSupport",
                                          nozzleSupportXSize,
                                          nozzleSupportYSize,
                                          nozzleSupportZSize);
    
    G4LogicalVolume* logicNozzleSupport = new G4LogicalVolume(solidNozzleSupport,
                                                              nozzleSupportMaterial,
                                                              "NozzleSupport");
    
    physiNozzleSupport = new G4PVPlacement(0, G4ThreeVector(nozzleSupportXPosition,0., 0.),
                                           "NozzleSupport",
                                           logicNozzleSupport,
                                           physicalTreatmentRoom,
                                           false,
                                           0);
    
    logicNozzleSupport -> SetVisAttributes(yellow);
    
    
    
    //------------------------------------//
    // HOLE IN THE SUPPORT                //
    //------------------------------------//
    const G4double innerRadiusHoleNozzleSupport = 0.*mm;
    const G4double outerRadiusHoleNozzleSupport = 21.5*mm;
    const G4double hightHoleNozzleSupport = 29.5 *mm;
    const G4double startAngleHoleNozzleSupport = 0.*deg;
    const G4double spanningAngleHoleNozzleSupport = 360.*deg;
    
    G4Tubs* solidHoleNozzleSupport = new G4Tubs("HoleNozzleSupport",
                                                innerRadiusHoleNozzleSupport,
                                                outerRadiusHoleNozzleSupport,
                                                hightHoleNozzleSupport,
                                                startAngleHoleNozzleSupport,
                                                spanningAngleHoleNozzleSupport);
    
    G4LogicalVolume* logicHoleNozzleSupport = new G4LogicalVolume(solidHoleNozzleSupport,
                                                                  holeNozzleSupportMaterial,
                                                                  "HoleNozzleSupport",
                                                                  0,
                                                                  0,
                                                                  0);
    
    
    physiHoleNozzleSupport = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector()),
                                               "HoleNozzleSupport",
                                               logicHoleNozzleSupport,
                                               physiNozzleSupport,
                                               false, 0);
    
    logicHoleNozzleSupport -> SetVisAttributes(darkOrange3);
    
    // ---------------------------------//
    //     BRASS TUBE 1 (phantom side)    //
    // ---------------------------------//
    const G4double innerRadiusBrassTube= 18.*mm;
    const G4double outerRadiusBrassTube = 21.5 *mm;
    const G4double hightBrassTube = 140.5*mm;
    const G4double startAngleBrassTube = 0.*deg;
    const G4double spanningAngleBrassTube = 360.*deg;
    
    const G4double brassTubeXPosition = -227.5 *mm;
    
    G4Tubs* solidBrassTube = new G4Tubs("BrassTube",
                                        innerRadiusBrassTube,
                                        outerRadiusBrassTube,
                                        hightBrassTube,
                                        startAngleBrassTube,
                                        spanningAngleBrassTube);
    
    G4LogicalVolume* logicBrassTube = new G4LogicalVolume(solidBrassTube,
                                                          brassTubeMaterial,
                                                          "BrassTube",
                                                          0, 0, 0);
    
    physiBrassTube = new G4PVPlacement(G4Transform3D(rm,
                                                     G4ThreeVector(brassTubeXPosition,
                                                                   0.,
                                                                   0.)),
                                       "BrassTube",
                                       logicBrassTube,
                                       physicalTreatmentRoom,
                                       false,
                                       0);
    
    logicBrassTube -> SetVisAttributes(darkOrange3);
    
    // ----------------------------------------------//
    //     BRASS TUBE 2 (inside the PMMA support)    //
    // ----------------------------------------------//
    const G4double innerRadiusBrassTube2= 18.*mm;
    const G4double outerRadiusBrassTube2 = 21.5 *mm;
    const G4double hightBrassTube2 = 29.5*mm;
    const G4double startAngleBrassTube2 = 0.*deg;
    const G4double spanningAngleBrassTube2 = 360.*deg;
    
    //  const G4double brassTube2XPosition = -227.5 *mm;
    
    G4Tubs* solidBrassTube2 = new G4Tubs("BrassTube2",
                                         innerRadiusBrassTube2,
                                         outerRadiusBrassTube2,
                                         hightBrassTube2,
                                         startAngleBrassTube2,
                                         spanningAngleBrassTube2);
    
    G4LogicalVolume* logicBrassTube2 = new G4LogicalVolume(solidBrassTube2,
                                                           brassTube2Material,
                                                           "BrassTube2",
                                                           0, 0, 0);
    
    physiBrassTube2 = new G4PVPlacement(G4Transform3D(rm,
                                                      G4ThreeVector(0,
                                                                    0.,
                                                                    0.)),
                                        "BrassTube2",
                                        logicBrassTube2,
                                        physiNozzleSupport,
                                        false,
                                        0);
    
    logicBrassTube2 -> SetVisAttributes(darkOrange3);
    
    
    // --------------------------------------//
    //     BRASS TUBE 3 (beam line side)    //
    // -------------------------------------//
    const G4double innerRadiusBrassTube3= 18.*mm;
    const G4double outerRadiusBrassTube3 = 21.5 *mm;
    const G4double hightBrassTube3 = 10.0 *mm;
    const G4double startAngleBrassTube3 = 0.*deg;
    const G4double spanningAngleBrassTube3 = 360.*deg;
    
    const G4double brassTube3XPosition = -437 *mm;
    
    G4Tubs* solidBrassTube3 = new G4Tubs("BrassTube3",
                                         innerRadiusBrassTube3,
                                         outerRadiusBrassTube3,
                                         hightBrassTube3,
                                         startAngleBrassTube3,
                                         spanningAngleBrassTube3);
    
    G4LogicalVolume* logicBrassTube3 = new G4LogicalVolume(solidBrassTube3,
                                                           brassTube3Material,
                                                           "BrassTube3",
                                                           0, 0, 0);
    
    physiBrassTube3 = new G4PVPlacement(G4Transform3D(rm,
                                                      G4ThreeVector(brassTube3XPosition,
                                                                    0.,
                                                                    0.)),
                                        "BrassTube3",
                                        logicBrassTube3,
                                        physicalTreatmentRoom,
                                        false,
                                        0);
    
    logicBrassTube3 -> SetVisAttributes(darkOrange3);
}

void PassiveProtonBeamLine::HadrontherapyBeamFinalCollimator()
{
    // -----------------------//
    //     FINAL COLLIMATOR   //
    //------------------------//
    const G4double outerRadiusFinalCollimator = 21.5*mm;
    const G4double hightFinalCollimator = 3.5*mm;
    const G4double startAngleFinalCollimator = 0.*deg;
    const G4double spanningAngleFinalCollimator = 360.*deg;
    const G4double finalCollimatorXPosition = -83.5 *mm;
    
    G4double phi = 90. *deg;
    
    // Matrix definition for a 90 deg rotation. Also used for other volumes
    G4RotationMatrix rm;
    rm.rotateY(phi);
    
    solidFinalCollimator = new G4Tubs("FinalCollimator",
                                      innerRadiusFinalCollimator,
                                      outerRadiusFinalCollimator,
                                      hightFinalCollimator,
                                      startAngleFinalCollimator,
                                      spanningAngleFinalCollimator);
    
    G4LogicalVolume* logicFinalCollimator = new G4LogicalVolume(solidFinalCollimator,
                                                                finalCollimatorMaterial,
                                                                "FinalCollimator",
                                                                0,
                                                                0,
                                                                0);
    
    physiFinalCollimator = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector(finalCollimatorXPosition,0.,0.)),
                                             "FinalCollimator", logicFinalCollimator, physicalTreatmentRoom, false, 0);
    
    logicFinalCollimator -> SetVisAttributes(yellow);
}
void PassiveProtonBeamLine::SetRangeShifterXPosition(G4double value)
{
    physiRangeShifterBox -> SetTranslation(G4ThreeVector(value, 0., 0.));
    G4RunManager::GetRunManager() -> GeometryHasBeenModified();
    G4cout << "The Range Shifter is translated to"<< value/mm <<"mm along the X axis" <<G4endl;
}

void PassiveProtonBeamLine::SetRangeShifterXSize(G4double value)
{
    solidRangeShifterBox -> SetXHalfLength(value) ;
    G4cout << "RangeShifter size X (mm): "<< ((solidRangeShifterBox -> GetXHalfLength())*2.)/mm
    << G4endl;
    G4RunManager::GetRunManager() -> GeometryHasBeenModified();
}

void PassiveProtonBeamLine::SetFirstScatteringFoilXSize(G4double value)
{
    firstScatteringFoil -> SetXHalfLength(value);
    G4RunManager::GetRunManager() -> GeometryHasBeenModified();
    G4cout <<"The X size of the first scattering foil is (mm):"<<
    ((firstScatteringFoil -> GetXHalfLength())*2.)/mm
    << G4endl;
}

void PassiveProtonBeamLine::SetSecondScatteringFoilXSize(G4double value)
{
    secondScatteringFoil -> SetXHalfLength(value);
    G4RunManager::GetRunManager() -> GeometryHasBeenModified();
    G4cout <<"The X size of the second scattering foil is (mm):"<<
    ((secondScatteringFoil -> GetXHalfLength())*2.)/mm
    << G4endl;
}

void PassiveProtonBeamLine::SetOuterRadiusStopper(G4double value)
{
    solidStopper -> SetOuterRadius(value);
    G4RunManager::GetRunManager() -> GeometryHasBeenModified();
    G4cout << "OuterRadius od the Stopper is (mm):"
    << solidStopper -> GetOuterRadius()/mm
    << G4endl;
}

void PassiveProtonBeamLine::SetInnerRadiusFinalCollimator(G4double value)
{
    solidFinalCollimator -> SetInnerRadius(value);
    G4RunManager::GetRunManager() -> GeometryHasBeenModified();
    G4cout<<"Inner Radius of the final collimator is (mm):"
        << solidFinalCollimator -> GetInnerRadius()/mm
        << G4endl;
}

void PassiveProtonBeamLine::SetRSMaterial(G4String materialChoice)
{
    if (G4Material* pttoMaterial = G4NistManager::Instance()->FindOrBuildMaterial(materialChoice, false) )
    {
        if (pttoMaterial)
        {
            rangeShifterMaterial  = pttoMaterial;
            logicRangeShifterBox -> SetMaterial(pttoMaterial);
            G4cout << "The material of the Range Shifter has been changed to " << materialChoice << G4endl;
        }
    }
    else
    {
        G4cout << "WARNING: material \"" << materialChoice << "\" doesn't exist in NIST elements/materials"
            " table [located in $G4INSTALL/source/materials/src/G4NistMaterialBuilder.cc]" << G4endl;
        G4cout << "Use command \"/parameter/nist\" to see full materials list!" << G4endl;
    }
}

void PassiveProtonBeamLine::SetModulatorAngle(G4double value)
{
    modulator -> SetModulatorAngle(value);
    //G4RunManager::GetRunManager() -> GeometryHasBeenModified();
}