Formation GEANT4 - LPC

1. Formation GEANT4 - LPC Emmanuel Delage, Loïc Lestand, Yann Perrot, Bogdan Vulpescu

2. Programme de la formation • Jeudi 15 novembre : introduction générale à Geant4 • Vendredi 23 novembre : exploitation + visualisation • Jeudi 29 novembre : matériaux + géométrie • Jeudi 6 décembre : physique et particules • Jeudi 13 décembre : récupération des données

3. EXTRAIRE DES DONNéES(tiré en partie des slides de M. Asai, Tutorial G4, SLAC 2012)

4. Extraire des données • Une expérience GEANT4 incluant la géométrie, la physique, la génération des particules primaires ne produit AUCUNE sortie • L’utilisateur DOIT donner du code pour enregistrer les données qu’il désire! • La récupération de données est le scoring

5. Niveaux de simulation GEANT4 Expérience Nombre total particules primaires Évènement Émission 1 particule primaire Interactions physiques Transport des particules Étape Énergie cinétique nulle primaire + secondaires Fin expérience

6. Méthodes pour extraire des données • Utilisation de scorers dans les volumes (primitive scorers) • Créer un score à chaqueévènement • Créersapropreclasse Run/Event pour accumuler les scores • Assigner un G4VSensitiveDetector à un volume pour générer des “hits”. • Implémenter des classes utilisateurs (G4UserEventAction, G4UserRunAction) pour récupérer au niveau de l’évènement, de la simulation • Utilisation des classes utilisateurs (G4UserTrackingAction, G4UserSteppingAction,…) • Accès à toutes les informations • Il faut tout écriresoimême

7. Volume sensible • Un volume sensible, représenté par un objet de la classeG4VSensitiveDetector (détecteur) peutêtreassigné à un volume logique G4LogicalVolume • Le rôle principal d’un volume sensible est de construire des objets de type G4Hit à partir des informationsobtenus au cours des étapes (steps) lors du suivi des particules. • Dans le casoù un step a lieu dans un volume logiqueayantétédéclaré volume sensible, l’objetG4VSensitiveDetectorestinvoqué avec l’objet courant G4Step • Vouspouvezimplémentervospropres classes définissantvosdétecteurs (sensitive detector) oubienutiliser les outilsdédiésofferts par GEANT4 (primitive scorers)

8. Collection de Hits • G4VHitsCollection est la classeabstraite des deux classes G4THitsCollection et G4THitsMap. • G4THitsCollectionest un vecteur de pointeursvers des objets de type G4Hit • G4THitsMapest un objet de type “map” qui enregistre des clés (numéro de copie du volume) associées à des pointeursd’objets • Cesobjets ne sont pas nécessairement des hits (G4double pour les primitive scorers de Geant4) • G4THitsMap peutnéanmoinsêtreutilisé par le volume sensible pour enregistrer des objets de type hit.

9. G4HCofThisEvent • A la fin d’un évènementréussi, un objet G4Event contient un objet du type G4HCofThisEvent. • Un objet G4HCofThisEventcontienttoutes les collections de hits généréesdurantl’évènement: • Pointeur(s) vers les collections NULL siaucune collection crééedurantl’évènement. • Chaque collection estréférencée par un index unique et inchangédurantune simulation. L’indexpeutêtreobtenu de la façonsuivante: G4SDManager::GetCollectionID(“detName/colName”); NB: la table des index estaussienregistréedans G4Run

10. Sensitive detector Vousdevezimplémentervospropres classes gérant le détecteur et les hits. Uneclassedéfinissant un hit peutcontenirplusieursquantités. Un hit peutcontenirchaque step individuellementouaccumuler des grandeurs. Une collection de hits estcréée pour chaquedétecteur. Primitive scorer Nombreux scorers fournis . Possibilitéd’ajouter les siens. Chaque scorer accumuleune grandeur pour un évènement. G4MultiFunctionalDetector crééune collection par scorer. Sensitive detector ou primitive scorer? En pratique: • Primitive scorers • Si vousn’êtes pas intéressé par l’enregistrement de chaque step maisplutôt par l’accumulation de grandeurs pour un évènementouune simulation et • Si vousn’avez pas besoin de trop de scorers. • Sinon, penser à implémenter son propre sensitive detector.

11. Primitive scorers

12. Liste des Primitive Scorers • Longueur du track • G4PSTrackLength, G4PSPassageTrackLength • Énergie déposée • G4PSEnergyDeposit, G4PSDoseDeposit, G4PSChargeDeposit • Courant/flux • G4PSFlatSurfaceCurrent, G4PSSphereSurfaceCurrent, G4PSPassageCurrent • G4PSFlatSurfaceFlux, G4PSCellFlux, G4PSPassageCellFlux • Autres • G4PSMinKinEAtGeneration, G4PSNofSecondary, G4PSNofStep

13. Les filtres • On peut être sélectifs quant aux informations récupérées via l’utilisation de filtres • charged * Filtreparticulechargée • neutral * Filtreparticule de charge neutre • Particle * Filtresur le nom de la particule • kineticEnergy * Filtresurunegammed’énergiecinétique • particleWithKineticEnergy * Filtreparticule + énergiecinétique

14. Définir un Primitive ScorerEn pratique… • Dans DetectorConstruction::Construct(): • Créer un détecteur, il s’agit d’un objet G4MultiFunctionalDetector • Attacher un volume logique à ce détecteur • Créer un scorer • Éventuellement, ajouter un filtre au scorer • Enregistrer le scorer au détecteur

15. Exemple B4d G4LogicalVolume* absorberLV = new G4LogicalVolume(…); // Scorers // declare Absorber as a MultiFunctionalDetector scorer G4MultiFunctionalDetector* absDetector = new G4MultiFunctionalDetector("Absorber"); G4SDManager::GetSDMpointer()->AddNewDetector(absDetector); absorberLV->SetSensitiveDetector(absDetector); G4VPrimitiveScorer* primitive; primitive = new G4PSEnergyDeposit("Edep"); absDetector->RegisterPrimitive(primitive); primitive = new G4PSTrackLength("TrackLength"); G4SDChargedFilter* charged = new G4SDChargedFilter("chargedFilter"); primitive ->SetFilter(charged); absDetector->RegisterPrimitive(primitive);

16. Accumulation des données • Principe: • Récupérer la collection de hits : G4HitsMap • Sommer la grandeur • Application: • Dans un code utilisateur, quel est le moment le plus adapté pour récupérer/accumuler les informations? • Copier dans votre répertoire de travail l’exemple /usr/local/geant4.9.5.p01/examples/basic/B4/B4d • Dans B4d, éditer le code source de la classe qui gèrera la récupération des hits et l’accumulation des données

17. Exemple B4dRécupération des hits voidB4dEventAction::EndOfEventAction(const G4Event* event) { // Get sum value from hits collections // G4double absoEdep = GetSum(GetHitsCollection("Absorber/Edep", event)); } Pour rappel dans B4dDetectorConstruction: G4MultiFunctionalDetector* absDetector = new G4MultiFunctionalDetector("Absorber"); G4VPrimitiveScorer* primitive; primitive = new G4PSEnergyDeposit("Edep"); absDetector->RegisterPrimitive(primitive);

18. Exemple B4dRécupération des hits GetHitsCollection("Absorber/Edep", event) G4THitsMap<G4double>* B4dEventAction::GetHitsCollection(const G4String& hcName, const G4Event* event) const { G4int hcID = G4SDManager::GetSDMpointer()->GetCollectionID(hcName); G4THitsMap<G4double>* hitsCollection = static_cast<G4THitsMap<G4double>*>( event->GetHCofThisEvent()->GetHC(hcID)); // Gestion des erreurs if ( ! hitsCollection ) { G4cerr << "Cannot access hitsCollection " << hcName << G4endl; exit(1); } return hitsCollection; }

19. Exemple B4dSomme de l’énergiedéposée GetSum(GetHitsCollection("Absorber/Edep", event)) G4double B4dEventAction::GetSum(G4THitsMap<G4double>* hitsMap) const { G4double sumValue= 0; std::map<G4int, G4double*>::iterator it; for ( it = hitsMap->GetMap()->begin(); it != hitsMap->GetMap()->end(); it++) { sumValue += *(it->second); } return sumValue; }

20. Exemple B4dExemple de fichier de sortie simple G4double absoEdep= GetSum(GetHitsCollection("Absorber/Edep", event)) FILE* myOutputFile; myOutputFile=fopen("output.txt", "a"); fprintf(myOutputFile, "#Event %i -> %d \n", event->GetEventID(), absoEdep); fclose(myOutputFile); Tester avec la macro run2.mac où vous aurez modifié le nombre de particules à simuler Comment avoir des sorties plus élaborées ??

21. Root

22. ROOT • Outil d’analyse de données: • Lecture/écriture de fichiers; • Manipulation d’histogrammes, de structures de données; • Graphiques, fits; • Calculs analytiques ou numériques; • Et tant d’autres choses! • http://root.cern.ch • Informations générales; • Téléchargement/installation; • Tutoriaux, guides utilisateur; • Forums.

23. ROOT et C++ • ROOT consiste en une collection de classes C++ • Code ROOT = code C++ interprété par CINT • Exemples d’objets ROOT: • Fichier ROOT: TFile • Histogramme 1D, 2D: TH1F, TH2F • Ntuple: TNtuple • Fonction 1D: TF1 • Graphique: TGraph • Pour utiliser ROOT: • De manière interactive • Par le biais de scripts (macro) root[0] .X script.C Hello World ! 4.79583 root[1]

24. Exemple B4dLecture d’un fichier ROOT 1/4 • Dans un terminal, lancer root: root • Lancer un explorateur: new TBrowser()

25. Exemple B4dLecture d’un fichier ROOT 2/4 • Ouvrir le fichier root File->Open • Le fichier sélectionné est ajouté • Pour visualiser le contenu: double-click, ici histogrammes et n-tuples

26. Exemple B4dLecture d’un fichier ROOT 3/4

27. Exemple B4dLecture d’un fichier ROOT 4/4

28. ROOT et GEANT4 • Geant4 dispose d’un manager d’analyse: interface uniforme pour manipuler les données dans différents formats (ROOT, AIDA, XML, CSV et HBOOK) • Idée : utiliser le manager pour créer, remplir et sauvegarder des histogrammes, des tables, … • Exemple B4d : déclaration de l’utilisation de ROOT dans le fichier d’entête B4Analysis.hh #include "g4analysis_defs.hh " using namespace G4Root;

29. 1- Créer le fichier ROOT En début de simulation: B4RunAction::BeginOfRunAction // Createanalysis manager // The choice of analysis technology is done via selectin of a namespace in B4Analysis.hh G4AnalysisManager* analysisManager = G4AnalysisManager::Instance(); // Open an output file G4String fileName = "B4"; analysisManager->OpenFile(fileName); analysisManager->SetFirstHistoId(1);

30. 2- Créer la structure du fichier ROOT En début de simulation: B4RunAction::BeginOfRunAction // Creatinghistograms analysisManager->CreateH1("1","Edep in absorber", 100, 0., 800*MeV); analysisManager->CreateH1("2","Edep in gap", 100, 0., 100*MeV); analysisManager->CreateH1("3","trackL in absorber", 100, 0., 1*m); analysisManager->CreateH1("4","trackL in gap", 100, 0., 50*cm); // Creatingntuple // analysisManager->CreateNtuple("B4", "Edep and TrackL"); analysisManager->CreateNtupleDColumn("Eabs"); analysisManager->CreateNtupleDColumn("Egap"); analysisManager->CreateNtupleDColumn("Labs"); analysisManager->CreateNtupleDColumn("Lgap"); analysisManager->FinishNtuple();

31. 3- Remplir le fichier ROOT A la fin de chaque évènement: B4dEventAction::EndOfEventAction // fillhistograms analysisManager->FillH1(1, absoEdep); analysisManager->FillH1(2, gapEdep); analysisManager->FillH1(3, absoTrackLength); analysisManager->FillH1(4, gapTrackLength); // fillntuple analysisManager->FillNtupleDColumn(0, absoEdep); analysisManager->FillNtupleDColumn(1, gapEdep); analysisManager->FillNtupleDColumn(2, absoTrackLength); analysisManager->FillNtupleDColumn(3, gapTrackLength); analysisManager->AddNtupleRow();

32. 4- Sauvegarder le fichier ROOT A la fin de la simulation: B4RunAction::EndOfRunAction G4AnalysisManager* analysisManager = G4AnalysisManager::Instance(); // savehistograms analysisManager->Write(); analysisManager->CloseFile();

33. Exercice d’application

34. Dénombrer le nombre de particules secondaires A partir de l’exampleB4d, dénombrer le nombre de photons et d’électrons secondaires créés dans le volume Absorber en utilisant des primitive scorers et des filtres. Stocker les données dans deux Ntuples ROOT. 1- Modifier la classe B4dDetectorConstruction pour y inclure les scorers. G4PSNofSecondary::G4PSNofSecondary(G4String name) G4SDParticleFilter::G4SDParticleFilter(G4String name_filter, const G4String& particleName) 2- Modifier la classe B4dEventAction pour dénombrer le nombre d’électrons et de photons secondaires à chaque évènement. 3- Modifier la classe B4RunAction pour créer les deux Ntuples 4- Modifier la classe B4dEventAction pour enregistrer à chaque évènement le nombre d’électrons et de photons secondaires dans les Ntuples ROOT.

35. Dénombrer le nombre de particules secondairesUne solution 1- Modifier la classe B4dDetectorConstruction pour y inclure les scorers. Déclarations #include "PSNofSecondary.hh" #include "G4SDParticleFilter.hh" Dans DefineVolumes() primitive= new G4PSNofSecondary("SecElec"); G4SDParticleFilter* elecFilter = new G4SDParticleFilter("elecFilter","e-"); primitive->SetFilter(elecFilter); absDetector->RegisterPrimitive(primitive); Idem avec G4PSNofSecondary("SecGamma"), et G4SDParticleFilter("elecGamma",”gamma");

36. Dénombrer le nombre de particules secondairesUne solution 2- Modifier la classe B4dEventAction pour dénombrer le nombre d’électrons et de photons secondaires à chaque évènement. Dans EndOfEventAction: G4int nSecElec = G4int( GetSum( GetHitsCollection("Absorber/SecElec",event) ) ); G4int nSecElec = G4int( GetSum( GetHitsCollection("Absorber/SecElec",event) ) );

37. Dénombrer le nombre de particules secondairesUne solution 3- Modifier la classe B4RunAction pour créer les deux Ntuples Dans BeginOfRunAction: //Creatingntuple … … analysisManager->CreateNtupleDColumn("SecGamma"); analysisManager->CreateNtupleDColumn("SecElec"); analysisManager->FinishNtuple();

38. Dénombrer le nombre de particules secondairesUne solution 4- Modifier la classe B4dEventAction pour enregistrer à chaque évènement le nombre d’électrons et de photons secondaires dans les Ntuples ROOT. Dans EndOfEventAction: //fillntuple … … analysisManager->FillNtupleDColumn(4,nSecGamma); analysisManager->FillNtupleDColumn(5,nSecElec); analysisManager->AddNtupleRow();

39. Utilisation avancéeSenSitive detector & hits

40. Volume sensible et hits • Chaque volume logiquepeutavoir un pointeurvers un détecteur (sensitive detector) • Ce volume devientalorssensible • Hit est un étatd’une trace dans un volume sensible • Un volume sensible créé des hits à partir des informationscontenusdansuneétape (objet G4Step). L’utilisateurdoitdonnerl’implémentation de la réponse du détecteur. • Les hits, qui sont des objetsdéfinis par l’utilisateur, sontcollectés à la fin de chaqueévènementdans un objet G4Event.

41. Les informations à partir de G4Step • Un step ->2 points • Un point = 1 position • 1 point = 1 volume physique, 1 volume logique, matériau, numéro de copie, etc… • 1 point = 1 processus qui a limité le step • Un step -> • Énergiedéposée • Longueurd’étape • … • Un step -> Une trace • Moment • Énergiecinétique.. • … • Et biend’autresinformations, voir les méthodesd’accès: http://geant4.web.cern.ch/geant4/UserDocumentation/UsersGuides/FAQ/html/ch01s04.html

42. La classe Hit • Un Hit estuneclasseutilisateurdérivée de G4VHit. • L’utilisateurstocke les donnéesqu’ildésire en implémentantsapropreclasse Hit. Par exemple: • Position et temps du step • Moment, énergie de la trace • Dépôtd’énergie du step • Information géométrique • Les objets Hit doiventêtrestockésdansune collection dédiéedérivée de la laclasse G4THitsCollection. • La collection estassociée à un évènement via G4HCofThisEvent. • La collection est accessible: • Par G4Event à la fin de l’évènement • Par G4SDManager lors de l’analyse d’un évènement (filtre par exemple)

43. Implémentation Hit #include "G4VHit.hh" class MyHit : public G4VHit { public: MyHit(some_arguments); virtual ~MyHit(); virtual void Draw(); virtual void Print(); private: // some data members public: // some set/get methods }; #include “G4THitsCollection.hh” typedef G4THitsCollection<MyHit> MyHitsCollection;

44. Implémentation Sensitive Detector #include "G4VSensitiveDetector.hh" #include "MyHit.hh" class G4Step; class G4HCofThisEvent; class MyDetector : public G4VSensitiveDetector { public: MyDetector(G4String name); virtual ~MyDetector(); virtual void Initialize(G4HCofThisEvent*HCE); virtual G4bool ProcessHits(G4Step*aStep, G4TouchableHistory*ROhist); virtual void EndOfEvent(G4HCofThisEvent*HCE); private: MyHitsCollection * hitsCollection; G4int collectionID; };

45. CLASSE UTILISATEURS

46. Principe Exécution 1- L’utilisateur doit créer si besoin ses propres classes pour décrire les différents niveaux. 2- Les informations sont recueillies au niveau de l’étape à partir d’un objet G4Step 3- Les niveaux de simulation ne communiquent pas les uns avec les autres, l’utilisateur doit inclure dans ses classes des moyens pour communiquer entre les différentes classes (pointeurs, méthodes lecture/écriture de données) Début expérience BeginOfRunAction RunAction:: G4UserRunAction EndOfRunAction Que faire au début de l’expérience? Début évènement BeginOfEventAction EventAction:: G4UserEventAction EndOfEventAction Étape Que faire au début d’un évènement? SteppingAction:: G4UserSteppingAction UserSteppingAction Fin évènement Que faire lors d’une interaction? Fin expérience Que faire à la fin d’un évènement? Que faire à la fin de l’expérience? Fermeture

47. PrincipeExemple simple EventAction.hh public: EventAction(); BeginOfEventAction(); EndOfEventAction(); CollectInfo(); private: G4double Info; SteppingAction.hh public: SteppingAction(EventAction *); UserSteppingAction(); private: EventAction* fEventAction; Au début de l’évènement BeginOfEventAction: Remettre à zéro la quantité A chaque StepUserSteppingAction: Collecter l’information stockée grâce à la méthode d’accès implémentée dans EventAction fEventAction->CollectInfo();

48. L’exemple B4aEventAction.hh class B4aEventAction : public G4UserEventAction { public: … voidAddAbs(G4double de, G4double dl); voidAddGap(G4double de, G4double dl); private: … G4double fEnergyAbs; G4double fEnergyGap; G4double fTrackLAbs; G4double fTrackLGap; };

49. L’exemple B4aEventAction.cc void B4aEventAction::BeginOfEventAction(const G4Event* evt) { … // initialisation per event fEnergyAbs = 0.; fEnergyGap = 0.; fTrackLAbs = 0.; fTrackLGap = 0.; }

50. L’exemple B4aSteppingAction.hh class B4aEventAction; class B4aSteppingAction : public G4UserSteppingAction { public: B4aSteppingAction(const B4DetectorConstruction* detectorConstruction, B4aEventAction*eventAction); virtual ~B4aSteppingAction(); virtual void UserSteppingAction(const G4Step* step); private: const B4DetectorConstruction* fDetConstruction; B4aEventAction* fEventAction; };