#include <iostream>
#include <cstdlib>
#define THICKNESS 0.001		// define the thickness of air and aluminum in each forward calorimeter

#include "FCALDetectorConstruction.hh"

#include "G4Element.hh"
#include "G4ElementTable.hh"
#include "G4Material.hh"
#include "G4MaterialTable.hh"
#include "G4NistManager.hh"
#include "G4Box.hh"
#include "G4Tubs.hh"
#include "G4LogicalVolume.hh"
#include "G4ThreeVector.hh"
#include "G4RotationMatrix.hh"
#include "G4Transform3D.hh"
#include "G4PVPlacement.hh"
#include "G4LogicalBorderSurface.hh"
#include "G4OpBoundaryProcess.hh"
#include "G4ios.hh"

FCALDetectorConstruction::FCALDetectorConstruction(G4int row, G4int column, G4double abs_Len)
: n_Rows(row), n_Columns(column), abs_Length(abs_Len)
{
	n_Cells = n_Rows * n_Columns;		// the number of the cells to be made
	front_Space = 11.5*cm;

//			-------- 1. LeadGlass: You must be carefule because all the length means half of the actual length.
	x_LeadGlass = 2.*cm;
	y_LeadGlass = 2.*cm;
	z_LeadGlass = 22.5*cm;
//	z_LeadGlass = 0.5*cm;

//			-------- 2. Air Gap and Foil
	x_AirGap = x_LeadGlass + 1.*THICKNESS*cm;
	y_AirGap = y_LeadGlass + 1.*THICKNESS*cm;
	z_AirGap = z_LeadGlass;

	x_FoilThick = x_AirGap + 1.*THICKNESS*cm;
	y_FoilThick = y_AirGap + 1.*THICKNESS*cm;
	z_FoilThick = z_AirGap;

//			-------- 3. PlexiGlass, PMT
	r_Min = 0.*cm, r_Max = 1.5*cm;
	angle_Min = 0.*deg, angle_Max = 360.*deg;
	r_cathode = 1.2*cm;
	r_AirMax = r_Max + 1.*THICKNESS*cm;
	r_FoilMax = r_Max + 2.*THICKNESS*cm;

	h_Cookie = 0.01*cm, h_PlexiGlass = 2.*cm, h_PMTGlass = 0.1*cm, h_Cathode = 0.01*cm;
	h_Total = h_Cookie + h_PlexiGlass + h_PMTGlass + h_Cathode;

//			-------- 4. World 
	x_World = 1.* n_Columns * x_FoilThick; 
	y_World = 1.* n_Rows * y_FoilThick;
	z_World = z_FoilThick + h_Total + front_Space;			

}

FCALDetectorConstruction::~FCALDetectorConstruction() {;}

G4VPhysicalVolume* FCALDetectorConstruction::Construct()
{
//			----------------------------------------------------------------------
//			-------- 1. Elements and Materials 
//			----------------------------------------------------------------------
	G4double a, z, density;			// a: mass of a mole, z: mean number of protons
	G4int n_Elements;			// n_Elements: the number of kinds of the elements
	
	G4NistManager* nist_Manager = G4NistManager::Instance();

//			-------- 1.1. Interface Cookie Material
	G4Element* O		= new G4Element("Oxygen", "O", z=8, a=16.00*g/mole);
	G4Element* Si		= new G4Element("Silicon", "Si", z=14, a=28.09*g/mole);
	G4Material* cookieMat 	= new G4Material("InterfaceCookie", density=2.32*g/cm3, n_Elements=2);					// Interface cookie
	cookieMat->AddElement(O, 2);
	cookieMat->AddElement(Si, 1);

//			-------- 1.2. Materials for various parts of the detector
	G4Material* air		= nist_Manager->FindOrBuildMaterial("G4_AIR");				// Air
	G4Material* Al		= nist_Manager->FindOrBuildMaterial("G4_Al");				// Aluminum
	G4Material* leadGlass	= nist_Manager->FindOrBuildMaterial("G4_GLASS_LEAD");			// Leadglass
	G4Material* plexiGlass	= nist_Manager->FindOrBuildMaterial("G4_PLEXIGLASS");			// Plexiglass
	G4Material* pmtGlass	= nist_Manager->FindOrBuildMaterial("G4_SILICON_DIOXIDE");		// PMTWindow


//			----------------------------------------------------------------------
//			-------- 2. Generate and Add Material Property Tables
//			----------------------------------------------------------------------
	const G4int num = 3;

	G4double eRange_Photon[num]		= {1.24*eV, 3.26*eV, 3.76*eV};	// the energy range of optical photons for which optical properties are described
	G4double absLength_LeadGlass[num]	= {abs_Length*cm, abs_Length*cm, 0.1*cm} ;	// Absorption length of Lead-Glass		absLength = absorption length
	G4double absLength_PlexiGlass[num]	= {200.*cm, 200.*cm, 0.1*cm};	// Absorption length of plexiglass

	G4double index_Air[num]			= {1., 1., 1.};		// Index of refraction of air
	G4double index_LeadGlass[num]		= {1.62, 1.62, 1.62};	// Index of refraction of LeadGlass
	G4double index_PMTGlass[num]		= {1.48, 1.48, 1.48};	// Index of refraction of PMTWindow
	G4double index_Cookie[num]		= {1.43, 1.43, 1.43};	// Index of refraction of the interface cookie
	G4double index_PlexiGlass[num]		= {1.48, 1.48, 1.48};	// Index of refraction of Plexiglass

	G4double eff_Cathode[num]		= {1., 1., 1.};		// eff = efficiency	
	G4double ref_Cathode[num]		= {0., 0., 0.};		// ref = reflectivity

//			-------- 2.1. Table - Air
	G4MaterialPropertiesTable* air_MPT = new G4MaterialPropertiesTable();		// MPT = Material Properties Table
		air_MPT->AddProperty("RINDEX", eRange_Photon, index_Air, num);
		air->SetMaterialPropertiesTable(air_MPT);

//			-------- 2.2. Table - LeadGlass
	G4MaterialPropertiesTable* leadGlass_MPT = new G4MaterialPropertiesTable();
		leadGlass_MPT->AddProperty("RINDEX", eRange_Photon, index_LeadGlass, num);
		leadGlass_MPT->AddProperty("ABSLENGTH", eRange_Photon, absLength_LeadGlass, num);
		leadGlass->SetMaterialPropertiesTable(leadGlass_MPT);

//			-------- 2.3. Table - interface cookie
	G4MaterialPropertiesTable* cookie_MPT = new G4MaterialPropertiesTable();
		cookie_MPT->AddProperty("RINDEX", eRange_Photon, index_Cookie, num);
		cookieMat->SetMaterialPropertiesTable(cookie_MPT);

//			-------- 2.4. Table - plexiglass
	G4MaterialPropertiesTable* plexiGlass_MPT = new G4MaterialPropertiesTable();
		plexiGlass_MPT->AddProperty("RINDEX", eRange_Photon, index_PlexiGlass, num);
		plexiGlass_MPT->AddProperty("ABSLENGTH", eRange_Photon, absLength_PlexiGlass, num);
		plexiGlass->SetMaterialPropertiesTable(plexiGlass_MPT);

//			-------- 2.5. Table - PMTglass
	G4MaterialPropertiesTable* pmtGlass_MPT = new G4MaterialPropertiesTable();
		pmtGlass_MPT->AddProperty("RINDEX", eRange_Photon, index_PMTGlass, num);
		pmtGlass->SetMaterialPropertiesTable(pmtGlass_MPT);

//			-------- 2.6. Table - PMTcathode
	G4MaterialPropertiesTable* cathode_MPT = new G4MaterialPropertiesTable();
		cathode_MPT->AddProperty("EFFICIENCY", eRange_Photon, eff_Cathode, num);
		cathode_MPT->AddProperty("REFLECTIVITY", eRange_Photon, ref_Cathode, num);
		Al->SetMaterialPropertiesTable(cathode_MPT);


//			----------------------------------------------------------------------
//			-------- 3. Volumes
//			----------------------------------------------------------------------
//			-------- 3.1. World Volume
//			----------------------------------------------------------------------
	G4Box* world_box			= new G4Box("World", x_World, y_World, z_World);		// Definiton for the world volume
	G4LogicalVolume* world_log		= new G4LogicalVolume(world_box, air, "World");			// Definiton of the logical volume of the world volume filled with the air
	G4VPhysicalVolume* world_phys		= new G4PVPlacement(0, G4ThreeVector(), world_log, "World", 0, false, 0);		// Placement of the world volume


//			----------------------------------------------------------------------
//			-------- 3.2. LeadGlass
//			----------------------------------------------------------------------
	G4Box* foilBox_box			= new G4Box("FoilBox", x_FoilThick, y_FoilThick, z_FoilThick);
	G4LogicalVolume* foilBox_log		= new G4LogicalVolume(foilBox_box, Al, "FoilBox");
	G4VPhysicalVolume* foilBox_phys;
	for (G4int i=0; i<n_Cells; i++)
	{
		G4double p_x, p_y, p_z;
		G4int column_Num, row_Num;
		column_Num = i % n_Columns;
		row_Num = (i - column_Num) / n_Columns;

		p_x = -x_World + x_FoilThick + column_Num*2.*x_FoilThick;
		p_y = -y_World + y_FoilThick + row_Num*2.*y_FoilThick;
		p_z = -z_World + z_FoilThick + 2.*front_Space;

		foilBox_phys			= new G4PVPlacement(0, G4ThreeVector(p_x, p_y, p_z), foilBox_log, "FoilBox", world_log, false, i);
	}

	G4Box* airGap_box			= new G4Box("AirGap", x_AirGap, y_AirGap, z_AirGap);
	G4LogicalVolume* airGap_log		= new G4LogicalVolume(airGap_box, air, "AirGap");
	G4VPhysicalVolume* airGap_phys		= new G4PVPlacement(0, G4ThreeVector(), airGap_log, "AirGap", foilBox_log, false, 0);

	G4Box* leadGlass_box			= new G4Box("LeadGlass", x_LeadGlass, y_LeadGlass, z_LeadGlass);
	G4LogicalVolume* leadGlass_log		= new G4LogicalVolume(leadGlass_box, leadGlass, "LeadGlass");
	G4VPhysicalVolume* leadGlass_phys	= new G4PVPlacement(0, G4ThreeVector(), leadGlass_log, "LeadGlass", airGap_log, false, 0);

//			----------------------------------------------------------------------
//			-------- 3.3. PlexiGlass
//			----------------------------------------------------------------------
	G4Tubs* foilTub_tubs			= new G4Tubs("FoilTub", r_Min, r_FoilMax, h_Total, angle_Min, angle_Max);
	G4LogicalVolume* foilTub_log		= new G4LogicalVolume(foilTub_tubs, Al, "FoilTub");
	G4VPhysicalVolume* foilTub_phys;
	for (G4int i=0; i<n_Cells; i++)
	{
		G4double p_x, p_y, p_z;
		G4int column_Num, row_Num;
		column_Num = i % n_Columns;
		row_Num = (i - column_Num) / n_Columns;

		p_x = -x_World + x_FoilThick + column_Num*2.*x_FoilThick;
		p_y = -y_World + y_FoilThick + row_Num*2.*y_FoilThick;
		p_z = -z_World + 2.*z_FoilThick + 2.*front_Space + h_Total;

		foilTub_phys			= new G4PVPlacement(0, G4ThreeVector(p_x, p_y, p_z), foilTub_log, "FoilTub", world_log, false, i);
	}

	G4Tubs* airTub_tubs			= new G4Tubs("AirTub", r_Min, r_AirMax, h_Total, angle_Min, angle_Max);
	G4LogicalVolume* airTub_log		= new G4LogicalVolume(airTub_tubs, air, "AirTub");
	G4VPhysicalVolume* airTub_phys		= new G4PVPlacement(0, G4ThreeVector(), airTub_log, "AirTub", foilTub_log, false, 0);

	G4Tubs* cookie_tubs			= new G4Tubs("InterfaceCookie", r_Min, r_Max, h_Cookie, angle_Min, angle_Max);				// interface cookie
	G4LogicalVolume* cookie_log		= new G4LogicalVolume(cookie_tubs, cookieMat, "InterfaceCookie");
	G4VPhysicalVolume* cookie_phys		= new G4PVPlacement(0, G4ThreeVector(0, 0, -h_Total+h_Cookie), cookie_log, "InterfaceCookie", airTub_log, false, 0);
	
	G4Tubs* plexiGlass_tubs			= new G4Tubs("PlexiGlass", r_Min, r_Max, h_PlexiGlass, angle_Min, angle_Max);					// plexiglass 
	G4LogicalVolume* plexiGlass_log		= new G4LogicalVolume(plexiGlass_tubs, plexiGlass, "PlexiGlass");
	G4VPhysicalVolume* plexiGlass_phys	= new G4PVPlacement(0, G4ThreeVector(0, 0, -h_Total+2.*h_Cookie+h_PlexiGlass), plexiGlass_log, "PlexiGlass", airTub_log, false, 0);

	G4Tubs* pmtGlass_tubs			= new G4Tubs("PMTGlass", r_Min, r_Max, h_PMTGlass, angle_Min, angle_Max);					// pmtglass 
	G4LogicalVolume* pmtGlass_log		= new G4LogicalVolume(pmtGlass_tubs, pmtGlass, "PMTGlass");
	G4VPhysicalVolume* pmtGlass_phys	= new G4PVPlacement(0, G4ThreeVector(0, 0, -h_Total+2.*h_Cookie+2.*h_PlexiGlass+h_PMTGlass), pmtGlass_log, "PMTGlass", airTub_log, false, 0);

	G4Tubs* pmtCathode_tubs			= new G4Tubs("PMTCathode", r_Min, r_cathode, h_Cathode, angle_Min, angle_Max);				// pmt cathode
	G4LogicalVolume* pmtCathode_log		= new G4LogicalVolume(pmtCathode_tubs, Al, "PMTCathode");
	G4VPhysicalVolume* pmtCathode_phys	= new G4PVPlacement(0, G4ThreeVector(0, 0, -h_Total+2.*h_Cookie+2.*h_PlexiGlass+2.*h_PMTGlass+h_Cathode), pmtCathode_log, "PMTCathode", airTub_log, false, 0);


//			----------------------------------------------------------------------
//			-------- 4. Surfaces
//			----------------------------------------------------------------------
//			-------- 4.1. Generate and Add Surface Property Tables
//			----------------------------------------------------------------------
	G4double ref_PlexiGlass[num]	= {0.95, 0.95, 0.95};		// the reflectivity of plexiglass
	G4double eff_PlexiGlass[num]	= {0.05, 0.05, 0.05};		// plexi_glass_reflectivity + plexi_glass_efficiency = 1
	G4double ref_LeadGlass[num]	= {0.85, 0.85, 0.85};
	G4double eff_LeadGlass[num]	= {0.15, 0.15, 0.15};

	G4double specular_Lobe[num]	= {0.045, 0.045, 0.045};
	G4double specular_Spike[num]	= {0.950, 0.950, 0.950};
	G4double back_Scatter[num]	= {0.005, 0.005, 0.005};		// specular_lobe + specular_spike + back_scatter = 1


//			----------------------------------------------------------------------
//			-------- 4.2. Surfaces
//			----------------------------------------------------------------------
//			-------- 4.2.1. Surface of the leadglass in contact with air
	G4OpticalSurface* leadGlassAir_Surf		= new G4OpticalSurface("LeadGlassAirSurface");
	leadGlassAir_Surf->SetType(dielectric_dielectric);
	leadGlassAir_Surf->SetFinish(ground);
	leadGlassAir_Surf->SetModel(unified);
	leadGlassAir_Surf->SetSigmaAlpha(0.12);
	G4LogicalBorderSurface* leadGlassAir_Surf_log	= new G4LogicalBorderSurface("LeadGlassAirSurface", leadGlass_phys, airGap_phys, leadGlassAir_Surf);

//			-------- 4.2.2 Surface of the plexiglass in contact with air
	G4OpticalSurface* plexiGlassAir_Surf		= new G4OpticalSurface("PlexiGlassAirSurface");
	plexiGlassAir_Surf->SetType(dielectric_dielectric);
	plexiGlassAir_Surf->SetFinish(polished);
	plexiGlassAir_Surf->SetModel(unified);
	G4LogicalBorderSurface* plexiGlassAir_Surf_log 	= new G4LogicalBorderSurface("PlexiGlassAirSurface", plexiGlass_phys, airTub_phys, plexiGlassAir_Surf);

//			-------- 4.2.3 Surface of air in contact with the Aluminum foil(LeadGlass)
	G4OpticalSurface* leadGlass_AirFoil_Surf	= new G4OpticalSurface("LeadGlassAirFoilSurface");
	leadGlass_AirFoil_Surf->SetType(dielectric_metal);
	leadGlass_AirFoil_Surf->SetFinish(polished);
	leadGlass_AirFoil_Surf->SetModel(unified);
	G4LogicalBorderSurface* leadGlass_AirFoil_Surf_log = new G4LogicalBorderSurface("LeadGlassAirFoilSurface", airGap_phys, foilBox_phys, leadGlass_AirFoil_Surf);

//			-------- 4.2.4 Surface of air in contact with the Aluminum foil(PlexiGlass)
	G4OpticalSurface* plexiGlass_AirFoil_Surf	= new G4OpticalSurface("PlexiGlassAirFoilSurface");
	plexiGlass_AirFoil_Surf->SetType(dielectric_metal);
	plexiGlass_AirFoil_Surf->SetFinish(polished);
	plexiGlass_AirFoil_Surf->SetModel(unified);
	G4LogicalBorderSurface* plexiGlass_AirFoil_Surf_log = new G4LogicalBorderSurface("PlexiGlassAirFoilSurface", airTub_phys, foilTub_phys, plexiGlass_AirFoil_Surf);

	
//			----------------------------------------------------------------------
//			-------- 4.3. Property Tables
//			----------------------------------------------------------------------
//			-------- 4.3.1. PlexiGlass(Air-Foil)
	G4MaterialPropertiesTable* plexiGlass_AirFoil_MPT = new G4MaterialPropertiesTable();
	plexiGlass_AirFoil_MPT->AddProperty("REFLECTIVITY", eRange_Photon, ref_PlexiGlass, num);
	plexiGlass_AirFoil_MPT->AddProperty("EFFICIENCY", eRange_Photon, eff_PlexiGlass, num);
	
//			-------- 4.1.2. LeadGlass(Air-Foil)
	G4MaterialPropertiesTable* leadGlass_AirFoil_MPT = new G4MaterialPropertiesTable();
	leadGlass_AirFoil_MPT->AddProperty("REFLECTIVITY", eRange_Photon, ref_LeadGlass, num);
	leadGlass_AirFoil_MPT->AddProperty("EFFICIENCY", eRange_Photon, eff_LeadGlass, num);

//			-------- 4.1.3. PlexiGlass(PlexiGlass-Air)
	G4MaterialPropertiesTable* plexiGlass_Air_MPT = new G4MaterialPropertiesTable();

//			-------- 4.1.4. LeadGlass(LeadGlass_Air)
	G4MaterialPropertiesTable* leadGlass_Air_MPT = new G4MaterialPropertiesTable();
	leadGlass_Air_MPT->AddProperty("SPECULARLOBECONSTANT", eRange_Photon, specular_Lobe, num);
	leadGlass_Air_MPT->AddProperty("SPECULARSPIKECONSTANT", eRange_Photon, specular_Spike, num);
	leadGlass_Air_MPT->AddProperty("BACKSCATTERCONSTANT", eRange_Photon, back_Scatter, num);

	
	leadGlassAir_Surf->SetMaterialPropertiesTable(leadGlass_Air_MPT);
	plexiGlassAir_Surf->SetMaterialPropertiesTable(plexiGlass_Air_MPT);
	leadGlass_AirFoil_Surf->SetMaterialPropertiesTable(leadGlass_AirFoil_MPT);
	plexiGlass_AirFoil_Surf->SetMaterialPropertiesTable(plexiGlass_AirFoil_MPT);

	G4cout << "End of the Geometry Construction" << G4endl;

	return world_phys;  // return the physical volume of the world volume
}