/** * @file Phase1WLSFibre.cc * @author Ibrahin Pinera * @date 2016 SoLid - University of Antwerp */ #include "Phase1WLSFibre.hh" #include "SolidMPPC.hh" #include "SolidMaterials.hh" #include "G4Material.hh" #include "G4MaterialTable.hh" #include "G4LogicalBorderSurface.hh" #include "G4OpBoundaryProcess.hh" #include "G4ThreeVector.hh" #include "G4VisAttributes.hh" #include "G4Transform3D.hh" #include "G4UnionSolid.hh" #include "G4Box.hh" #include "G4SubtractionSolid.hh" #include "G4Tubs.hh" #include "G4PVPlacement.hh" #include "G4UnitsTable.hh" #include "G4LogicalBorderSurface.hh" #include "G4SystemOfUnits.hh" #include "G4PhysicalConstants.hh" #include Phase1WLSFibre::Phase1WLSFibre(G4int verboseLevel, G4int numberMPPCs, G4int numberCladdings) { m_verboseLevel = verboseLevel; m_numberMPPCs = numberMPPCs; m_numberCladdings = numberCladdings; m_mirrorThickness = 25.*um; m_greaseThickness = 50.*um; m_coreLength = 922.4*mm; // fiber total length m_fiberSection = 3.*mm; // fiber total section m_clad1Thickness = 0.04*m_fiberSection; m_clad2Thickness = 0.0; if(m_numberCladdings==2) m_clad2Thickness = 0.02*m_fiberSection; m_clad1Section = m_fiberSection - 2.*m_clad2Thickness; m_coreSection = m_clad1Section - 2.*m_clad1Thickness; mppcConstruction = new SolidMPPC(m_verboseLevel); m_mppcThickness = mppcConstruction->GetThickness(); m_fiberLength = m_coreLength + m_mppcThickness + m_greaseThickness; if(m_numberMPPCs==2) // the fiber have 2 SiPMs m_fiberLength += m_mppcThickness + m_greaseThickness; else // the fiber have only 1 MPPC m_fiberLength += m_mirrorThickness; if( m_verboseLevel > 0 ) { G4cout << "==============================================" << G4endl; G4cout << "number of claddings: " << m_numberCladdings << G4endl; G4cout << "number of SiPM/fiber: " << m_numberMPPCs << G4endl; G4cout << "external size of the WLS fiber: "; G4cout << m_fiberLength << " x " << m_fiberSection << " x " << m_fiberSection << G4endl; G4cout << "size of the fiber core: "; G4cout << m_coreLength << " x " << m_coreSection << " x " << m_coreSection << G4endl; G4cout << "cladding1 (internal) thickness: " << m_clad1Thickness << G4endl; G4cout << "cladding2 (external) thickness: " << m_clad2Thickness << G4endl; G4cout << "mirror thickness: " << m_mirrorThickness << G4endl; G4cout << "grease thickness: " << m_greaseThickness << G4endl; G4cout << "SiPM thickness: " << m_mppcThickness << G4endl; G4cout << "==============================================" << G4endl; G4cout << "section of fiber:" << G4endl; G4cout << " ------------------------------" << G4endl; G4cout << " | cladding2 |" << G4endl; G4cout << " | ---------------------- |" << G4endl; G4cout << " | | cladding1 | |" << G4endl; G4cout << " | | ---------- | |" << G4endl; G4cout << " | | | | | |" << G4endl; G4cout << " | | | core | | |" << G4endl; G4cout << " | | | | | |" << G4endl; G4cout << " | | ---------- | |" << G4endl; G4cout << " | | | |" << G4endl; G4cout << " | ---------------------- |" << G4endl; G4cout << " | |" << G4endl; G4cout << " ------------------------------" << G4endl; G4cout << "==============================================" << G4endl; } } Phase1WLSFibre::~Phase1WLSFibre() { } // numberMPPC: number of MPPCs per fiber // numberCladding: number of claddings per fiber G4LogicalVolume *Phase1WLSFibre::WLSConstruction() { G4Material * mate; // make colours ********************************************************* G4Colour white (1.0, 1.0, 1.0) ; G4Colour white_t (1.0, 1.0, 1.0, .95) ; G4Colour grey (0.5, 0.5, 0.5) ; G4Colour lgrey (.75, .75, .75) ; G4Colour red (1.0, 0.0, 0.0) ; G4Colour blue (0.0, 0.0, 1.0) ; G4Colour blue_t (0.0, 0.0, 1.0, .75) ; G4Colour cyan (0.0, 1.0, 1.0, 0.5) ; G4Colour magenta (1.0, 0.0, 1.0) ; G4Colour yellow (1.0, 1.0, 0.0) ; G4Colour lblue (0.0, 0.0, .75) ; G4Colour black (0.0, 0.0, 0.0) ; G4Colour green (0.0, 1.0, 0.0) ; G4Colour lgreen (0.0, .75, 0.0) ; G4Colour lgreen2 (0.0, .85, 0.0) ; //Build MPPC mate = G4Material::GetMaterial("MPPCFilm"); logFibreMPPC = mppcConstruction->MPPCConstruction(); G4VisAttributes* vatMPPC = new G4VisAttributes(yellow); vatMPPC->SetVisibility(true); vatMPPC->SetForceSolid(true); logFibreMPPC->SetVisAttributes(vatMPPC); //Build fiber total along the Y axe (fiber+MPPC+Almirror+cladding+grease) G4Box * FibBox = new G4Box("FibBox", 0.5*m_fiberSection, 0.5*m_fiberLength, 0.5*m_fiberSection); mate = G4Material::GetMaterial("Air"); logFibre = new G4LogicalVolume(FibBox, mate, "logFibre", 0, 0, 0); // build the core ************************************************************* G4Box * CoreBox = new G4Box("CoreBox", 0.5*m_coreSection, 0.5*m_coreLength, 0.5*m_coreSection); mate = G4Material::GetMaterial("Polystyrene"); logFibreCore = new G4LogicalVolume(CoreBox, mate, "logFibreCore", 0, 0, 0); //Color Settings G4VisAttributes* vatCore = new G4VisAttributes(green); vatCore->SetVisibility(true); vatCore->SetForceSolid(true); logFibreCore->SetVisAttributes(vatCore); // build the first cladding (internal) G4Box * CladBoxOuter = new G4Box("CladBoxOuter", 0.5*m_clad1Section, 0.5*m_coreLength, 0.5*m_clad1Section); G4Box * CladBoxInner = new G4Box("CladBoxInner", 0.5*m_coreSection, 0.5*m_coreLength + 0.5*mm, 0.5*m_coreSection); G4SubtractionSolid * CladBox = new G4SubtractionSolid("CladBox", CladBoxOuter, CladBoxInner, 0, G4ThreeVector(0,0,0)); mate = G4Material::GetMaterial("PMMA"); logFibreClad = new G4LogicalVolume(CladBox, mate, "logFibreClad", 0, 0, 0); G4VisAttributes* vatClad = new G4VisAttributes(lgreen); vatClad->SetVisibility(true); vatClad->SetForceSolid(true); logFibreClad->SetVisAttributes(vatClad); // build the external 2nd cladding if exist G4Box * CladBoxOuter2 = new G4Box("CladBoxOuter2", 0.5*m_fiberSection, 0.5*m_coreLength, 0.5*m_fiberSection); G4Box * CladBoxInner2 = new G4Box("CladBoxInner2", 0.5*m_clad1Section, 0.5*m_coreLength + 0.5*mm, 0.5*m_clad1Section); G4SubtractionSolid * CladBox2 = new G4SubtractionSolid("CladBox2", CladBoxOuter2, CladBoxInner2, 0, G4ThreeVector(0,0,0)); if( m_numberCladdings == 2 ) { mate = G4Material::GetMaterial("PMMA"); logFibreClad2 = new G4LogicalVolume(CladBox2, mate, "logFibreClad2", 0, 0, 0); G4VisAttributes* vatClad2 = new G4VisAttributes(lgreen2); vatClad2->SetVisibility(true); vatClad2->SetForceSolid(true); logFibreClad2->SetVisAttributes(vatClad2); } //Build Grease G4Box * GreaseBox = new G4Box("GreaseBox", 0.5*m_fiberSection, 0.5*m_greaseThickness, 0.5*m_fiberSection); mate = G4Material::GetMaterial("Grease"); logFibreGrease = new G4LogicalVolume(GreaseBox, mate, "logFibreGrease", 0, 0, 0); G4VisAttributes* vatGrease = new G4VisAttributes(lgrey); vatGrease->SetVisibility(true); vatGrease->SetForceSolid(true); logFibreGrease->SetVisAttributes(vatGrease); //Build Mirror G4Box * MirrorBox = new G4Box("MirrorBox", 0.5*m_fiberSection, 0.5*m_mirrorThickness, 0.5*m_fiberSection); mate = G4Material::GetMaterial("Mirror"); logFibreMirror = new G4LogicalVolume(MirrorBox, mate, "logFibreMirror", 0, 0, 0); G4VisAttributes* vatMirror = new G4VisAttributes(white); vatMirror->SetVisibility(true); vatMirror->SetForceSolid(true); logFibreMirror->SetVisAttributes(vatMirror); // placing physical volumes physFibreMPPC = new G4PVPlacement(0, G4ThreeVector(0., 0.5*(m_fiberLength - m_mppcThickness), 0.), logFibreMPPC, "physFibreMPPC", logFibre, false, 0); physFibreGrease = new G4PVPlacement(0, G4ThreeVector(0., 0.5*(m_fiberLength - m_greaseThickness) - m_mppcThickness, 0.), logFibreGrease, "physFibreGrease", logFibre, false, 0); G4double core_offSet = 0.5*(m_mirrorThickness - m_greaseThickness - m_mppcThickness); if( m_numberMPPCs == 1 ) { // the fiber have only 1 MPPC, and then an Al mirror on the other extreme physFibreMirror = new G4PVPlacement(0, G4ThreeVector(0., -0.5*m_fiberLength + 0.5*m_mirrorThickness, 0.), logFibreMirror, "physFibreMirror", logFibre, false, 0); }else{ // the fiber have 2 MPPCs core_offSet = 0.; physFibreMPPC = new G4PVPlacement(0, G4ThreeVector(0., -0.5*(m_fiberLength - m_mppcThickness), 0.), logFibreMPPC, "physFibreMPPC", logFibre, false, 1); physFibreGrease = new G4PVPlacement(0, G4ThreeVector(0., -0.5*m_coreLength - 0.5*m_greaseThickness, 0.), logFibreGrease, "physFibreGrease", logFibre, false, 1); } physFibreCore = new G4PVPlacement(0, G4ThreeVector(0., core_offSet, 0.), logFibreCore, "physFibreCore", logFibre, false, 0); physFibreClad = new G4PVPlacement(0, G4ThreeVector(0., core_offSet, 0.), logFibreClad, "physFibreClad", logFibre, false, 0); if( m_numberCladdings == 2 ) physFibreClad2 = new G4PVPlacement(0, G4ThreeVector(0., core_offSet, 0.), logFibreClad2, "physFibreClad2", logFibre, false, 0); if( m_verboseLevel > 0 ) { G4float fiberCore_mass = logFibreCore->GetMass()/g; G4float fiberCore_vol = CoreBox->GetCubicVolume()/cm3; G4float fiberClad1_mass = logFibreClad->GetMass()/g; G4float fiberClad1_vol = CladBox->GetCubicVolume()/cm3; G4float grease_mass = logFibreGrease->GetMass()/g; G4float grease_vol = GreaseBox->GetCubicVolume()/cm3; G4float mirror_mass = logFibreMirror->GetMass()/g; G4float mirror_vol = MirrorBox->GetCubicVolume()/cm3; G4float mppc_mass = logFibreMPPC->GetMass()/g; G4float mppc_vol = m_mppcThickness*m_fiberSection*m_fiberSection/cm3; G4float total_mass = fiberCore_mass + fiberClad1_mass + m_numberMPPCs*grease_mass + (2.-m_numberMPPCs)*mirror_mass + m_numberMPPCs*mppc_mass; G4float total_vol = fiberCore_vol + fiberClad1_vol + m_numberMPPCs*grease_vol + (2.-m_numberMPPCs)*mirror_vol + m_numberMPPCs*mppc_vol; G4cout << "**************VOLUME MASS****************" << G4endl; G4cout << " part: mass (g)\t volume (cm^3)" << G4endl; G4cout << " fiber core: " << fiberCore_mass << ", \t" << fiberCore_vol << G4endl; G4cout << " fiber cladding1: " << fiberClad1_mass << ", \t" << fiberClad1_vol << G4endl; if( m_numberCladdings == 2 ) { G4float fiberClad2_mass = logFibreClad2->GetMass()/g; G4float fiberClad2_vol = CladBox2->GetCubicVolume()/cm3; G4cout << " fiber cladding2: " << fiberClad2_mass << ", \t" << fiberClad2_vol << G4endl; total_mass += fiberClad2_mass; total_vol += fiberClad2_vol; } G4cout << " grease: " << grease_mass << ", \t" << grease_vol << G4endl; G4cout << " mirror: " << mirror_mass << ", \t" << mirror_vol << G4endl; G4cout << " SiPM: " << mppc_mass << ", \t" << mppc_vol << G4endl; G4cout << " ALL fiber total: " << total_mass << ", \t" << total_vol << G4endl; G4cout << " entire Phase1: " << 64.*50.*total_mass << ", \t" << 64.*50.*total_vol << G4endl; G4cout << "**********END VOLUME MASS****************" << G4endl; } return logFibre; }