diff --git a/NPSimulation/src/GaspardTrackerDummyShape.cc b/NPSimulation/src/GaspardTrackerDummyShape.cc index e4dcf07c007ae19acd419c3172edefb1ba5fe622..1867e685e4793e3754e509def0597614f2514212 100644 --- a/NPSimulation/src/GaspardTrackerDummyShape.cc +++ b/NPSimulation/src/GaspardTrackerDummyShape.cc @@ -212,21 +212,23 @@ void GaspardTrackerDummyShape::VolumeMaker(G4int TelescopeNumber, //Uncomment to help debugging or if you want to understand the way the code work. //I should recommand to Comment it during simulation to avoid perturbation of simulation //Remember G4 is limitationg step on geometry constraints. - G4ThreeVector positionMarkerU = MMCenter*0.8 + MMu*FirstStageFace/4; - G4Box* solidMarkerU = new G4Box("solidMarkerU", FirstStageFace/4, 1*mm, 1*mm); - G4LogicalVolume* logicMarkerU = new G4LogicalVolume(solidMarkerU, Vacuum, "logicMarkerU", 0, 0, 0); - PVPBuffer = new G4PVPlacement(G4Transform3D(*MMrot,positionMarkerU), logicMarkerU, "MarkerU", world, false, 0); + /* + G4ThreeVector positionMarkerU = CT*0.98 + MMu*SiliconFace/4; + G4Box* solidMarkerU = new G4Box( "solidMarkerU" , SiliconFace/4 , 1*mm , 1*mm ) ; + G4LogicalVolume* logicMarkerU = new G4LogicalVolume( solidMarkerU , Vacuum , "logicMarkerU",0,0,0) ; + PVPBuffer = new G4PVPlacement(G4Transform3D(*MMrot,positionMarkerU),logicMarkerU,"MarkerU",world,false,0) ; - G4VisAttributes* MarkerUVisAtt= new G4VisAttributes(G4Colour(0.,0.,0.5)); //blue - logicMarkerU->SetVisAttributes(MarkerUVisAtt); + G4VisAttributes* MarkerUVisAtt= new G4VisAttributes(G4Colour(0.,0.,0.5));//blue + logicMarkerU->SetVisAttributes(MarkerUVisAtt); - G4ThreeVector positionMarkerV = MMCenter*0.8 + MMv*FirstStageFace/4; - G4Box* solidMarkerV = new G4Box("solidMarkerU", 1*mm, FirstStageFace/4, 1*mm); - G4LogicalVolume* logicMarkerV = new G4LogicalVolume(solidMarkerV, Vacuum, "logicMarkerV", 0, 0, 0); - PVPBuffer = new G4PVPlacement(G4Transform3D(*MMrot,positionMarkerV), logicMarkerV, "MarkerV", world, false, 0); + G4ThreeVector positionMarkerV = CT*0.98 + MMv*SiliconFace/4; + G4Box* solidMarkerV = new G4Box( "solidMarkerU" , 1*mm , SiliconFace/4 , 1*mm ) ; + G4LogicalVolume* logicMarkerV = new G4LogicalVolume( solidMarkerV , Vacuum , "logicMarkerV",0,0,0) ; + PVPBuffer = new G4PVPlacement(G4Transform3D(*MMrot,positionMarkerV),logicMarkerV,"MarkerV",world,false,0) ; - G4VisAttributes* MarkerVVisAtt= new G4VisAttributes(G4Colour(0.,0.5,0.)); //green - logicMarkerV->SetVisAttributes(MarkerVVisAtt); + G4VisAttributes* MarkerVVisAtt= new G4VisAttributes(G4Colour(0.,0.5,0.5));//green + logicMarkerV->SetVisAttributes(MarkerVVisAtt); + */ //////////////////////////////////////////////////////////////// ///////////////// First Stage Construction ///////////////////// @@ -241,6 +243,7 @@ void GaspardTrackerDummyShape::VolumeMaker(G4int TelescopeNumber, PVPBuffer = new G4PVPlacement(0, positionFirstStage, logicFirstStage, +// "G" + DetectorNumber + "FirstStage", Name + "_FirstStage", logicGPDDummyShape, false, @@ -250,7 +253,7 @@ void GaspardTrackerDummyShape::VolumeMaker(G4int TelescopeNumber, logicFirstStage->SetSensitiveDetector(m_FirstStageScorer); ///Visualisation of FirstStage Strip - G4VisAttributes* FirstStageVisAtt = new G4VisAttributes(G4Colour(0.0, 0.0, 0.9)); // blue + G4VisAttributes* FirstStageVisAtt = new G4VisAttributes(G4Colour(0.2, 0.2, 0.2)); logicFirstStage->SetVisAttributes(FirstStageVisAtt); } @@ -267,6 +270,7 @@ void GaspardTrackerDummyShape::VolumeMaker(G4int TelescopeNumber, PVPBuffer = new G4PVPlacement(0, positionSecondStage, logicSecondStage, +// "G" + DetectorNumber + "SecondStage", Name + "_SecondStage", logicGPDDummyShape, false, @@ -276,8 +280,8 @@ void GaspardTrackerDummyShape::VolumeMaker(G4int TelescopeNumber, logicSecondStage->SetSensitiveDetector(m_SecondStageScorer); ///Visualisation of SecondStage Strip - G4VisAttributes* SecondStageVisAtt = new G4VisAttributes(G4Colour(0.5, 0.5, 0.5)); - logicSecondStage->SetVisAttributes(SecondStageVisAtt); + G4VisAttributes* SecondStageVisAtt = new G4VisAttributes(G4Colour(0.5, 0.5, 0.5)) ; + logicSecondStage->SetVisAttributes(SecondStageVisAtt) ; } //////////////////////////////////////////////////////////////// @@ -293,6 +297,7 @@ void GaspardTrackerDummyShape::VolumeMaker(G4int TelescopeNumber, PVPBuffer = new G4PVPlacement(0, positionThirdStage, logicThirdStage, +// "G" + DetectorNumber + "ThirdStage", Name + "_ThirdStage", logicGPDDummyShape, false, @@ -302,7 +307,7 @@ void GaspardTrackerDummyShape::VolumeMaker(G4int TelescopeNumber, logicThirdStage->SetSensitiveDetector(m_ThirdStageScorer); ///Visualisation of Third Stage - G4VisAttributes* ThirdStageVisAtt = new G4VisAttributes(G4Colour(0.0, 0.9, 0.0)); // green + G4VisAttributes* ThirdStageVisAtt = new G4VisAttributes(G4Colour(0.7, 0.7, 0.7)); logicThirdStage->SetVisAttributes(ThirdStageVisAtt); } } @@ -534,10 +539,10 @@ void GaspardTrackerDummyShape::ConstructDetector(G4LogicalVolume* world) { G4RotationMatrix* MMrot = NULL ; G4ThreeVector MMpos = G4ThreeVector(0, 0, 0) ; -// G4ThreeVector MMu = G4ThreeVector(0, 0, 0) ; -// G4ThreeVector MMv = G4ThreeVector(0, 0, 0) ; + G4ThreeVector MMu = G4ThreeVector(0, 0, 0) ; + G4ThreeVector MMv = G4ThreeVector(0, 0, 0) ; G4ThreeVector MMw = G4ThreeVector(0, 0, 0) ; -// G4ThreeVector MMCenter = G4ThreeVector(0, 0, 0) ; + G4ThreeVector MMCenter = G4ThreeVector(0, 0, 0) ; bool FirstStage = true ; bool SecondStage = true ; bool ThirdStage = true ; @@ -555,12 +560,17 @@ void GaspardTrackerDummyShape::ConstructDetector(G4LogicalVolume* world) MMv = m_X1_Y128[i] - m_X1_Y1[i]; MMv = MMv.unit(); + + G4ThreeVector MMscal = MMu.dot(MMv); + MMw = MMu.cross(MMv); +// if (MMw.z() > 0) MMw = MMv.cross(MMu) ; MMw = MMw.unit(); MMCenter = (m_X1_Y1[i] + m_X1_Y128[i] + m_X128_Y1[i] + m_X128_Y128[i]) / 4; // Passage Matrix from Lab Referential to Telescope Referential + // MUST2 MMrot = new G4RotationMatrix(MMu, MMv, MMw); // translation to place Telescope MMpos = MMw * Length * 0.5 + MMCenter; @@ -582,7 +592,7 @@ void GaspardTrackerDummyShape::ConstructDetector(G4LogicalVolume* world) MMw = G4ThreeVector(wX, wY, wZ); // vector corresponding to the center of the module - MMCenter = MMw; + G4ThreeVector CT = MMw; // vector parallel to one axis of silicon plane G4double ii = cos(Theta / rad) * cos(Phi / rad); @@ -597,20 +607,21 @@ void GaspardTrackerDummyShape::ConstructDetector(G4LogicalVolume* world) MMu = MMu.unit(); // Passage Matrix from Lab Referential to Telescope Referential + // MUST2 MMrot = new G4RotationMatrix(MMu, MMv, MMw); // Telescope is rotate of Beta angle around MMv axis. MMrot->rotate(m_beta_u[i], MMu); MMrot->rotate(m_beta_v[i], MMv); MMrot->rotate(m_beta_w[i], MMw); // translation to place Telescope - MMpos = MMw * Length * 0.5 + MMCenter; + MMpos = MMw * Length * 0.5 + CT ; } FirstStage = m_wFirstStage[i] ; SecondStage = m_wSecondStage[i] ; ThirdStage = m_wThirdStage[i] ; - VolumeMaker(i + 1, MMpos, MMrot, FirstStage, SecondStage, ThirdStage, world); + VolumeMaker(i + 1, MMpos, MMrot, FirstStage, SecondStage, ThirdStage , world); } delete MMrot ; diff --git a/NPSimulation/src/GaspardTrackerSquare.cc b/NPSimulation/src/GaspardTrackerSquare.cc index 99f3974a02a9bd8ceeb260967feb075a00dd113c..71e021a8e7654bc974534bc8036ef05129d243cb 100644 --- a/NPSimulation/src/GaspardTrackerSquare.cc +++ b/NPSimulation/src/GaspardTrackerSquare.cc @@ -82,12 +82,12 @@ GaspardTrackerSquare::~GaspardTrackerSquare() //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void GaspardTrackerSquare::AddModule(G4ThreeVector X1_Y1 , - G4ThreeVector X128_Y1 , - G4ThreeVector X1_Y128 , - G4ThreeVector X128_Y128 , - bool wFirstStage , - bool wSecondStage , - bool wThirdStage) + G4ThreeVector X128_Y1 , + G4ThreeVector X1_Y128 , + G4ThreeVector X128_Y128 , + bool wFirstStage , + bool wSecondStage , + bool wThirdStage) { m_DefinitionType.push_back(true) ; @@ -111,14 +111,14 @@ void GaspardTrackerSquare::AddModule(G4ThreeVector X1_Y1 , //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void GaspardTrackerSquare::AddModule(G4double R , - G4double Theta , - G4double Phi , - G4double beta_u , - G4double beta_v , - G4double beta_w , - bool wFirstStage , - bool wSecondStage , - bool wThirdStage) + G4double Theta , + G4double Phi , + G4double beta_u , + G4double beta_v , + G4double beta_w , + bool wFirstStage , + bool wSecondStage , + bool wThirdStage) { G4ThreeVector empty = G4ThreeVector(0, 0, 0); @@ -188,7 +188,7 @@ void GaspardTrackerSquare::VolumeMaker(G4int TelescopeNumber, // Al density = 2.702 * g / cm3; a = 26.98 * g / mole; -// G4Material* Aluminium = new G4Material("Aluminium", z = 13., a, density); + G4Material* Aluminium = new G4Material("Aluminium", z = 13., a, density); // Iron // density = 7.874 * g / cm3; @@ -228,15 +228,37 @@ void GaspardTrackerSquare::VolumeMaker(G4int TelescopeNumber, // Little trick to avoid warning in compilation: Use a PVPlacement "buffer". // If don't you will have a Warning unused variable 'myPVP' G4PVPlacement* PVPBuffer ; - G4String Name = "GPDSquare" + DetectorNumber; - G4Box* solidGPDSquare = new G4Box(Name, 0.5*FaceFront, 0.5*FaceFront, 0.5*Length); - G4LogicalVolume* logicGPDSquare = new G4LogicalVolume(solidGPDSquare, Vacuum, Name, 0, 0, 0); + G4Trd* solidMM = new G4Trd("GPDSquare" + DetectorNumber, 0.5*FaceFront, 0.5*FaceFront, 0.5*FaceFront, 0.5*FaceFront, 0.5*Length); +// G4LogicalVolume* logicMM = new G4LogicalVolume(solidMM, Iron, "GPDSquare" + DetectorNumber, 0, 0, 0) ; + G4LogicalVolume* logicMM = new G4LogicalVolume(solidMM, Vacuum, "GPDSquare" + DetectorNumber, 0, 0, 0) ; - PVPBuffer = new G4PVPlacement(G4Transform3D(*MMrot, MMpos), logicGPDSquare, Name, world, false, 0); + G4String Name = "GPDSquare" + DetectorNumber ; + PVPBuffer = new G4PVPlacement(G4Transform3D(*MMrot, MMpos) , + logicMM , + Name , + world , + false , + 0); - logicGPDSquare->SetVisAttributes(G4VisAttributes::Invisible); - if (m_non_sensitive_part_visiualisation) logicGPDSquare->SetVisAttributes(G4VisAttributes(G4Colour(0.90, 0.90, 0.90))); + logicMM->SetVisAttributes(G4VisAttributes::Invisible); + if (m_non_sensitive_part_visiualisation) logicMM->SetVisAttributes(G4VisAttributes(G4Colour(0.90, 0.90, 0.90))); + + G4ThreeVector positionVacBox = G4ThreeVector(0, 0, VacBox_PosZ); + + G4Trd* solidVacBox = new G4Trd("solidVacBox", 0.5*SiliconFace, 0.5*SiliconFace, 0.5*SiliconFace, 0.5*SiliconFace, 0.5*VacBoxThickness); + G4LogicalVolume* logicVacBox = new G4LogicalVolume(solidVacBox, Vacuum, "logicVacBox", 0, 0, 0); + + PVPBuffer = new G4PVPlacement(0, positionVacBox, logicVacBox, "G" + DetectorNumber + "VacBox", logicMM, false, 0); + + logicVacBox->SetVisAttributes(G4VisAttributes::Invisible); + + // Add a degrader plate between Si and CsI: + /* + G4Box* Degrader = new G4Box("Degrader" , 50*mm , 50*mm , 0.1*mm ); + G4LogicalVolume* logicDegrader = new G4LogicalVolume( Degrader , Harvar, "logicDegrader",0,0,0); + PVPBuffer = new G4PVPlacement(0,G4ThreeVector(0,0,0),logicDegrader,"Degrader",logicVacBox,false,0) ; + */ //Place two marker to identify the u and v axis on silicon face: //marker are placed a bit before the silicon itself so they don't perturbate simulation @@ -262,55 +284,166 @@ void GaspardTrackerSquare::VolumeMaker(G4int TelescopeNumber, */ //////////////////////////////////////////////////////////////// - //////////////// First Stage Construction ////////////////////// + /////////////////Si Strip Construction////////////////////////// //////////////////////////////////////////////////////////////// if (wFirstStage) { + // Aluminium dead layers + G4ThreeVector positionAluStripFront = G4ThreeVector(0, 0, AluStripFront_PosZ); + G4ThreeVector positionAluStripBack = G4ThreeVector(0, 0, AluStripBack_PosZ); + + G4Box* solidAluStrip = new G4Box("AluBox", 0.5*SiliconFace, 0.5*SiliconFace, 0.5*AluStripThickness); +// G4LogicalVolume* logicAluStrip = new G4LogicalVolume(solidAluStrip, Aluminium, "logicAluStrip", 0, 0, 0); + G4LogicalVolume* logicAluStrip = new G4LogicalVolume(solidAluStrip, Vacuum, "logicAluStrip", 0, 0, 0); + + PVPBuffer = new G4PVPlacement(0, positionAluStripFront, logicAluStrip, "G" + DetectorNumber + "AluStripFront", logicMM, false, 0); + PVPBuffer = new G4PVPlacement(0, positionAluStripBack, logicAluStrip, "G" + DetectorNumber + "AluStripBack", logicMM, false, 0); + + logicAluStrip->SetVisAttributes(G4VisAttributes::Invisible); + // Silicon detector itself - G4ThreeVector positionFirstStage = G4ThreeVector(0, 0, FirstStage_PosZ); + G4ThreeVector positionSilicon = G4ThreeVector(0, 0, Silicon_PosZ); + + G4Box* solidSilicon = new G4Box("solidSilicon", 0.5*SiliconFace, 0.5*SiliconFace, 0.5*SiliconThickness); + G4LogicalVolume* logicSilicon = new G4LogicalVolume(solidSilicon, Silicon, "logicSilicon", 0, 0, 0); - G4Box* solidFirstStage = new G4Box("solidFirstStage", 0.5*FirstStageFace, 0.5*FirstStageFace, 0.5*FirstStageThickness); - G4LogicalVolume* logicFirstStage = new G4LogicalVolume(solidFirstStage, Silicon, "logicFirstStage", 0, 0, 0); + PVPBuffer = new G4PVPlacement(0, positionSilicon, logicSilicon, Name + "_Silicon", logicMM, false, 0); - PVPBuffer = new G4PVPlacement(0, - positionFirstStage, - logicFirstStage, - Name + "_FirstStage", - logicGPDSquare, - false, - 0); // Set First Stage sensible - logicFirstStage->SetSensitiveDetector(m_FirstStageScorer); + logicSilicon->SetSensitiveDetector(m_FirstStageScorer); - ///Visualisation of FirstStage Strip - G4VisAttributes* FirstStageVisAtt = new G4VisAttributes(G4Colour(0.0, 0.0, 0.9)); // blue - logicFirstStage->SetVisAttributes(FirstStageVisAtt); + ///Visualisation of Silicon Strip + G4VisAttributes* SiliconVisAtt = new G4VisAttributes(G4Colour(0.5, 0.5, 0.5)) ; + logicSilicon->SetVisAttributes(SiliconVisAtt) ; } //////////////////////////////////////////////////////////////// - //////////////////// Second Stage Construction //////////////// + //////////////////// SiLi Construction //////////////////////// //////////////////////////////////////////////////////////////// if (wSecondStage) { - // Second stage silicon detector - G4ThreeVector positionSecondStage = G4ThreeVector(0, 0, SecondStage_PosZ); - - G4Box* solidSecondStage = new G4Box("solidSecondStage", 0.5*SecondStageFace, 0.5*SecondStageFace, 0.5*SecondStageThickness); - G4LogicalVolume* logicSecondStage = new G4LogicalVolume(solidSecondStage, Silicon, "logicSecondStage", 0, 0, 0); - - PVPBuffer = new G4PVPlacement(0, - positionSecondStage, - logicSecondStage, - Name + "_SecondStage", - logicGPDSquare, - false, - 0); + G4double SiLiSpace = 8 * mm; + + G4RotationMatrix* rotSiLi = Rotation(0., 0., 0.); + + G4ThreeVector positionSiLi = G4ThreeVector(-0.25 * SiliconFace - 0.5 * SiLiSpace, 0, 0); + G4ThreeVector positionSiLi2 = G4ThreeVector(0.25 * SiliconFace + 0.5 * SiLiSpace, 0, 0); + + G4Box* solidSiLi = new G4Box("SiLi", 0.5*SiLiFaceX, 0.5*SiLiFaceY, 0.5*SiLiThickness); + + G4LogicalVolume* logicSiLi = new G4LogicalVolume(solidSiLi, Aluminium, "SiLi" + DetectorNumber, 0, 0, 0); + + // First Si(Li) 2 time 4 detectore + PVPBuffer = new G4PVPlacement(G4Transform3D(*rotSiLi, positionSiLi) , + logicSiLi , + "SiLi" + DetectorNumber , + logicVacBox , + false , + 0); + + // Second Si(Li) 2 time 4 detectore + PVPBuffer = new G4PVPlacement(G4Transform3D(*rotSiLi, positionSiLi2) , + logicSiLi , + "SiLi" + DetectorNumber , + logicVacBox , + false , + 1); + + // SiLi are placed inside of the VacBox... + // Left/Right define when looking to detector from Si to CsI + G4double SiLi_HighY_Upper = 19.86 * mm; + G4double SiLi_HighY_Center = 25.39 * mm ; + G4double SiLi_WidthX_Left = 22.85 * mm; + G4double SiLi_WidthX_Right = 24.9 * mm ; + G4double SiLi_ShiftX = 0.775 * mm ; + + // SiLi : left side of SiLi detector + G4Box* solidSiLi_LT = new G4Box("SiLi_LT" , 0.5*SiLi_WidthX_Left , 0.5*SiLi_HighY_Upper , 0.5*SiLiThickness); + G4Box* solidSiLi_RT = new G4Box("SiLi_RT" , 0.5*SiLi_WidthX_Right , 0.5*SiLi_HighY_Upper , 0.5*SiLiThickness); + G4Box* solidSiLi_LC1 = new G4Box("SiLi_LC1" , 0.5*SiLi_WidthX_Left , 0.5*SiLi_HighY_Center , 0.5*SiLiThickness); + G4Box* solidSiLi_RC1 = new G4Box("SiLi_RC1" , 0.5*SiLi_WidthX_Right , 0.5*SiLi_HighY_Center , 0.5*SiLiThickness); + G4Box* solidSiLi_LB = new G4Box("SiLi_LB" , 0.5*SiLi_WidthX_Left , 0.5*SiLi_HighY_Upper , 0.5*SiLiThickness); + G4Box* solidSiLi_RB = new G4Box("SiLi_RB" , 0.5*SiLi_WidthX_Right , 0.5*SiLi_HighY_Upper , 0.5*SiLiThickness); + G4Box* solidSiLi_LC2 = new G4Box("SiLi_LC2" , 0.5*SiLi_WidthX_Left , 0.5*SiLi_HighY_Center , 0.5*SiLiThickness); + G4Box* solidSiLi_RC2 = new G4Box("SiLi_RC2" , 0.5*SiLi_WidthX_Right , 0.5*SiLi_HighY_Center , 0.5*SiLiThickness); + + G4LogicalVolume* logicSiLi_LT = new G4LogicalVolume(solidSiLi_LT , Silicon , "SiLi_LT" , 0 , 0 , 0); + G4LogicalVolume* logicSiLi_RT = new G4LogicalVolume(solidSiLi_RT , Silicon , "SiLi_RT" , 0 , 0 , 0); + G4LogicalVolume* logicSiLi_LC1 = new G4LogicalVolume(solidSiLi_LC1 , Silicon , "SiLi_LC1" , 0 , 0 , 0); + G4LogicalVolume* logicSiLi_RC1 = new G4LogicalVolume(solidSiLi_RC1 , Silicon , "SiLi_RC1" , 0 , 0 , 0); + G4LogicalVolume* logicSiLi_LB = new G4LogicalVolume(solidSiLi_LB , Silicon , "SiLi_LB" , 0 , 0 , 0); + G4LogicalVolume* logicSiLi_RB = new G4LogicalVolume(solidSiLi_RB , Silicon , "SiLi_RB" , 0 , 0 , 0); + G4LogicalVolume* logicSiLi_LC2 = new G4LogicalVolume(solidSiLi_LC2 , Silicon , "SiLi_LC2" , 0 , 0 , 0); + G4LogicalVolume* logicSiLi_RC2 = new G4LogicalVolume(solidSiLi_RC2 , Silicon , "SiLi_RC2" , 0 , 0 , 0); + + G4double interSiLi = 0.5 * mm; + + // Top + G4ThreeVector positionSiLi_LT = G4ThreeVector(-0.5 * SiLi_WidthX_Left - interSiLi - SiLi_ShiftX , + 0.5 * SiLi_HighY_Upper + SiLi_HighY_Center + 1.5 * interSiLi , + 0); + + G4ThreeVector positionSiLi_RT = G4ThreeVector(0.5 * SiLi_WidthX_Right - SiLi_ShiftX , + 0.5 * SiLi_HighY_Upper + SiLi_HighY_Center + 1.5 * interSiLi , + 0); + + G4ThreeVector positionSiLi_LC1 = G4ThreeVector(-0.5 * SiLi_WidthX_Left - interSiLi - SiLi_ShiftX , + 0.5 * SiLi_HighY_Center + 0.5 * interSiLi , + 0); + + G4ThreeVector positionSiLi_RC1 = G4ThreeVector(0.5 * SiLi_WidthX_Right - SiLi_ShiftX , + 0.5 * SiLi_HighY_Center + 0.5 * interSiLi , + 0); + + // Bottom + G4ThreeVector positionSiLi_LB = G4ThreeVector(-0.5 * SiLi_WidthX_Left - interSiLi - SiLi_ShiftX , + -0.5 * SiLi_HighY_Upper - SiLi_HighY_Center - 1.5 * interSiLi , + 0); + + G4ThreeVector positionSiLi_RB = G4ThreeVector(0.5 * SiLi_WidthX_Right - SiLi_ShiftX , + -0.5 * SiLi_HighY_Upper - SiLi_HighY_Center - 1.5 * interSiLi , + 0); + + G4ThreeVector positionSiLi_LC2 = G4ThreeVector(-0.5 * SiLi_WidthX_Left - interSiLi - SiLi_ShiftX , + -0.5 * SiLi_HighY_Center - 0.5 * interSiLi , + 0); + + G4ThreeVector positionSiLi_RC2 = G4ThreeVector(0.5 * SiLi_WidthX_Right - SiLi_ShiftX , + -0.5 * SiLi_HighY_Center - 0.5 * interSiLi , + 0); + + PVPBuffer = new G4PVPlacement(0 , positionSiLi_LT , logicSiLi_LT , Name + "_SiLi_LT" , logicSiLi , false , 0) ; + PVPBuffer = new G4PVPlacement(0 , positionSiLi_RT , logicSiLi_RT , Name + "_SiLi_RT" , logicSiLi , false , 0) ; + PVPBuffer = new G4PVPlacement(0 , positionSiLi_LC1 , logicSiLi_LC1 , Name + "_SiLi_LC1" , logicSiLi , false , 0) ; + PVPBuffer = new G4PVPlacement(0 , positionSiLi_RC1 , logicSiLi_RC1 , Name + "_SiLi_RC1" , logicSiLi , false , 0) ; + + PVPBuffer = new G4PVPlacement(0 , positionSiLi_LB , logicSiLi_LB , Name + "_SiLi_LB" , logicSiLi , false , 0) ; + PVPBuffer = new G4PVPlacement(0 , positionSiLi_RB , logicSiLi_RB , Name + "_SiLi_RB" , logicSiLi , false , 0) ; + PVPBuffer = new G4PVPlacement(0 , positionSiLi_LC2 , logicSiLi_LC2 , Name + "_SiLi_LC2" , logicSiLi , false , 0) ; + PVPBuffer = new G4PVPlacement(0 , positionSiLi_RC2 , logicSiLi_RC2 , Name + "_SiLi_RC2" , logicSiLi , false , 0) ; + + logicSiLi->SetVisAttributes(G4VisAttributes(G4Colour(1, 1., 1.))); // Set Second Stage sensible - logicSecondStage->SetSensitiveDetector(m_SecondStageScorer); - - ///Visualisation of SecondStage Strip - G4VisAttributes* SecondStageVisAtt = new G4VisAttributes(G4Colour(0.5, 0.5, 0.5)) ; - logicSecondStage->SetVisAttributes(SecondStageVisAtt) ; + logicSiLi_LT->SetSensitiveDetector(m_SecondStageScorer); + logicSiLi_RT->SetSensitiveDetector(m_SecondStageScorer); + logicSiLi_LC1->SetSensitiveDetector(m_SecondStageScorer); + logicSiLi_RC1->SetSensitiveDetector(m_SecondStageScorer); + + logicSiLi_LB->SetSensitiveDetector(m_SecondStageScorer); + logicSiLi_RB->SetSensitiveDetector(m_SecondStageScorer); + logicSiLi_LC2->SetSensitiveDetector(m_SecondStageScorer); + logicSiLi_RC2->SetSensitiveDetector(m_SecondStageScorer); + + // Mark blue a SiLi to see telescope orientation + logicSiLi_LT->SetVisAttributes(G4VisAttributes(G4Colour(0, 0., 1.))); + logicSiLi_RT->SetVisAttributes(G4VisAttributes(G4Colour(0, 1., 0))); + logicSiLi_LC1->SetVisAttributes(G4VisAttributes(G4Colour(0, 1., 0))); + logicSiLi_RC1->SetVisAttributes(G4VisAttributes(G4Colour(0, 1., 0))); + + logicSiLi_LB->SetVisAttributes(G4VisAttributes(G4Colour(0, 1., 0))); + logicSiLi_RB->SetVisAttributes(G4VisAttributes(G4Colour(0, 1., 0))); + logicSiLi_LC2->SetVisAttributes(G4VisAttributes(G4Colour(0, 1., 0))); + logicSiLi_RC2->SetVisAttributes(G4VisAttributes(G4Colour(0, 1., 0))); } //////////////////////////////////////////////////////////////// @@ -320,23 +453,19 @@ void GaspardTrackerSquare::VolumeMaker(G4int TelescopeNumber, // Third stage silicon detector G4ThreeVector positionThirdStage = G4ThreeVector(0, 0, ThirdStage_PosZ); - G4Box* solidThirdStage = new G4Box("solidThirdStage", 0.5*ThirdStageFace, 0.5*ThirdStageFace, 0.5*ThirdStageThickness); +// G4Box* solidThirdStage = new G4Box("solidThirdStage", 0.5*SiliconFace, 0.5*SiliconFace, 0.5*ThirdStageThickness); + G4Box* solidThirdStage = new G4Box("solidThirdStage", 0.5*FaceFront, 0.5*FaceFront, 0.5*ThirdStageThickness); G4LogicalVolume* logicThirdStage = new G4LogicalVolume(solidThirdStage, Silicon, "logicThirdStage", 0, 0, 0); - PVPBuffer = new G4PVPlacement(0, - positionThirdStage, - logicThirdStage, - Name + "_ThirdStage", - logicGPDSquare, - false, - 0); + PVPBuffer = new G4PVPlacement(0, positionThirdStage, logicThirdStage, Name + "_ThirdStage", logicMM, false, 0); + + ///Visualisation of Third Stage + G4VisAttributes* ThirdStageVisAtt = new G4VisAttributes(G4Colour(0.7, 0.7, 0.7)) ; + logicThirdStage->SetVisAttributes(ThirdStageVisAtt) ; +// logicThirdStage->SetVisAttributes(G4VisAttributes::Invisible); // Set Third Stage sensible logicThirdStage->SetSensitiveDetector(m_ThirdStageScorer); - - ///Visualisation of Third Stage - G4VisAttributes* ThirdStageVisAtt = new G4VisAttributes(G4Colour(0.0, 0.9, 0.0)); // green - logicThirdStage->SetVisAttributes(ThirdStageVisAtt); } } @@ -582,14 +711,18 @@ void GaspardTrackerSquare::ReadConfiguration(string Path) void GaspardTrackerSquare::ConstructDetector(G4LogicalVolume* world) { G4RotationMatrix* MMrot = NULL ; - G4ThreeVector MMpos = G4ThreeVector(0, 0, 0) ; +/* G4ThreeVector MMpos = G4ThreeVector(0, 0, 0) ; G4ThreeVector MMu = G4ThreeVector(0, 0, 0) ; G4ThreeVector MMv = G4ThreeVector(0, 0, 0) ; - G4ThreeVector MMw = G4ThreeVector(0, 0, 0) ; + G4ThreeVector MMw = G4ThreeVector(0, 0, 0) ;*/ + MMpos = G4ThreeVector(0, 0, 0) ; + MMu = G4ThreeVector(0, 0, 0) ; + MMv = G4ThreeVector(0, 0, 0) ; + MMw = G4ThreeVector(0, 0, 0) ; G4ThreeVector MMCenter = G4ThreeVector(0, 0, 0) ; - bool FirstStage = true; - bool SecondStage = true; - bool ThirdStage = true; + bool FirstStage = true ; + bool SecondStage = true ; + bool ThirdStage = true ; G4int NumberOfTelescope = m_DefinitionType.size() ; @@ -599,21 +732,31 @@ void GaspardTrackerSquare::ConstructDetector(G4LogicalVolume* world) // (u,v,w) unitary vector associated to telescope referencial // (u,v) // to silicon plan // w perpendicular to (u,v) plan and pointing ThirdStage - MMu = m_X128_Y1[i] - m_X1_Y1[i]; - MMu = MMu.unit(); + G4cout << "############ Gaspard " << i << " #############" << G4endl; + MMu = m_X128_Y1[i] - m_X1_Y1[i] ; + G4cout << "MMu: X = " << MMu(0) << " , Y = " << MMu(1) << " , Z = " << MMu(2) << G4endl; + MMu = MMu.unit() ; + G4cout << "Norm MMu: X = " << MMu(0) << " , Y = " << MMu(1) << " , Z = " << MMu(2) << G4endl; - MMv = m_X1_Y128[i] - m_X1_Y1[i]; - MMv = MMv.unit(); + MMv = m_X1_Y128[i] - m_X1_Y1[i] ; + G4cout << "MMv X = " << MMv(0) << " , Y = " << MMv(1) << " , Z = " << MMv(2) << G4endl; + MMv = MMv.unit() ; + G4cout << "Norm MMv X = " << MMv(0) << " , Y = " << MMv(1) << " , Z = " << MMv(2) << G4endl; - MMw = MMu.cross(MMv); - MMw = MMw.unit(); + G4ThreeVector MMscal = MMu.dot(MMv); + G4cout << "Norm MMu.MMv X = " << MMv(0) << " , Y = " << MMv(1) << " , Z = " << MMv(2) << G4endl; - MMCenter = (m_X1_Y1[i] + m_X1_Y128[i] + m_X128_Y1[i] + m_X128_Y128[i]) / 4; + MMw = MMu.cross(MMv) ; +// if (MMw.z() > 0) MMw = MMv.cross(MMu) ; + MMw = MMw.unit() ; + + MMCenter = (m_X1_Y1[i] + m_X1_Y128[i] + m_X128_Y1[i] + m_X128_Y128[i]) / 4 ; // Passage Matrix from Lab Referential to Telescope Referential - MMrot = new G4RotationMatrix(MMu, MMv, MMw); + // MUST2 + MMrot = new G4RotationMatrix(MMu, MMv, MMw) ; // translation to place Telescope - MMpos = MMw * Length * 0.5 + MMCenter; + MMpos = MMw * Length * 0.5 + MMCenter ; } // By Angle @@ -626,13 +769,13 @@ void GaspardTrackerSquare::ConstructDetector(G4LogicalVolume* world) // w perpendicular to (u,v) plan and pointing ThirdStage // Phi is angle between X axis and projection in (X,Y) plan // Theta is angle between position vector and z axis - G4double wX = m_R[i] * sin(Theta / rad) * cos(Phi / rad); - G4double wY = m_R[i] * sin(Theta / rad) * sin(Phi / rad); - G4double wZ = m_R[i] * cos(Theta / rad); - MMw = G4ThreeVector(wX, wY, wZ); + G4double wX = m_R[i] * sin(Theta / rad) * cos(Phi / rad) ; + G4double wY = m_R[i] * sin(Theta / rad) * sin(Phi / rad) ; + G4double wZ = m_R[i] * cos(Theta / rad) ; + MMw = G4ThreeVector(wX, wY, wZ) ; // vector corresponding to the center of the module - MMCenter = MMw; + CT = MMw; // vector parallel to one axis of silicon plane G4double ii = cos(Theta / rad) * cos(Phi / rad); @@ -654,7 +797,7 @@ void GaspardTrackerSquare::ConstructDetector(G4LogicalVolume* world) MMrot->rotate(m_beta_v[i], MMv); MMrot->rotate(m_beta_w[i], MMw); // translation to place Telescope - MMpos = MMw * Length * 0.5 + MMCenter; + MMpos = MMw * Length * 0.5 + CT ; } FirstStage = m_wFirstStage[i] ; @@ -717,8 +860,8 @@ void GaspardTrackerSquare::ReadSensitive(const G4Event* event) // NULL pointer are given to avoid warning at compilation // Si(Li) - std::map<G4int, G4double*>::iterator SecondStageEnergy_itr ; - G4THitsMap<G4double>* SecondStageEnergyHitMap = NULL ; + std::map<G4int, G4double*>::iterator SiLiEnergy_itr ; + G4THitsMap<G4double>* SiLiEnergyHitMap = NULL ; // Third Stage std::map<G4int, G4double*>::iterator ThirdStageEnergy_itr ; G4THitsMap<G4double>* ThirdStageEnergyHitMap = NULL ; @@ -776,14 +919,14 @@ void GaspardTrackerSquare::ReadSensitive(const G4Event* event) Ang_Phi_itr = AngPhiHitMap->GetMap()->begin() ; - // Read the Scorer associate to the SecondStage + // Read the Scorer associate to the SiLi //Energy - G4int SecondStageEnergyCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("SecondStageScorerGPDSquare/SecondStageEnergy") ; - SecondStageEnergyHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(SecondStageEnergyCollectionID)) ; - SecondStageEnergy_itr = SecondStageEnergyHitMap->GetMap()->begin() ; + G4int SiLiEnergyCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("SecondStageScorerGPDSquare/SecondStageEnergy") ; + SiLiEnergyHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(SiLiEnergyCollectionID)) ; + SiLiEnergy_itr = SiLiEnergyHitMap->GetMap()->begin() ; - // Read the Scorer associate to the ThirdStage + // Read the Scorer associate to the CsI crystal //Energy G4int ThirdStageEnergyCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThirdStageScorerGPDSquare/ThirdStageEnergy"); ThirdStageEnergyHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(ThirdStageEnergyCollectionID)); @@ -920,13 +1063,13 @@ void GaspardTrackerSquare::ReadSensitive(const G4Event* event) } // Second Stage - SecondStageEnergy_itr = SecondStageEnergyHitMap->GetMap()->begin() ; - for (G4int h = 0 ; h < SecondStageEnergyHitMap->entries() ; h++) { - G4int SecondStageEnergyTrackID = SecondStageEnergy_itr->first - N; - G4double SecondStageEnergy = *(SecondStageEnergy_itr->second); + SiLiEnergy_itr = SiLiEnergyHitMap->GetMap()->begin() ; + for (G4int h = 0 ; h < SiLiEnergyHitMap->entries() ; h++) { + G4int SiLiEnergyTrackID = SiLiEnergy_itr->first - N; + G4double SiLiEnergy = *(SiLiEnergy_itr->second); - if (SecondStageEnergyTrackID == NTrackID) { - ms_Event->SetGPDTrkSecondStageEEnergy(RandGauss::shoot(SecondStageEnergy, ResoSecondStage)) ; + if (SiLiEnergyTrackID == NTrackID) { + ms_Event->SetGPDTrkSecondStageEEnergy(RandGauss::shoot(SiLiEnergy, ResoSecondStage)) ; ms_Event->SetGPDTrkSecondStageEPadNbr(1); ms_Event->SetGPDTrkSecondStageTPadNbr(1); ms_Event->SetGPDTrkSecondStageTTime(1); @@ -934,7 +1077,7 @@ void GaspardTrackerSquare::ReadSensitive(const G4Event* event) ms_Event->SetGPDTrkSecondStageEDetectorNbr(m_index["Square"] + N); } - SecondStageEnergy_itr++; + SiLiEnergy_itr++; } // Third Stage @@ -969,7 +1112,7 @@ void GaspardTrackerSquare::ReadSensitive(const G4Event* event) PosZHitMap ->clear(); AngThetaHitMap ->clear(); AngPhiHitMap ->clear(); - SecondStageEnergyHitMap ->clear() ; + SiLiEnergyHitMap ->clear() ; ThirdStageEnergyHitMap ->clear() ; } } diff --git a/NPSimulation/src/GaspardTrackerTrapezoid.cc b/NPSimulation/src/GaspardTrackerTrapezoid.cc index 259923c65a62c1d8a58a5cdc3f7199d89969f43c..1a501934aee355b17b4a2ea32127491e28cef858 100644 --- a/NPSimulation/src/GaspardTrackerTrapezoid.cc +++ b/NPSimulation/src/GaspardTrackerTrapezoid.cc @@ -224,19 +224,47 @@ void GaspardTrackerTrapezoid::VolumeMaker(G4int TelescopeNumber , // Little trick to avoid warning in compilation: Use a PVPlacement "buffer". // If don't you will have a Warning unused variable 'myPVP' G4PVPlacement* PVPBuffer ; - G4String Name = "GPDTrapezoid" + DetectorNumber ; // Definition of the volume containing the sensitive detector - G4Trap* solidGPDTrapezoid = new G4Trap(Name, - Length/2, 0*deg, 0*deg, - Height/2, BaseLarge/2, BaseSmall/2, 0*deg, - Height/2, BaseLarge/2, BaseSmall/2, 0*deg); - G4LogicalVolume* logicGPDTrapezoid = new G4LogicalVolume(solidGPDTrapezoid, Vacuum, Name, 0, 0, 0); + G4Trap* solidMM = new G4Trap("GPDTrapezoid" + DetectorNumber, + Length/2, 0*deg, 0*deg, + Height/2, BaseSmall/2, BaseLarge/2, 0*deg, + Height/2, BaseSmall/2, BaseLarge/2, 0*deg); + +// G4LogicalVolume* logicMM = new G4LogicalVolume(solidMM, Iron, "GPDTrapezoid" + DetectorNumber, 0, 0, 0) ; + G4LogicalVolume* logicMM = new G4LogicalVolume(solidMM, Vacuum, "GPDTrapezoid" + DetectorNumber, 0, 0, 0) ; + + G4String Name = "GPDTrapezoid" + DetectorNumber ; + PVPBuffer = new G4PVPlacement(G4Transform3D(*MMrot, MMpos) , + logicMM , + Name , + world , + false , + 0); + + logicMM->SetVisAttributes(G4VisAttributes::Invisible); + if (m_non_sensitive_part_visiualisation) logicMM->SetVisAttributes(G4VisAttributes(G4Colour(0.90, 0.90, 0.90))); + + // Definition of a vaccuum volume + G4ThreeVector positionVacBox = G4ThreeVector(0, 0, VacBox_PosZ); - PVPBuffer = new G4PVPlacement(G4Transform3D(*MMrot, MMpos), logicGPDTrapezoid, Name, world, false, 0); + G4Trap* solidVacBox = new G4Trap("solidVacBox", + VacBoxThickness/2, 0*deg, 0*deg, + FirstStageHeight/2, FirstStageBaseSmall/2, FirstStageBaseLarge/2, 0*deg, + FirstStageHeight/2, FirstStageBaseSmall/2, FirstStageBaseLarge/2, 0*deg); - logicGPDTrapezoid->SetVisAttributes(G4VisAttributes::Invisible); - if (m_non_sensitive_part_visiualisation) logicGPDTrapezoid->SetVisAttributes(G4VisAttributes(G4Colour(0.90, 0.90, 0.90))); + G4LogicalVolume* logicVacBox = new G4LogicalVolume(solidVacBox, Vacuum, "logicVacBox", 0, 0, 0); + + PVPBuffer = new G4PVPlacement(0, positionVacBox, logicVacBox, "G" + DetectorNumber + "VacBox", logicMM, false, 0); + + logicVacBox->SetVisAttributes(G4VisAttributes::Invisible); + + // Add a degrader plate between Si and CsI: + /* + G4Box* Degrader = new G4Box("Degrader" , 50*mm , 50*mm , 0.1*mm ); + G4LogicalVolume* logicDegrader = new G4LogicalVolume( Degrader , Harvar, "logicDegrader",0,0,0); + PVPBuffer = new G4PVPlacement(0,G4ThreeVector(0,0,0),logicDegrader,"Degrader",logicVacBox,false,0) ; + */ //Place two marker to identify the u and v axis on silicon face: //marker are placed a bit before the silicon itself so they don't perturbate simulation @@ -265,58 +293,46 @@ void GaspardTrackerTrapezoid::VolumeMaker(G4int TelescopeNumber , /////////////////// First Stage Construction//////////////////// //////////////////////////////////////////////////////////////// if (wFirstStage) { + // Aluminium dead layers + G4ThreeVector positionAluStripFront = G4ThreeVector(0, 0, AluStripFront_PosZ); + G4ThreeVector positionAluStripBack = G4ThreeVector(0, 0, AluStripBack_PosZ); + + G4Trap* solidAluStrip = new G4Trap("AluBox", + AluStripThickness/2, 0*deg, 0*deg, + FirstStageHeight/2, FirstStageBaseSmall/2, FirstStageBaseLarge/2, 0*deg, + FirstStageHeight/2, FirstStageBaseSmall/2, FirstStageBaseLarge/2, 0*deg); + +// G4LogicalVolume* logicAluStrip = new G4LogicalVolume(solidAluStrip, Aluminium, "logicAluStrip", 0, 0, 0); + G4LogicalVolume* logicAluStrip = new G4LogicalVolume(solidAluStrip, Vacuum, "logicAluStrip", 0, 0, 0); + + PVPBuffer = new G4PVPlacement(0, positionAluStripFront, logicAluStrip, "G" + DetectorNumber + "AluStripFront", logicMM, false, 0); + PVPBuffer = new G4PVPlacement(0, positionAluStripBack, logicAluStrip, "G" + DetectorNumber + "AluStripBack", logicMM, false, 0); + + logicAluStrip->SetVisAttributes(G4VisAttributes::Invisible); + // Silicon detector itself - G4ThreeVector positionFirstStage = G4ThreeVector(0, 0, FirstStage_PosZ); - - G4Trap* solidFirstStage = new G4Trap("solidFirstStage", - FirstStageThickness/2, 0*deg, 0*deg, - FirstStageHeight/2, FirstStageBaseLarge/2, FirstStageBaseSmall/2, 0*deg, - FirstStageHeight/2, FirstStageBaseLarge/2, FirstStageBaseSmall/2, 0*deg); - G4LogicalVolume* logicFirstStage = new G4LogicalVolume(solidFirstStage, Silicon, "logicFirstStage", 0, 0, 0); - - PVPBuffer = new G4PVPlacement(0, - positionFirstStage, - logicFirstStage, - Name + "_FirstStage", - logicGPDTrapezoid, - false, - 0); + G4ThreeVector positionSilicon = G4ThreeVector(0, 0, Silicon_PosZ); + + G4Trap* solidSilicon = new G4Trap("solidSilicon", + FirstStageThickness/2, 0*deg, 0*deg, + FirstStageHeight/2, FirstStageBaseSmall/2, FirstStageBaseLarge/2, 0*deg, + FirstStageHeight/2, FirstStageBaseSmall/2, FirstStageBaseLarge/2, 0*deg); + G4LogicalVolume* logicSilicon = new G4LogicalVolume(solidSilicon, Silicon, "logicSilicon", 0, 0, 0); + + PVPBuffer = new G4PVPlacement(0, positionSilicon, logicSilicon, Name + "_Silicon", logicMM, false, 0); // Set First Stage sensible - logicFirstStage->SetSensitiveDetector(m_FirstStageScorer); + logicSilicon->SetSensitiveDetector(m_FirstStageScorer); - ///Visualisation of FirstStage Strip - G4VisAttributes* FirstStageVisAtt = new G4VisAttributes(G4Colour(0.0, 0.0, 0.9)); // blue - logicFirstStage->SetVisAttributes(FirstStageVisAtt); + // Visualisation of Silicon Strip + G4VisAttributes* SiliconVisAtt = new G4VisAttributes(G4Colour(0.5, 0.5, 0.5)) ; + logicSilicon->SetVisAttributes(SiliconVisAtt) ; } //////////////////////////////////////////////////////////////// //////////////// Second Stage Construction //////////////////// //////////////////////////////////////////////////////////////// if (wSecondStage) { - // Second stage silicon detector - G4ThreeVector positionSecondStage = G4ThreeVector(0, 0, SecondStage_PosZ); - - G4Trap* solidSecondStage = new G4Trap("solidSecondStage", - SecondStageThickness/2, 0*deg, 0*deg, - FirstStageHeight/2, FirstStageBaseLarge/2, FirstStageBaseSmall/2, 0*deg, - FirstStageHeight/2, FirstStageBaseLarge/2, FirstStageBaseSmall/2, 0*deg); - G4LogicalVolume* logicSecondStage = new G4LogicalVolume(solidSecondStage, Silicon, "logicSecondStage", 0, 0, 0); - - PVPBuffer = new G4PVPlacement(0, - positionSecondStage, - logicSecondStage, - Name + "_SecondStage", - logicGPDTrapezoid, - false, - 0); - - // Set Second Stage sensible - logicSecondStage->SetSensitiveDetector(m_SecondStageScorer); - - ///Visualisation of SecondStage Strip - G4VisAttributes* SecondStageVisAtt = new G4VisAttributes(G4Colour(0.5, 0.5, 0.5)); - logicSecondStage->SetVisAttributes(SecondStageVisAtt); } //////////////////////////////////////////////////////////////// @@ -328,24 +344,20 @@ void GaspardTrackerTrapezoid::VolumeMaker(G4int TelescopeNumber , G4Trap* solidThirdStage = new G4Trap("solidThirdStage", ThirdStageThickness/2, 0*deg, 0*deg, - FirstStageHeight/2, FirstStageBaseLarge/2, FirstStageBaseSmall/2, 0*deg, - FirstStageHeight/2, FirstStageBaseLarge/2, FirstStageBaseSmall/2, 0*deg); + Height/2, BaseSmall/2, BaseLarge/2, 0*deg, + Height/2, BaseSmall/2, BaseLarge/2, 0*deg); + G4LogicalVolume* logicThirdStage = new G4LogicalVolume(solidThirdStage, Silicon, "logicThirdStage", 0, 0, 0); - PVPBuffer = new G4PVPlacement(0, - positionThirdStage, - logicThirdStage, - Name + "_ThirdStage", - logicGPDTrapezoid, - false, - 0); + PVPBuffer = new G4PVPlacement(0, positionThirdStage, logicThirdStage, Name + "_ThirdStage", logicMM, false, 0); + + // Visualisation of Third Stage + G4VisAttributes* ThirdStageVisAtt = new G4VisAttributes(G4Colour(0.7, 0.7, 0.7)) ; + logicThirdStage->SetVisAttributes(ThirdStageVisAtt) ; +// logicThirdStage->SetVisAttributes(G4VisAttributes::Invisible); // Set Third Stage sensible logicThirdStage->SetSensitiveDetector(m_ThirdStageScorer); - - ///Visualisation of Third Stage - G4VisAttributes* ThirdStageVisAtt = new G4VisAttributes(G4Colour(0.0, 0.9, 0.0)); // red - logicThirdStage->SetVisAttributes(ThirdStageVisAtt); } } @@ -590,13 +602,16 @@ void GaspardTrackerTrapezoid::ReadConfiguration(string Path) // Called After DetecorConstruction::AddDetector Method void GaspardTrackerTrapezoid::ConstructDetector(G4LogicalVolume* world) { - G4RotationMatrix* MMrot = NULL; - G4ThreeVector MMpos = G4ThreeVector(0, 0, 0); - G4ThreeVector MMu = G4ThreeVector(0, 0, 0); - G4ThreeVector MMv = G4ThreeVector(0, 0, 0); - G4ThreeVector MMw = G4ThreeVector(0, 0, 0); - G4ThreeVector MMCenter = G4ThreeVector(0, 0, 0); - + G4RotationMatrix* MMrot = NULL ; +/* G4ThreeVector MMpos = G4ThreeVector(0, 0, 0) ; + G4ThreeVector MMu = G4ThreeVector(0, 0, 0) ; + G4ThreeVector MMv = G4ThreeVector(0, 0, 0) ; + G4ThreeVector MMw = G4ThreeVector(0, 0, 0) ;*/ + MMpos = G4ThreeVector(0, 0, 0) ; + MMu = G4ThreeVector(0, 0, 0) ; + MMv = G4ThreeVector(0, 0, 0) ; + MMw = G4ThreeVector(0, 0, 0) ; + G4ThreeVector MMCenter = G4ThreeVector(0, 0, 0) ; bool FirstStage = true ; bool SecondStage = true ; bool ThirdStage = true ; @@ -606,97 +621,86 @@ void GaspardTrackerTrapezoid::ConstructDetector(G4LogicalVolume* world) for (G4int i = 0; i < NumberOfModule; i++) { // By Point if (m_DefinitionType[i]) { - // (u,v,w) unitary vector associated to trapezoidal referencial + // (u,v,w) unitary vector associated to telescope referencial // (u,v) // to silicon plan - // ------- - // / \ ^ - // / \ | v - // / \ | - // --------------- <------ - // u // w perpendicular to (u,v) plan and pointing ThirdStage - G4cout << "XXXXXXXXXXXX Trapezoid " << i << " XXXXXXXXXXXXX" << G4endl; - MMu = m_X128_Y1[i] - m_X1_Y1[i]; - MMu = MMu.unit(); + G4cout << "############ Gaspard Trapezoid " << i << " #############" << G4endl; + MMu = m_X128_Y1[i] - m_X1_Y128[i] ; G4cout << "MMu: " << MMu << G4endl; + MMu = MMu.unit() ; + G4cout << "Norm MMu: " << MMu << G4endl; - MMv = 0.5 * (m_X1_Y128[i] + m_X128_Y128[i] - m_X1_Y1[i] - m_X128_Y1[i]); - MMv = MMv.unit(); + MMv = -0.5 * (m_X1_Y1[i] + m_X128_Y128[i] - m_X1_Y128[i] - m_X128_Y1[i]); G4cout << "MMv: " << MMv << G4endl; + MMv = MMv.unit() ; + G4cout << "Norm MMv: " << MMv << G4endl; - MMw = MMu.cross(MMv); - MMw = MMw.unit(); - G4cout << "MMw: " << MMw << G4endl; + G4ThreeVector MMscal = MMu.dot(MMv); + G4cout << "Norm MMu.MMv: " << MMscal << G4endl; + + MMw = MMu.cross(MMv) ; +// if (MMw.z() > 0) MMw = MMv.cross(MMu) ; + MMw = MMw.unit() ; + G4cout << "Norm MMw: " << MMw << G4endl; // Center of the module - MMCenter = (m_X1_Y1[i] + m_X1_Y128[i] + m_X128_Y1[i] + m_X128_Y128[i]) / 4; + MMCenter = (m_X1_Y1[i] + m_X1_Y128[i] + m_X128_Y1[i] + m_X128_Y128[i]) / 4 ; // Passage Matrix from Lab Referential to Module Referential - MMrot = new G4RotationMatrix(MMu, MMv, MMw); + MMrot = new G4RotationMatrix(MMu, MMv, MMw) ; // translation to place Module - MMpos = MMw * Length * 0.5 + MMCenter; + MMpos = MMw * Length * 0.5 + MMCenter ; } // By Angle else { - G4double Theta = m_Theta[i]; - G4double Phi = m_Phi[i]; + G4double Theta = m_Theta[i] ; + G4double Phi = m_Phi[i] ; + //This part because if Phi and Theta = 0 equation are false + if (Theta == 0) Theta = 0.0001 ; + if (Theta == 2*cos(0)) Theta = 2 * acos(0) - 0.00001 ; + if (Phi == 0) Phi = 0.0001 ; // (u,v,w) unitary vector associated to telescope referencial // (u,v) // to silicon plan - // ------- - // / \ ^ - // / \ | v - // / \ | - // --------------- <------ - // u // w perpendicular to (u,v) plan and pointing ThirdStage // Phi is angle between X axis and projection in (X,Y) plan // Theta is angle between position vector and z axis - G4double wX = m_R[i] * sin(Theta / rad) * cos(Phi / rad); - G4double wY = m_R[i] * sin(Theta / rad) * sin(Phi / rad); - G4double wZ = m_R[i] * cos(Theta / rad); - MMw = G4ThreeVector(wX, wY, wZ); - - // vector corresponding to the center of the module - MMCenter = MMw; - - // vector parallel to one axis of silicon plane - // in fact, this is vector u - G4double ii = cos(Theta / rad) * cos(Phi / rad); - G4double jj = cos(Theta / rad) * sin(Phi / rad); - G4double kk = -sin(Theta / rad); - G4ThreeVector Y = G4ThreeVector(ii, jj, kk); - - MMw = MMw.unit(); - MMv = MMw.cross(Y); - MMu = MMv.cross(MMw); - MMv = MMv.unit(); - MMu = MMu.unit(); + G4double wX = m_R[i] * sin(Theta / rad) * cos(Phi / rad) ; + G4double wY = m_R[i] * sin(Theta / rad) * sin(Phi / rad) ; + G4double wZ = m_R[i] * cos(Theta / rad) ; - G4cout << "XXXXXXXXXXXX Trapezoid " << i << " XXXXXXXXXXXXX" << G4endl; - G4cout << "MMu: " << MMu << G4endl; - G4cout << "MMv: " << MMv << G4endl; - G4cout << "MMw: " << MMw << G4endl; + MMw = G4ThreeVector(wX, wY, wZ) ; +// G4ThreeVector CT = MMw ; + CT = MMw ; + MMw = MMw.unit() ; + + G4ThreeVector Y = G4ThreeVector(0 , 1 , 0) ; + + MMu = MMw.cross(Y) ; + MMv = MMw.cross(MMu) ; + MMv = MMv.unit(); + MMu = MMu.unit(); // Passage Matrix from Lab Referential to Telescope Referential + // MUST2 MMrot = new G4RotationMatrix(MMu, MMv, MMw); // Telescope is rotate of Beta angle around MMv axis. MMrot->rotate(m_beta_u[i], MMu); MMrot->rotate(m_beta_v[i], MMv); MMrot->rotate(m_beta_w[i], MMw); // translation to place Telescope - MMpos = MMw * Length * 0.5 + MMCenter; + MMpos = MMw * Length * 0.5 + CT ; } - FirstStage = m_wFirstStage[i]; - SecondStage = m_wSecondStage[i]; - ThirdStage = m_wThirdStage[i]; + FirstStage = m_wFirstStage[i] ; + SecondStage = m_wSecondStage[i] ; + ThirdStage = m_wThirdStage[i] ; - VolumeMaker(i + 1, MMpos, MMrot, FirstStage, SecondStage, ThirdStage, world); + VolumeMaker(i + 1, MMpos, MMrot, FirstStage, SecondStage, ThirdStage , world); } - delete MMrot; + delete MMrot ; } @@ -748,13 +752,12 @@ void GaspardTrackerTrapezoid::ReadSensitive(const G4Event* event) G4THitsMap<G4double>* AngPhiHitMap; // NULL pointer are given to avoid warning at compilation - // Second Stage - std::map<G4int, G4double*>::iterator SecondStageEnergy_itr; - G4THitsMap<G4double>* SecondStageEnergyHitMap = NULL; + // Third Stage std::map<G4int, G4double*>::iterator ThirdStageEnergy_itr; G4THitsMap<G4double>* ThirdStageEnergyHitMap = NULL; + // Read the Scorer associated to the first Stage //Detector Number G4int StripDetCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("FirstStageScorerGPDTrapezoid/DetectorNumber") ; @@ -806,12 +809,8 @@ void GaspardTrackerTrapezoid::ReadSensitive(const G4Event* event) AngPhiHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(InterCoordAngPhiCollectionID)) ; Ang_Phi_itr = AngPhiHitMap->GetMap()->begin() ; - // Read the Scorer associated to the Second and Third Stage - // Energy second stage - G4int SecondStageEnergyCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("SecondStageScorerGPDTrapezoid/SecondStageEnergy") ; - SecondStageEnergyHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(SecondStageEnergyCollectionID)) ; - SecondStageEnergy_itr = SecondStageEnergyHitMap->GetMap()->begin() ; - // Energy third stage + // Read the Scorer associated to the Third Stage + //Energy G4int ThirdStageEnergyCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThirdStageScorerGPDTrapezoid/ThirdStageEnergy") ; ThirdStageEnergyHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(ThirdStageEnergyCollectionID)) ; ThirdStageEnergy_itr = ThirdStageEnergyHitMap->GetMap()->begin() ; @@ -948,41 +947,25 @@ void GaspardTrackerTrapezoid::ReadSensitive(const G4Event* event) } // Second Stage - SecondStageEnergy_itr = SecondStageEnergyHitMap->GetMap()->begin() ; - for (G4int h = 0 ; h < SecondStageEnergyHitMap->entries() ; h++) { - G4int SecondStageEnergyTrackID = SecondStageEnergy_itr->first - N; - G4double SecondStageEnergy = *(SecondStageEnergy_itr->second); - - if (SecondStageEnergyTrackID == NTrackID) { - ms_Event->SetGPDTrkSecondStageEEnergy(RandGauss::shoot(SecondStageEnergy, ResoSecondStage)); - ms_Event->SetGPDTrkSecondStageEPadNbr(1); - ms_Event->SetGPDTrkSecondStageTPadNbr(1); - ms_Event->SetGPDTrkSecondStageTTime(1); - ms_Event->SetGPDTrkSecondStageTDetectorNbr(m_index["Trapezoid"] + N); - ms_Event->SetGPDTrkSecondStageEDetectorNbr(m_index["Trapezoid"] + N); - } - - SecondStageEnergy_itr++; - } // Third Stage - ThirdStageEnergy_itr = ThirdStageEnergyHitMap->GetMap()->begin() ; - for (G4int h = 0 ; h < ThirdStageEnergyHitMap->entries() ; h++) { - G4int ThirdStageEnergyTrackID = ThirdStageEnergy_itr->first - N; - G4double ThirdStageEnergy = *(ThirdStageEnergy_itr->second); - - if (ThirdStageEnergyTrackID == NTrackID) { - ms_Event->SetGPDTrkThirdStageEEnergy(RandGauss::shoot(ThirdStageEnergy, ResoThirdStage)); - ms_Event->SetGPDTrkThirdStageEPadNbr(1); - ms_Event->SetGPDTrkThirdStageTPadNbr(1); - ms_Event->SetGPDTrkThirdStageTTime(1); - ms_Event->SetGPDTrkThirdStageTDetectorNbr(m_index["Trapezoid"] + N); - ms_Event->SetGPDTrkThirdStageEDetectorNbr(m_index["Trapezoid"] + N); + ThirdStageEnergy_itr = ThirdStageEnergyHitMap->GetMap()->begin() ; + for (G4int h = 0 ; h < ThirdStageEnergyHitMap->entries() ; h++) { + G4int ThirdStageEnergyTrackID = ThirdStageEnergy_itr->first - N; + G4double ThirdStageEnergy = *(ThirdStageEnergy_itr->second); + + if (ThirdStageEnergyTrackID == NTrackID) { + ms_Event->SetGPDTrkThirdStageEEnergy(RandGauss::shoot(ThirdStageEnergy, ResoThirdStage)); + ms_Event->SetGPDTrkThirdStageEPadNbr(1); + ms_Event->SetGPDTrkThirdStageTPadNbr(1); + ms_Event->SetGPDTrkThirdStageTTime(1); + ms_Event->SetGPDTrkThirdStageTDetectorNbr(m_index["Trapezoid"] + N); + ms_Event->SetGPDTrkThirdStageEDetectorNbr(m_index["Trapezoid"] + N); + } + + ThirdStageEnergy_itr++; } - ThirdStageEnergy_itr++; - } - DetectorNumber_itr++; } @@ -997,7 +980,6 @@ void GaspardTrackerTrapezoid::ReadSensitive(const G4Event* event) PosZHitMap ->clear(); AngThetaHitMap ->clear(); AngPhiHitMap ->clear(); - SecondStageEnergyHitMap ->clear(); ThirdStageEnergyHitMap ->clear(); } } diff --git a/NPSimulation/src/ParisCluster.cc b/NPSimulation/src/ParisCluster.cc index 3ff1f1c10c6c8c3816e8b96f0033b6d584c9d305..464b611fe9489e03c7f70e49c1b6ed76d9a8da9c 100644 --- a/NPSimulation/src/ParisCluster.cc +++ b/NPSimulation/src/ParisCluster.cc @@ -564,6 +564,8 @@ void ParisCluster::ConstructDetector(G4LogicalVolume* world) MMv = m_X1_Y128[i] - m_X1_Y1[i]; MMv = MMv.unit(); + G4ThreeVector MMscal = MMu.dot(MMv); + MMw = MMu.cross(MMv); // if (MMw.z() > 0) MMw = MMv.cross(MMu) ; MMw = MMw.unit(); @@ -772,8 +774,6 @@ void ParisCluster::ReadSensitive(const G4Event* event) G4int sizeNCsI= CsIDetectorNumberHitMap->entries(); G4int sizeECsI= CsIStageEnergyHitMap->entries(); - sizeC *= 1; // remove warning at compilation - sizeECsI *= 1; // remove warning at compilation //G4cout <<"SizeN=" << sizeN << endl; //G4cout <<"SizeC=" << sizeC << endl; //G4cout <<"SizeN CsI =" << sizeNCsI << endl; @@ -806,7 +806,7 @@ void ParisCluster::ReadSensitive(const G4Event* event) G4double T = *(Time_itr->second); G4int NCryst= *(CrystalNumber_itr->second); - NCryst *= 1; //remove warning at compilation + //G4cout <<"NTrackID=" << NTrackID << G4endl; //G4cout <<"N_first=" << N_first << G4endl; //G4cout <<"CrystalNumber_first=" << NCryst << G4endl; diff --git a/NPSimulation/src/ParisPhoswich.cc b/NPSimulation/src/ParisPhoswich.cc index a0a956a8114f1989945a5223c3f295ad594f4c87..cd7a6a4f4342f650dab49819f8af4c21ee3ff4ab 100644 --- a/NPSimulation/src/ParisPhoswich.cc +++ b/NPSimulation/src/ParisPhoswich.cc @@ -475,6 +475,8 @@ void ParisPhoswich::ConstructDetector(G4LogicalVolume* world) MMv = m_X1_Y128[i] - m_X1_Y1[i]; MMv = MMv.unit(); + G4ThreeVector MMscal = MMu.dot(MMv); + MMw = MMu.cross(MMv); // if (MMw.z() > 0) MMw = MMv.cross(MMu) ; MMw = MMw.unit(); @@ -684,8 +686,6 @@ void ParisPhoswich::ReadSensitive(const G4Event* event) G4int sizeNCsI= CsIDetectorNumberHitMap->entries(); G4int sizeECsI= CsIStageEnergyHitMap->entries(); - sizeC *= 1; // remove warning at compilation - sizeECsI *= 1; // remove warning at compilation //G4cout <<"SizeN=" << sizeN << endl; //G4cout <<"SizeC=" << sizeC << endl; //G4cout <<"SizeN CsI =" << sizeNCsI << endl; @@ -717,7 +717,6 @@ void ParisPhoswich::ReadSensitive(const G4Event* event) G4double E = *(Energy_itr->second); G4double T = *(Time_itr->second); G4int NCryst= *(CrystalNumber_itr->second); - NCryst *= 1; // remove warning at compilation //G4cout <<"NTrackID=" << NTrackID << G4endl; diff --git a/NPSimulation/src/ShieldClParis.cc b/NPSimulation/src/ShieldClParis.cc index f3fe6a464dddf0450840d05d4fa8e398faea9632..2b77daa1872f600c4d441b9d0551885fcc843bc7 100644 --- a/NPSimulation/src/ShieldClParis.cc +++ b/NPSimulation/src/ShieldClParis.cc @@ -463,6 +463,8 @@ void ShieldClParis::ConstructDetector(G4LogicalVolume* world) //MMv = -0.5 * (m_X1_Y1[i] + m_X128_Y128[i] - m_X1_Y128[i] - m_X128_Y1[i]); MMv = MMv.unit(); + G4ThreeVector MMscal = MMu.dot(MMv); + MMw = MMu.cross(MMv); // if (MMw.z() > 0) MMw = MMv.cross(MMu) ; MMw = MMw.unit(); diff --git a/NPSimulation/src/ShieldPhParis.cc b/NPSimulation/src/ShieldPhParis.cc index 8d8edc55db700584e415c66ac2189649a034cdba..0f6b140dcf54dbd26f18d4ed95396aea0ea61cc1 100644 --- a/NPSimulation/src/ShieldPhParis.cc +++ b/NPSimulation/src/ShieldPhParis.cc @@ -473,6 +473,8 @@ void ShieldPhParis::ConstructDetector(G4LogicalVolume* world) //MMv = -0.5 * (m_X1_Y1[i] + m_X128_Y128[i] - m_X1_Y128[i] - m_X128_Y1[i]); MMv = MMv.unit(); + G4ThreeVector MMscal = MMu.dot(MMv); + MMw = MMu.cross(MMv); // if (MMw.z() > 0) MMw = MMv.cross(MMu) ; MMw = MMw.unit(); diff --git a/NPSimulation/vis.mac b/NPSimulation/vis.mac index 486ab51efdeaa869f5b0e15e31204710f8e4baa0..a5d707fb4bdc1aa0e81ebe7c24fcf0af0c750fda 100644 --- a/NPSimulation/vis.mac +++ b/NPSimulation/vis.mac @@ -9,19 +9,19 @@ # choose a graphic system ##/vis/open OGLIX ##/vis/open OGLSX -/vis/open VRML2FILE -/vis/scene/create -/vis/drawVolume -/vis/viewer/set/viewpointThetaPhi 0 0 deg -/vis/viewer/zoom 7 +#/vis/open VRML2FILE +#/vis/scene/create +#/vis/drawVolume +#/vis/viewer/set/viewpointThetaPhi 0 0 deg +#/vis/viewer/zoom 7 ## options to draw trajectories -/vis/scene/endOfEventAction accumulate -/vis/scene/add/trajectories 1 -/tracking/storeTrajectory 1 -/vis/scene/add/axes 0 0 0 20 cm -/vis/viewer/refresh +#/vis/scene/endOfEventAction accumulate +#/vis/scene/add/trajectories 1 +#/tracking/storeTrajectory 1 +#/vis/scene/add/axes 0 0 0 20 cm +#/vis/viewer/refresh # run event #/run/beamOn 0 -#/run/beamOn 1 +#/run/beamOn 10000