Monte carlo and analysis tool study

1. Monte carlo and analysis tool study Shi-hong yao

2. motivation Geant4 -- interest, because it’s seem to very useful for many science domain. -- to understandhow to design a big program? -- work need…if I need to add/change something. ROOT, Mysql -- for analysis, histogram, data store.

3. In your Main Required G4VUserDetectorConstruction User must be implement his own classes drivered from the three abstract classes respective and register them to the RunManager. G4VUserPhysicsList G4VUserPrimaryGeneratorAction Optional G4RunManager G4UserRunAction You can collect data from these classes or even change/interfere the process in order to reach some purpose. G4UserEventAction G4UserStackingAction G4RunManager is the only manager class in Geant4 kernel. SetUserInitialization(…) SetUserAction(…) Initialize(); beamOn(numberOfEvent); G4UserTrackingAction G4UserSteppingAction

4. In your Main G4UIManager G4UIterminal In interactive mode, G4UIsession* session = new G4UIterminal; session->SessionStart(); // Work… work… work…your work… delete session; It’s hard-coded batch mode, G4UImanager* UI =G4UImanager::GetUIpointer(); UI->ApplyCommand(“/run/verbose 1"); UI->ApplyCommand(“/control/execute macrofile"); #macro file #set verbose level /run/verbose 2 #start run /run/beamOn 100 Idle> your command

5. Geant4 as a state machine PreInit Initialized Idle beamOn exit GeomClose Quit EventProc Abort

6. Basic - Units • Length(meter) km, m, cm, mm, um, nm, fm … • Surface(meter2) km2, m2, cm2, mm2, barn, mbarn, … • Angle(radian) rad, mrad, sr, deg … …and a lot of extend unit… • Energy(electron volt) • eV, keV, MeV, …, PeV J(joule) • Mass(gram) • mg, g, kg • Time(second) • s, ms, mus, ns, ps The units in blue are defined by one. global

7. When you read-in a value, it is recommend to set the units. • If the units are not specified, it implicitly use the internal units. • You can output your data with the units you wish. • By divide the unit. • You can get the list of units by static function G4UnitDefinition::PrintUnitsTable() or UI command /units/list . density = 12.0*g/cm3; length = 1*cm; G4cout << density/(g/cm3); G4cout << length/cm;

8. Call from G4RunManager Construct Detector G4VUserDetectorConstruction myDetectorConstruction Public your own concrete class Construct(); //pure virtual method which is invoked by G4RunManager when it's Initialize() method is invoked. This method must return the G4VPhysicalVolume pointer which represents the world volume.

9. Material – molecule or mixture It include: Name, atomic number Z, number of nucleons N, atomic mass A, Or, Isotopes. It include: Name, Density, State, Temperature, Pressure, Or Components G4Element H x2 G4Material H2O G4Element O x1 G4Material Other … G4Element* elH = new G4Element(name="Hydrogen", symbol="H" , z= 1., a = 1.01*g/mole); G4Element* elO = new G4Element(name="Oxygen" ,symbol="O" , z= 8., a 16.00*g/mole); G4Material* H2O = new G4Material(name="Water", density = 1.000*g/cm3, ncomponents=2); H2O->AddElement (elH, natoms=2); H2O->AddElement (elO, natoms=1); Or get material from Geant4 database G4NistManager* man = G4NistManager::Instance(); G4Material* H2O = man->FindOrBuildMaterial("G4_WATER"); See the list: Geant4 Material Database

10. Detector construct Return to G4RunManager SD World volume G4PVPlacement PhysicsWorld G4LogicalVolume LogiVolDaughter G4Tube G4PVPlacement PhysicsVolume G4LogicalVolume LogiVolMother Air G4PVPlacement other G4BOX G4LogicalVolume LogiVolOther … SolidBox = new G4Box(name, xLangth, yLangth, zLangth); LogiVolMother = new G4LogicalVolume( SolidBox, Air, name2); LogiVolDaughter = new G4LogicalVolume( SolidTube, H2O, name3); PhysicsVolume = new G4PVPlacement( Rot, trans, LogiVolDaughter , name?, LogiVolMother, 0, copyNo); PhysicsWorld = new G4PVPlacement(0, G4ThreeVector(), LogicalVolume, “WORLD”, 0, false, 0);

11. Call from G4RunManager Physics list G4VUserPhysicsList myPhysicsList Public your own concrete class Three pure virtual method must be implemented by the user. ConstructParticle(); // construction of particles ConstructProcess(); // construct processes and register them to particles SetCuts(); // setting a range cut value for all particles Other method for using: AddTransportation(); //register the G4Transportation class which describes the particle motion in space and time with all particles. SetCutsWithDefault();//

12. G4BosonConstructor pure virtual method. User must define "All PARTICLE TYPES" include not only primary particles, but also all secondaries. G4Proton G4LeptonConstructor G4Gamma ConstructParticle() G4MesonConstructor G4Electron G4BarionConstructor G4MuonPlus G4IonConstructor … G4ShortlivedConstructor. More than 100 types of particles are provided by default. Define particles one by one or use this six utility classes, corresponding to each of the particle categories. In constructParticle(), explicitly invoke static methods of all particle classes G4Proton::ProtonDefinition(); G4Gamma::GammaDefinition(); G4MuonMinus::MuonMinusDefinition(); … #include “G4LeptonConstructor.hh” G4LeptonConstructor pConstructor; pConstructor.ConstructParticle();

13. What is process Physics processes describe how particles interact with materials. G4VProcess G4Transportation The GPIL method gives the step length from the current space-time point to the next space-time point. It does this by calculating the probability of interaction based on the process's cross section information. Abstract base class for all physics processes. GPIL method: PostStepGetPhysicalInteractionLength(…); AlongStepGetPhysicalInteractionLength(…); AtRestGetPhysicalInteractionLength(…); DoIt method: AlongStepDoIt(..); PostStepDoIt(…); AtRestDoIt(…); G4Cerenkov G4ComptonScattering The DoIt method implements the details of the interaction, changing the particle's energy, momentum, direction and position, and producing secondary tracks if required. G4eplusAnnihilation … seven major categories of processes, electromagnetic, hadronic, transportation, decay, optical, photolepton_hadron, and parameterisation.

14. ConstructProcess G4Gamma G4PhotoElectricEffect G4ProcessManager G4ComptonScattering Register related physics process to each particle type. This maybe a big project. … AddTransportation(); particle = G4Gamma::GammaDefinition(); G4ProcessManager* pmanager = particle->GetProcessManager(); G4PhotoElectricEffect * thePhotoElectricEffect = new G4PhotoElectricEffect(); pmanager->AddDiscreteProcess(thePhotoElectricEffect); pmanager->AddDiscreteProcess(theComptonEffect); pmanager->AddDiscreteProcess(thePairProduction); … a lot ...

15. Range Cut SetCuts pure virtual method. A "unique cut value in range“ should be defined as a distance which is internally converted to an energy for individual materials. Geant4 recommend cut value is 1.0 mm. SetCutsWithDefault(); or defaultCutValue = 1.0*mm; //set cut value for different particle types. SetCutValue(cutForGamma, "gamma"); SetCutValue(cutForElectron, "e-"); …

16. Call from G4RunManager Primary Generation G4VUserPrimaryGeneratorAction myPrimaryGeneratorAction Public your own concrete class GeneratePrimaries(G4Event* anEvent ); Pure virtual method which is invoked at the beginning of each event. This class generates primary particle(s) with a given momentum and position. G4VPrimaryGenerator G4ParticleGun Public G4GeneralParticleSource For spatial distributions http://reat.space.qinetiq.com/gps/ This is an abstract base class of all of primary generators. GeneratePrimaryVertex() = 0; G4HEPEvtInterface G4int n_particle = 1; particleGun = new G4ParticleGun(n_particle); particleGun->SetParticleDefinition(G4Proton::ProtonDefinition()); particleGun->SetParticleEnergy(120.0*GeV); particleGun->SetParticlePosition(G4ThreeVector(0.0*m, 0.0*m, -5.0*m)); particleGun->GeneratePrimaryVertex(anEvent);

17. PYTHIA8 PYTHIA is a program for the generation of high-energy physics events, i.e. for the description of collisions at high energies between elementary particles such as e+, e-, p and pbar in various combinations. ……Excerpt from PYTHIA’s site

18. Three step of PYTHIA Pythia pythia; pythia.readString(“ProcessGroup:ProcessName = ON|OFF”); //List see appendix 1 pythia.readFile(FileName); … pythia.init(idA, idB, eA, eB); //pythia.init(idA, idB, eCM); //Currently the program only works with pp, pbar p, e+e-and μ+μ-incoming beams. • Only PYTHIA itself: • Initialization • Event loop • Finish, statistics while(!pythia.next()); //event record found in pythia.event pythia.event.list(); //do something…. pythia.statistics();

20. PYTHIA play a generator role myPrimaryGeneratorAction PYTHIA G4ParticleGun Particle remnants sent to the particle gun. while(!pythia.next()); for (int i = 1; i < pythia.event.size(); i++) { Pythia8::Particle par = pythia.event[i]; if (par.status() > 0) { ParticleGun -> SetParticlePosition(G4ThreeVector(par.xProd()*mm,par.yProd()*mm,par.zProd()*mm)); ParticleGun -> SetParticleDefinition(ParticleTable->FindParticle(par.id())); ParticleGun -> SetParticleMomentumDirection(G4ThreeVector(par.px(),par.py(),par.pz())*GeV); ParticleGun -> SetParticleEnergy(par.pAbs()*GeV); ParticleGun -> GeneratePrimaryVertex(anEvent); } }

21. visualization HepRep/WIRED What you can visualize -Detector components -A hierarchical structure of physical volumes -Particle trajectories and tracking steps -Hits of particles in detector components. OpenGL DAWN ???

22. OpenGL #Create an empty scene ("world" is default): /vis/scene/create #Create a scene handler for a specific graphics system /vis/open OGLIX #/vis/open DAWNFILE #/vis/open HepRepFile # for drawing the tracks,hits /vis/scene/add/trajectories /vis/scene/add/hits /tracking/storeTrajectory 1 #too many tracks may cause core dump /vis/scene/endOfEventAction accumulate | refresh /run/beamOn number View directly from Geant4 when you are running, it can zoom, rotate, translate, but poor for graphics. Addition on follow two line before session start. G4VisManager* visManager = new G4VisExecutive; visManager->Initialize(); Then execute your program and stop at idle state.

23. WIRED Replace OGLIX by HepRepFile beamOn and produce G4DataX.heprep by Geant4, X is the event number. View in the WIRED Event Display, you need to install java runtime environment. And download WIRED from here. >java -jar HepRApp.jar -file G4DataX.heprep It can zoom, rotate, translate,… click to show attributes, control visibility from hierarchical (tree) view of data Export to many vector graphic formats (PostScript, PDF, etc.)

24. Data store -- Mysql • Benefits It is a well-designed database manager system. It is convenient for query specify data with some condition. Some calculate capability, average, maximum, minimum, and so on. Run two or more program parallelize without worry about the file conflict. • Short You need one computer to do the server. You can’t run simulation without the server. Efficiency depend on the network. You still need other analysis tool for complex calculate.

25. Using Mysql API Get the data you want to store from Geant4 or PYTHIA. According to the Mysql syntax, combine the data into string. Through function mysql_query() or mysql_ #include “mysql.hh” if you can’t find this file, please contact your system manager…. MYSQL *conn = mysql_init(NULL); mysql_real_connect(conn, host,user,pass,database, 0, NULL, 0) mysql_query(conn, querystring); //receive returned data… result = mysql_store_result(conn); for(num_rows=0;row = mysql_fetch_row(result);num_rows++){ eventlist[num_rows]=atoi(row[0]); } mysql_free_result(result); mysql_close(conn);

26. Data store– ROOT • Benefits ROOT is a very efficient and complete frameworks for analysing large amounts of data. Root file can store many type of object include TTree. Structure of TTree is similar to mysql but root provide many fancy function for drawing histogram. The data directly store at local disk – fast. • Short If data has some problem or error, it’s inconvenient to find out and correct.

27. User Action //at the beginning of the BeamOn() method virtual void BeginOfRunAction( const G4Run* aRun); //at the end of the BeamOn() method virtual void EndOfRunAction( const G4Run* aRun); G4UserRunAction virtual void BeginOfEventAction( const G4Event* anEvent); virtual void EndOfEventAction( const G4Event* anEvent); G4UserEventAction G4UserStackingAction void PrepareNewEvent (); G4ClassificationOfNewTrack ClassifyNewTrack(const G4Track*) void NewStage () NewStage G4UserTrackingAction void PreUserTrackingAction( const G4Track*); void PostUserTrackingAction( const G4Track*) G4UserSteppingAction void UserSteppingAction(const G4Step*)

28. Feature work • I have not yet understand all the detail/concept of Geant4. • Try to understand the physics process. Thanks for Listening.

29. Appendix 1 Currently implemented processes, complete with respect to groups, but with some individual processes missing for lack of space (represented by “...”). In the names, a “2” separates initial and final state, an “(s:X)”, “(t:X)” or “(l:X)” occasionally appends info on an s- or t-channel- or loop-exchanged particle X.

30. back