Diapositiva 1 - GEANT4 at LNS

Download Report

Transcript Diapositiva 1 - GEANT4 at LNS

Utilizzo del codice GEANT4 in campo medico LNS
(Mon 12-Oct-2009 ---- Wed 14-Oct-2009)
Gestione ed analisi dei dati prodotti
Francesco Di Rosa
14 - 10 - 2009
[email protected]
1
Outline
I.
Overview of basic Analysis
a)
b)
II.
Overview of advanced Analysis
a)
b)
c)
d)
e)
III.
Downloading the AIDA utility package
Turning on the analysis package
Loading the package
Producing/Analyzing first histograms
The packages: JAS, ROOT, PAW, OpenScientist
Few examples of Analysis
a)
b)
c)
IV.
Generic and simple Output in Geant4
Analysis of ASCII files
JAS
ROOT
PAW
Conclusion
14 - 10 - 2009
[email protected]
2
I. a) Generic and simple Output in Geant4
Geant4 does not attempt to provide its own data analysis tools, focusing instead on its
central mission as a simulation toolkit
The Geant4 user is expected to provide his own code to output results to an
appropriate analysis format
The Geant4 simple output: the simplest possible output format here
G4Cout
ASCII
*. TXT (Row, Matrices, Line)
*. CSV (Comma-Separated Values),
*. PIPPO (What you want!!!!!!!!! Exe: mixed formats)
The resulting files can be analyzed by tools such as:
Excel, OpenOffice, Gnuplot, Matlab, OriginLab
14 - 10 - 2009
[email protected]
3
I. a) Generic and simple Output in Geant4
An Example of G4COUT (good for test and verify)
G4cout
<< std::setw( 6) << "Energia depositata--->" << " " << std::setw( 6) << edep << " "
<< std::setw( 6) << "Carica--->" << " " << std::setw( 2) << particleCharge << " "
<< std::setw( 6) << "Energia Cinetica--->" << " " << std::setw( 2) << kineticEnergy << " "
<< G4endl;
14 - 10 - 2009
[email protected]
4
I. a) Generic and simple Output in Geant4
An Example of ASCII File (Save As……)
std::ofstream myfile(“File_di_OUTPUT”, std::ios::app);
myfile << edep << '\t' <<particleCharge<<'\t' << kineticEnergy << '\n' << G4endl;
controlloTally=evtNb;
14 - 10 - 2009
[email protected]
5
I. a) Generic and simple Output in Geant4
14 - 10 - 2009
[email protected]
6
I. b) Analysis of ASCII files
Excel
•
Comes from a small company called MicroSoft
•
Part of Microsoft Office (very common software)
From the “Data” menu, select “Get External Data”…“Import Data…”
 Select the file using the extension (*.*)
 The Text Import Wizard will come up
 Check the “Delimited” box
 Hit “Next”
 Check the “Comma” box
 The bottom part of the Wizard should now show your data nicely separated into
appropriate columns.
 Hit “Finish”
 You will see a dialog titled “Import Data” which will ask where to put the data.
 Accept the defaults:

Existing worksheet:

=$A$1
 Hit “OK”
The data should now show up in the spreadsheet (see the next slide)
14 - 10 - 2009
[email protected]
7
I. b) Analysis of ASCII files
Importing
data in
Excel
14 - 10 - 2009
[email protected]
8
I. b) Analysis of ASCII files
Plotting with Excel
Now that the data is in the
spreadsheet, you can plot it as
follows:
Select at least one item in the filled
data area.
From the “Insert” menu, select
“Chart…”
Select Chart type: “Histogram”
Hit “Finish”
You’re done 
14 - 10 - 2009
[email protected]
9
I. b) Analysis of ASCII files
OpenOffice
OpenOffice.org is a multiplatform and multilingual office suite
and an open-source project
Compatible with all other major office suites, the product is
free to download, use, and distribute
http://www.openoffice.org/
From the “File” menu,
Select “Open…”
In the “File type” box, select “Text CSV”
Then select your file Typology
……………………………………….
14 - 10 - 2009
[email protected]
10
I. b) Analysis of ASCII files
Importing Data from a .csv File Into OpenOffice
Once you select your file:
• The Text Import dialog will come up
• Check the “Comma” box
• The bottom part of the dialog should now show your data nicely
separated into appropriate columns, then Hit “OK”
The data should
now show up in the
spreadsheet 
(maybe)
14 - 10 - 2009
[email protected]
11
I. b) Analysis of ASCII files
Plotting with OpenOffice
Now that the data is in the spreadsheet, you can plot it as follows:
Select at least one item in the filled data area.
From the “Insert” menu, select “Chart…”
The defaults should be fine
Hit “Next”
• Select the Icon
for “Line” chart
• Hit “Create”
• You’re done! 
14 - 10 - 2009
[email protected]
12
I. b) Analysis of ASCII files
Gnuplot
•
“Gnuplot is a portable command-line driven interactive data and function plotting
utility for UNIX, IBM OS/2, MS Windows, DOS, Macintosh, VMS, Atari and many
other platforms.”
•
“The software is copyrighted but freely distributed (i.e., you don't have to pay for it)”
•
http://www.gnuplot.info/
•
See the demo page (http://gnuplot.sourceforge.net/demo_4.2/)
14 - 10 - 2009
[email protected]
13
I. b) Analysis of ASCII files
Plotting with MatLab
1) Choose the data Folder
4. One-click with right mouse bottom on
2) LOAD the
output file
Workspace
3) In the
appearand
your
datathe
5. workspace
Select Plot Catalog
prefer
Chart type
14 - 10 - 2009
[email protected]
6. You’re
done…
:-/
14
Outline
I.
Overview of basic Analysis
a)
b)
II.
Overview of advanced Analysis
a)
b)
c)
d)
e)
III.
Downloading the AIDA utility package
Turning on the analysis package
Loading the package
Producing/Analyzing first histograms
The packages: JAS, ROOT, PAW, OpenScientist
Few examples of Analysis
a)
b)
c)
IV.
Generic and simple Output in Geant4
Analysis of ASCII files
JAS
ROOT
PAW
Conclusion
14 - 10 - 2009
[email protected]
15
II. Overview of advanced Analysis (AIDA)
While the Geant4 collaboration does not maintain any data analysis system, many
data analysis tool developers have interfaced their systems to Geant4 and each of
the following has an active Geant4 user community
AIDA (Analysis Interface for Data Analysis) is an interface standard used by
several advanced data analysis tools. It provides a single standard way to
create your data
Commands to: instantiate/fill an nTuple
instantiate/fill a Histogram
You can switch from one AIDA-compliant tool to another without having to change
how you book and fill your histograms (http://aida.freehep.org/)
ROOT is an Object Oriented Data Analysis Framework. It is heavily
used in High Energy Physics. Freely available at http://root.cern.ch/
14 - 10 - 2009
[email protected]
16
II. Overview of advanced Analysis (AIDA)
AIDA Example
Use same code with any AIDA-compliant analysis tool
A
User code
(e.g.
GEANT4)
14 - 10 - 2009
I
Analysis tool 1
Analysis tool 2
D
A
[email protected]
Analysis tool 3
17
II. a) Downloading & Setting the AIDA utility package
Now, you have just AIDA installed……
But…..
If you DO NOT have AIDA installed, then (will take ~5-10 min):
1.
Create AIDA directory under "geant4/work" directory
mkdir geant4/work/AIDA
2.
Download "osc_batch-v16r0-XX-i386-gcc_401.zip" into "AIDA" folder from
http://aida.freehep.org
– Windows
XX = CYGWIN
– Linux
XX = Linux
– MACOS
XX = Darwin
3.
Decompress the file
4.
Go into “osc_batch/v16r2” directory
cd osc_batch/v16r0
5.
Setup the AIDA environment
source aida-setup.csh or source aida-setup.sh
6.
Done!
14 - 10 - 2009
[email protected]
18
II. b) Turning on the Analysis Package
• G4ANALYSIS_USE
• Set variable
– setenv G4ANALYSIS_USE 1 (csh)
– export G4ANALYSIS_USE=1 (bash o setup)
• Is it on? (linux commands)
– 1 (on) or 0 (off)
– “echo”
– “env”
14 - 10 - 2009
echo $G4ANALYSIS_USE
env then look for G4ANALYSIS_USE
[email protected]
19
II. b) Loading the analysis package I
Required files
#include <AIDA/AIDA.h>

Define
AnalysisManager.cc
AnalysisManager.hh

Create & Fill
14 - 10 - 2009
[email protected]
20
II. b) Loading the analysis package II: AnalysisManager.cc (create)
From Hadrontherapy Example
14 - 10 - 2009
[email protected]
21
II. b) Loading the analysis package II: AnalysisManager.cc (create)
From Extended Example A01
………….
A01AnalysisManager::A01AnalysisManager()
:analysisFactory(0), hFactory(0), tFactory(0)
{
// Hooking an AIDA compliant analysis system.
analysisFactory = AIDA_createAnalysisFactory();
if(analysisFactory)
{
ITreeFactory* treeFactory = analysisFactory->createTreeFactory();
tree = treeFactory->create("A01.aida","xml",false,true,"compress=yes");
hFactory = analysisFactory->createHistogramFactory(*tree);
tFactory = analysisFactory->createTupleFactory(*tree);
delete treeFactory; // Will not delete the ITree.
}
}
……..
14 - 10 - 2009
[email protected]
22
II. b) Loading the analysis package II: AnalysisManager.cc (create)
From Extended Example AnaEx01: create histograms
…
std::string opts = "compress=no";
fTree = treeFactory->create("AnaEx01.aida","xml",false,true,opts);
…..
// Create an histo factory that will create histo in the tree :
AIDA::IHistogramFactory* histoFactory = fAIDA>createHistogramFactory(*fTree);
if (histoFactory)
{
fEAbs = histoFactory->createHistogram1D("EAbs",100,0,100);
fLAbs = histoFactory->createHistogram1D("LAbs",100,0,100);
…
delete histoFactory;
}
…
14 - 10 - 2009
[email protected]
23
II. b) Loading the analysis package III: AnalysisManager.cc (create)
From Extended Example AnaEx01: create ntuples
…
// Get a tuple factory :
AIDA::ITupleFactory* tupleFactory =
fAIDA->createTupleFactory(*fTree);
if (tupleFactory)
{
// Create a tuple :
fTuple = tupleFactory->create("AnaEx01","AnaEx01",
"double EAbs,double LAbs,double EGap,double LGap");
delete tupleFactory;
}
…
14 - 10 - 2009
[email protected]
24
II. b) Loading the analysis package IV: AnalysisManager.cc (fill)
From Hadrontherapy Example
14 - 10 - 2009
[email protected]
25
II. b) Loading the analysis package IV: AnalysisManager.cc (fill)
From Extened Example AnaEx01: fill histograms and ntuples
…
if (CHC)
{
G4int n_hit = CHC->entries();
for (G4int i=0;i<n_hit;i++) {
G4double EAbs = (*CHC)[i]->GetEdepAbs();
G4double LAbs = (*CHC)[i]->GetTrakAbs();
…
fEAbs->fill(EAbs);
fLAbs->fill(LAbs);
…
fTuple->fill(0,EAbs);
fTuple->fill(1,LAbs);
…
fTuple->addRow();
}
}…
14 - 10 - 2009
[email protected]
26
II. c) Producing/Analyzing first histograms
Homework assignment
 Copy example into your work area (hadrontherapy)
 Go into example directory: cd hadrontherapy
 Compile: gmake or make
 Run: ~/geant4/work/bin/Linux-g++/……
 A new file should appear: hadronterapy.root (.hbk, .aida)
ANALYSIS

Start JAS (…ROOT,AIDA,HBOOK)
 Open the file: File  “Open File” then select XX.root
 Double click on “Blue Folder”
 Double click on histogram of interest
14 - 10 - 2009
[email protected]
27
Outline
I.
Overview of basic Analysis
a)
b)
II.
Overview of advanced Analysis
a)
b)
c)
d)
e)
III.
Downloading the AIDA utility package
Turning on the analysis package
Loading the package
Producing/Analyzing first histograms
The packages: JAS, ROOT, PAW, OpenScientist
Few examples of Analysis
a)
b)
c)
IV.
Generic and simple Output in Geant4
Analysis of ASCII files
JAS
ROOT
PAW
Conclusion
14 - 10 - 2009
[email protected]
28
III. Few examples of advanced analysis
• Java Analysis Studio:
JAS (all files)
• The Root of Everything:
ROOT (root files)
• Physics Analysis Workstation:
PAW (hbook files)
• OpenScientist
(all files)
14 - 10 - 2009
[email protected]
29
III a) Java Analysis Studio (JAS, http://jas.freehep.org/jas3)
Java Analysis Studio (JAS)
– Written entirely in Java (really portable)
– Originally targeted at offline analysis (but also used
extensively for online monitoring)
– Rich, easy to use GUI
• Can view trees, plots, HTML pages with embedded active objects
– Built-in editor/compiler for writing analysis code
• Based on open source JEdit editor (pure Java editor, syntax
highlighting, auto indent, parenthesis matching, etc.)
– Local and Client-Server Operation
• Access local or remote analysis objects
• Can send analysis code to the data and ship back only the results
14 - 10 - 2009
[email protected]
30
III a) Java Analysis Studio (Example)
• Start jas
• Open the file: File 
“Open File” then select
OUTPUT.aida
• Double click on “Blue
Folder” (OUTPUT.aida)
• Double click on histogram
and then the ones of
interest
• Double click on tuples
then “Yellow Folder”
14 - 10 - 2009
[email protected]
31
III a) Java Analysis Studio (GUI)
Pages can display histograms. User can
control layout, add, remove plots.
Built-in HTML viewer
with embeddable
“objects” (plots,
buttons, etc.).
Tree area
provides
access to
analysis
objects. Each
item has
popup menu
and doubleclick action.
Plot Widget,
shows data in
real time,
optimized for
fast refresh
performance.
14 - 10 - 2009
[email protected]
32
III b) ROOT (http://root.cern.ch/ - http://www.mi.infn.it/~palombo/lazzaro/ )
•
In AnaEx01AnalysisManager.cc
 Verify or type the following lines
std::string opts = “export=root";
fTree = treeFactory->create("AnaEx01.root",“ROOT",false,true,opts);
•
•
•
•
•
•
•
Start root
Launch browser: new TBrowser()
Open the file: File  “Open” then select nameFile.root
Double click on “Yellow Folder” (ROOT Files)
Double click on “Yellow Folder” (file.root)
Double click on “Yellow Folder” (histograms) and then the ones of interest
Double click on “Yellow Folder” (tuples) then right-mouse on “nameFile” and select Start
Viewer
Drag variables of interest into (X,Y,Z) (Ex: X=EAbs and Y=Labs)
Select Draw Icon (bottom left corner)
Create a cut: double click on E() then type: EGapcut (alias) and EGap<20 (expression)
Drag cut into Scissor
Select Draw Icon (bottom left corner)
Can enable or disable the cut by double clicking on it
•
•
•
•
•
•
14 - 10 - 2009
[email protected]
33
III b) ROOT (Hist Browser)
14 - 10 - 2009
[email protected]
34
III b) ROOT (TGhtml web browser plug-in)
URL
You can
browse a root
file
You can execute a
script
14 - 10 - 2009
[email protected]
35
III b) ROOT (Macro Manager/Editor plug-in)
Click on button to
execute script
14 - 10 - 2009
[email protected]
36
III b) ROOT (Math libraries)
14 - 10 - 2009
[email protected]
37
III c) Physics Analysis Workstation (PAW, cernlib)
PAW is no longer supported: you are on your own!!!
•
In AnaEx01AnalysisManager.cc
 Verify or type the following lines
std::string opts = “export=hbook";
fTree = treeFactory->create("AnaEx01.hbook",“HBOOK",false,true,opts);
•
•
•
•
•
•
•
•
•
•
Should get: AnaEx01.hbook
Start paw++
Open file: h/file 0 [filepath]/AnaEx01.hbook
Double click on logical unit 1 (lun1)
Double click on histograms or ntuple
Select a variable and double lick in (X,Y,Z) (Ex: X=EAbs and Y=Labs)
Click on Plot
Click on Cut Editor …
Left box = variable, middle box = expression, right box = value (Ex: EGap<20)
Click on Plot
14 - 10 - 2009
[email protected]
38
Conclusion
Geant4 does not attempt to provide its own data
analysis tools, focusing instead on its central mission as
a simulation toolkit
Several analysis tools are available
The Geant4 user is expected to provide his own code
to output results to an appropriate analysis format
Use as it suits your needs
Real
Conclusion
14 - 10 - 2009
Use Whatever Analysis Tool You Like
[email protected]
39
MORE
INFORMATIONS
ABOUT AIDA
14 - 10 - 2009
[email protected]
40
n -Tuple (The term originated as an abstraction of the sequence: single, double,
triple, quadruple, quintuple, n-tuple)
Let us use Example
Let us consider the following example A Data Summary Tape contains 10000
events.
Each event consists of many variables (200) for which we would like to make some
distributions according to various selection criteria
One possibility is to create and fill 200 histograms on an event-by-event basis
while reading the DST
An alternative solution is to create one Ntuple
Instead of histograming the 200 variables directly, and therefore losing the exact
values of the variables for each event and their correlations, the variables are first
stored in an Ntuple.
An Ntuple is like a table where the 200 variables mentioned above are the columns
and each event is a row
14 - 10 - 2009
[email protected]
41
AIDA Factories – IAnalysisFactory
• Factory – “a building or set of buildings with
facilities for manufacturing” (Webster Dictionary).
• Objects in AIDA are created by Factories or internally by
other objects. User should newer utilize “new”.
 static IAnalysisFactory create ()







IHistogramFactory createHistogramFactory (ITree tree)
IFunctionFactory createFunctionFactory (ITree tree)
ITupleFactory createTupleFactory (ITree tree)
ITreeFactory createTreeFactory ()
IPlotterFactory createPlotterFactory ()
IDataPointSetFactory createDataPointSetFactory (ITree tree)
IFitFactory createFitFactory ()
14 - 10 - 2009
[email protected]
42
AIDA Factories – IHistogramFactory
• Creates 1D, 2D, 3D IHistograms, IClouds, IProfiles:
– IHistogram1D createHistogram1D (String name, String title,
int nBins, double lowerEdge, double upperEdge, String options )
– IHistogram1D createHistogram1D (String name, String title,
int nBins, double[] binEdges, String options )
– Options: “type=efficiency” – create efficiency histogram.
– Options: “autoconvert=true” – to enable auto conversion for Clouds.
• Does operations on those objects:
– Copy and Arithmetic operations add/subtract/divide/multiply
IHistogram1D add (String name, IHistogram1D h1, IHistogram1D h2)
– Projections and Slices of 2D and 3D Objects: 2D -> 1D, 3D -> 2D
IHistogram1D projectionX (String name, IHistogram2D h2D)
IHistogram1D sliceX (String name, IHistogram2D h2D, int index)
14 - 10 - 2009
[email protected]
43
IHistogram (1D-3D)
• Binned histogram: IHistogram1D, 2D,
3D
– “fill” methods (with/without weight)
– Histogram info: entries, mean, rms, axis
– Bin info: centre, entries, height, error
– Histogram arithmetic: add, multiply, divide
– Convenience methods, like coordinate-toindex conversion
14 - 10 - 2009
[email protected]
44
ITuple
• ITuple - interface to the Data
– “get/set” methods for double, float, int, …
– Information about columns: min, max, mean, rms
– Navigating: start(), next(), skip(int nRows)
– Project ITuple into 1D, 2D, 3D histogram
– New features for AIDA 3.0:
• Support for complex internal structures (subfolders)
• Chaining of ITuples
14 - 10 - 2009
[email protected]
45
ICloud
• Unbinned collection of points: ICloud1D, 2D,
3D
– Can represent scatter plot, dynamically
rebinnable histogram
– Can be converted to a binned histogram
– Standard “get/set” methods for entries
– Collection info: lower, upper, mean, rms
14 - 10 - 2009
[email protected]
46
IFunction and Fitting
• Fitting: IFunction, IFitFunction
– IFunction – simple interface, allows to set
parameters and get function value
– IFitFunction – fit function to a histogram
•
•
•
•
Extends IFunction
Various fit control methods: step size, bounds, etc.
Allows to perform fit and get results
AIDA 2.2 fitting functionality fairly limited
– AIDA 3.0 (Under discussion) extended functionality
14 - 10 - 2009
[email protected]
47
ITree
• ITree
– directory-like structure (Unix directory convention)
• Methods like: cd, ls, mkdir, etc.
– AIDA analysis objects (tuples, histograms, clouds,
ets.) exist within ITree directories
– “save/restore” functionality, hides storage details
from the user
• Compatible with database or file storage
• Can support multiple file formats
• Mount/Unmount functionality (like unix) allows multiple
stores to be seamlessly merged
• AIDA XML format is defined for data interchange
14 - 10 - 2009
[email protected]
48