/** * @file Phase1Plane.cc * @author Ibrahin Pinera * @date 2017 SoLid - University of Antwerp */ #include "Phase1Plane.hh" #include "SolidMaterials.hh" #include "Phase1AlFrame.hh" #include "Phase1ElectronicBox.hh" #include "G4Material.hh" #include "G4MaterialTable.hh" #include "G4ThreeVector.hh" #include "G4VisAttributes.hh" #include "G4Transform3D.hh" #include "G4Box.hh" #include "G4SubtractionSolid.hh" #include "G4Tubs.hh" #include "G4PVPlacement.hh" #include "G4UnitsTable.hh" #include "G4SystemOfUnits.hh" #include "G4PhysicalConstants.hh" #include Phase1Plane::Phase1Plane(G4int verboseLevel, G4int nCubesX, G4int nCubesY, G4double gapX, G4double gapY) { m_verboseLevel = verboseLevel; m_nCubesX = nCubesX; m_nCubesY = nCubesY; m_gapX = gapX; m_gapY = gapY; //DetMod (cube + Li + Tyvek) m_scintConstruction = new Phase1Scintillator(m_verboseLevel); logDetMod = m_scintConstruction->ScintConstruction(0); logDetModExternal = m_scintConstruction->ScintConstruction(1); m_detModSizeX = m_scintConstruction->GetDetModX(); m_detModSizeY = m_scintConstruction->GetDetModY(); m_detModSizeZ = m_scintConstruction->GetDetModZ(); m_pvtOffsetZ = m_scintConstruction->GetPVTOffsetZ(); m_grooveOffsetZ1 = m_scintConstruction->GetGrooveOffsetZ1(); m_grooveOffsetZ2 = m_scintConstruction->GetGrooveOffsetZ2(); m_grooveOffsetX1 = m_scintConstruction->GetGrooveOffsetX1(); m_grooveOffsetX2 = m_scintConstruction->GetGrooveOffsetX2(); m_grooveOffsetY1 = m_scintConstruction->GetGrooveOffsetY1(); m_grooveOffsetY2 = m_scintConstruction->GetGrooveOffsetY2(); // Construct logical WLS fibers G4cout << G4endl << "**************WLS fibers construction****************" << G4endl << G4endl; m_fiberConstruction = new Phase1WLSFibre(m_verboseLevel, 1, 2); //fiber with 1 SiPM and double cladding logFibre = m_fiberConstruction->WLSConstruction(); m_fiberLength = m_fiberConstruction->GetFiberLength(); // fiber length G4cout << G4endl << "**************END WLS fibers construction****************" << G4endl << G4endl; m_fibersGapX = m_detModSizeX + m_gapX; // gap between holes in the X direction m_fibersGapY = m_detModSizeY + m_gapY; // gap between holes in the Y direction m_holesSeparationX = m_grooveOffsetX1 - m_grooveOffsetX2; // between fibers in same cube in X direction m_holesSeparationY = m_grooveOffsetY1 - m_grooveOffsetY2; // between fibers in same cube in Y direction //x-dimension --> Width //y-dimension --> Length //z-dimension --> Thickness // Al Frame m_AlFrameConstruction = new Phase1AlFrame(m_verboseLevel, m_nCubesX, m_nCubesY, m_fibersGapX, m_fibersGapY, m_holesSeparationX, m_holesSeparationY); m_AlFrameWidth = m_AlFrameConstruction->GetWidth(); m_AlFrameFullLength = m_AlFrameConstruction->GetFullLength(); m_AlFrameFullHeight = m_AlFrameConstruction->GetFullHeight(); m_AlFrameInternalLength = m_AlFrameConstruction->GetInternalLength(); m_AlFrameInternalHeight = m_AlFrameConstruction->GetInternalHeight(); m_AlFrameOffSetX = m_AlFrameConstruction->GetOffSetX(); m_AlFrameOffSetY = m_AlFrameConstruction->GetOffSetY(); G4double holesOffsetZ1 = m_AlFrameConstruction->GetHolesOffSetZ1(); G4double holesOffsetZ2 = m_AlFrameConstruction->GetHolesOffSetZ2(); //PE tight sealing reflector bars G4cout << G4endl << "**************PE sealing reflector bars construction****************" << G4endl << G4endl; m_horPEReflectorWidth = 0.5*(m_AlFrameInternalHeight - 16.*m_fibersGapY); m_verPEReflectorWidth = 0.5*(m_AlFrameInternalLength - 16.*m_fibersGapX); // from external +X of Al frame (150.6*mm from design) m_holesOffsetX = 2.*m_AlFrameWidth + m_verPEReflectorWidth + 0.5*m_fibersGapX - m_grooveOffsetX1; // from external bottom of Al frame (102.7*mm from design) m_holesOffsetY = m_AlFrameWidth + m_horPEReflectorWidth + 0.5*m_fibersGapY + m_grooveOffsetY1; m_PEreflectorConstruction = new Phase1PEReflector(m_verboseLevel); logPlanePEReflector = m_PEreflectorConstruction->PEReflectorConstruction(m_horPEReflectorWidth, m_verPEReflectorWidth, m_nCubesX, m_nCubesY, m_fibersGapX, m_fibersGapY, m_holesSeparationX, m_holesSeparationY, holesOffsetZ1, holesOffsetZ2); G4cout << G4endl << "**************END PE sealing reflector bars construction****************" << G4endl << G4endl; logAlFrame = m_AlFrameConstruction->AlFrameConstruction(m_holesOffsetX, m_holesOffsetY); //Electronic Box G4cout << G4endl << "**************Electronic box construction****************" << G4endl << G4endl; m_electBoxConstruction = new Phase1ElectronicBox(m_verboseLevel); logElectronicBox = m_electBoxConstruction->ConstructElectronicBox(); m_electronicBoxLength = m_electBoxConstruction->GetFullLength(); // m_electronicBoxOffSetX = 97.955*mm; m_electronicBoxOffSetY = 0.; G4cout << G4endl << "**************END Electronic box construction****************" << G4endl << G4endl; m_planeFullWidth = 51.5*mm; // fix width for Planes (Al frame + Tyvek covers) m_planeFullLength = m_AlFrameFullLength + m_electronicBoxLength; // outer horizontal measure m_planeFullHeight = m_AlFrameFullHeight; // outer vertical measure m_planeOffSetX = m_AlFrameOffSetX - 0.5*m_electronicBoxLength; m_planeOffSetY = m_AlFrameOffSetY; m_coverLength = m_AlFrameInternalLength + 40.*mm; m_coverHeight = m_AlFrameInternalHeight + 40.*mm; m_coverThickness = 0.3*mm; if( m_verboseLevel > 0 ) { G4cout << "=== dimensions of the detector plane ===" << G4endl; G4cout << m_planeFullLength << " x " << m_planeFullHeight << " x " << m_planeFullWidth << G4endl; G4cout << "========================================" << G4endl; G4cout << "=== dimensions of the Tyvek covers for plane ===" << G4endl; G4cout << m_coverLength << " x " << m_coverHeight << " x " << m_coverThickness << G4endl; G4cout << "========================================" << G4endl; } } Phase1Plane::~Phase1Plane() { } G4LogicalVolume *Phase1Plane::PlaneConstruction() { 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, 0.5) ; 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, 0.5) ; 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) ; // G4bool forceSolid = false; G4bool forceSolid = true; // construct full Plane box G4Box* planeBox = new G4Box("planeBox", 0.5*m_planeFullLength, 0.5*m_planeFullHeight, 0.5*m_planeFullWidth); mate = G4Material::GetMaterial("Air"); logPlane = new G4LogicalVolume(planeBox, mate, "logPlane", 0, 0, 0); logPlane->SetVisAttributes(G4VisAttributes::Invisible); //==================================================== //locating the Al Frame //==================================================== physAlFrame = new G4PVPlacement(0, G4ThreeVector(-0.5*m_planeFullLength + 0.5*m_AlFrameFullLength, 0., 0.), logAlFrame, "physAlFrame", logPlane, false, 0); //==================================================== //locating the Al electronic box //==================================================== physElectronicBox = new G4PVPlacement(0, G4ThreeVector(0.5*m_planeFullLength - 0.5*m_electronicBoxLength, m_electronicBoxOffSetY, 0.), logElectronicBox, "physElectronicBox", logPlane, false, 0); //==================================================== // Locating the cube module (PVT + Tyvek + LiF:ZnS) //==================================================== //cube ID (XX-YY) //1st-2nd digits are the horizontal position in frame (0-15) //3rd-4th digits are the vertical position in frame (0-15) G4double posX = 0.; G4double posY = 0.; G4double posZ = 0.; int cubeID = 0; //Calculate Offset to make origin of ref system respect to detector 0,0,0 G4double offSetX = -0.5*(m_nCubesX - 1)*(m_detModSizeX + m_gapX); G4double offSetY = -0.5*(m_nCubesY - 1)*(m_detModSizeY + m_gapY); G4double offSetZ = 0.; for(int y = 0; y < m_nCubesY; y++) { for(int x = 0; x < m_nCubesX; x++) { cubeID = 100*x + y; //XX-YY posX = x*m_fibersGapX + m_planeOffSetX; posY = y*m_fibersGapY + m_planeOffSetY; posZ = 0; if( x == 0 || x == 15 || y == 0 || y == 15 ) // is external cube (no LiF:ZnS backing) physDetMod = new G4PVPlacement(0, G4ThreeVector(posX + offSetX, posY + offSetY, posZ + offSetZ), logDetModExternal, "physDetMod", logPlane, false, cubeID); else physDetMod = new G4PVPlacement(0, G4ThreeVector(posX + offSetX, posY + offSetY, posZ + offSetZ), logDetMod, "physDetMod", logPlane, false, cubeID); } } //==================================================== // Place physical WLS fibers in the frame //==================================================== //fiber ID (XXorYY-o/l) //1st-2nd digits are fiber # through the horizontal/vertical direction (YY/XX) (0-15) //3rd digit indicates the orientation of fibers (0 for horizontal, 1 for vertical) //4th digit indicates the location of the SiPM in the fiber (0 for bottom/left, 1 for top/right) int fibID; G4RotationMatrix* rotFibre = new G4RotationMatrix(); // Populate X direction Fibres (horizontal) int fiberOffSetX = 2.*mm; // shift even fibers to the right posX = (m_nCubesX - 1)*(0.5*m_fibersGapX) + m_planeOffSetX; for(int y = 0; y < m_nCubesY; y++) { // horizontal fiber on top -Z (MPPC on -X) (id -> YY-0/0) rotFibre->rotateZ(90.*deg); posY = y*m_fibersGapY + m_grooveOffsetY2 + m_planeOffSetY; posZ = m_grooveOffsetZ2 + m_pvtOffsetZ; fibID = 100*(y) + 0; //YY-0/0 physFibreX = new G4PVPlacement(rotFibre, G4ThreeVector(posX + offSetX - fiberOffSetX, posY + offSetY, posZ + offSetZ), logFibre, "physFibreX", logPlane, false, fibID); // horizontal fiber on bottom +Z (MPPC on +X) (id -> YY-0/1) posY = y*m_fibersGapY + m_grooveOffsetY1 + m_planeOffSetY; posZ = -m_grooveOffsetZ1 + m_pvtOffsetZ; fibID = 100*(y) + 1; //YY-0/1 rotFibre = new G4RotationMatrix(); rotFibre->rotateZ(-90.*deg); physFibreX = new G4PVPlacement(rotFibre, G4ThreeVector(posX + offSetX + fiberOffSetX, posY + offSetY, posZ + offSetZ), logFibre, "physFibreX", logPlane, false, fibID); rotFibre = new G4RotationMatrix(); } // Populate Y direction Fibres (vertical) int fiberOffSetY = 2.*mm; // shift odd fibers up posY = (m_nCubesY - 1)*(0.5*m_fibersGapY) + m_planeOffSetY; for(int x = 0; x < m_nCubesX; x++) { // vertical fiber on +X (MPPC on +Y) (id -> XX-1/1) posX = x*m_fibersGapX + m_grooveOffsetX1 + m_planeOffSetX; posZ = m_grooveOffsetZ1 + m_pvtOffsetZ; fibID = 100*(x) + 11; //XX-1/1 physFibreY = new G4PVPlacement(rotFibre, G4ThreeVector(posX + offSetX, posY + offSetY - fiberOffSetY, posZ + offSetZ), logFibre, "physFibreY", logPlane, false, fibID); // vertical fiber on -X (MPPC on -Y) (id -> XX-1/0) posX = x*m_fibersGapX + m_grooveOffsetX2 + m_planeOffSetX; posZ = -m_grooveOffsetZ2 + m_pvtOffsetZ; fibID = 100*(x) + 10; //XX-1/0 rotFibre = new G4RotationMatrix(); rotFibre->rotateX(180.*deg); physFibreY = new G4PVPlacement(rotFibre, G4ThreeVector(posX + offSetX, posY + offSetY + fiberOffSetY, posZ + offSetZ), logFibre, "physFibreY", logPlane, false, fibID); rotFibre = new G4RotationMatrix(); } //==================================================== // Construct Polyethylene tight sealing reflector bars //==================================================== physPlanePEReflector = new G4PVPlacement(0, G4ThreeVector(m_planeOffSetX, -0.5*m_AlFrameInternalHeight + 0.5*m_horPEReflectorWidth, 0.), logPlanePEReflector, "physPlanePEReflector", logPlane, false, 0); //==================================================== // Construct Tyvek cover sheets //==================================================== G4Box* TyvekCoverSheetBox = new G4Box("TyvekCoverSheetBox", 0.5*m_coverLength, 0.5*m_coverHeight, 0.5*m_coverThickness); mate = G4Material::GetMaterial("Tyvek"); logPlaneCover = new G4LogicalVolume(TyvekCoverSheetBox, mate, "logPlaneCover", 0, 0, 0); //Set Color G4VisAttributes* vatCoverSheet = new G4VisAttributes(white_t); vatCoverSheet->SetVisibility(true); vatCoverSheet->SetForceSolid(forceSolid); logPlaneCover->SetVisAttributes(vatCoverSheet); physPlaneCoverFront = new G4PVPlacement(0, G4ThreeVector(m_planeOffSetX, m_planeOffSetY, 0.5*m_planeFullWidth - 0.5*m_coverThickness), logPlaneCover, "physPlaneCoverFront", logPlane, false, 0); physPlaneCoverBack = new G4PVPlacement(0, G4ThreeVector(m_planeOffSetX, m_planeOffSetY, -0.5*m_planeFullWidth + 0.5*m_coverThickness), logPlaneCover, "physPlaneCoverBack", logPlane, false, 0); if( m_verboseLevel > 0 ) { G4float fibers_mass = 64.*logFibre->GetMass()/g; G4float detMods_mass = 196.*logDetMod->GetMass()/g + 60.*logDetModExternal->GetMass()/g; G4float AlFrame_mass = logAlFrame->GetMass()/g; G4float electBox_mass = logElectronicBox->GetMass()/g; G4float PEbars_mass = logPlanePEReflector->GetMass()/g; G4float cover_mass = 2.*logPlaneCover->GetMass()/g; G4float cover_vol = 2.*TyvekCoverSheetBox->GetCubicVolume()/cm3; G4float total_mass = fibers_mass + detMods_mass + AlFrame_mass + electBox_mass + PEbars_mass + cover_mass; G4float total_vol = planeBox->GetCubicVolume()/cm3; G4cout << "**************VOLUME MASS****************" << G4endl; G4cout << " part: mass (g)\t volume (cm^3)" << G4endl; G4cout << " Al frame: " << AlFrame_mass << G4endl; G4cout << " PE bars: " << PEbars_mass << G4endl; G4cout << " 64 fibers: " << fibers_mass << G4endl; G4cout << " 196(int)+60(ext) DetMods: " << detMods_mass << G4endl; G4cout << " Electronic box: " << electBox_mass << G4endl; G4cout << " 2 Tyvek covers: " << cover_mass << ", \t" << cover_vol << G4endl; G4cout << " entire plane sum: " << total_mass << ", \t" << total_vol << G4endl; G4cout << " entire plane total: " << logPlane->GetMass()/g << G4endl; G4cout << " 50 planes total: " << 50.*total_mass << ", \t" << 50.*total_vol << G4endl; G4cout << " 60 planes total: " << 60.*total_mass << ", \t" << 60.*total_vol << G4endl; G4cout << "**********END VOLUME MASS****************" << G4endl; } return logPlane; }