Title: GALPROP: under the hood
1GALPROPunder the hood
Igor V. Moskalenko (stanford/kipac)
IAS, Princeton/Apr.23, 2009
2GALPROP on the Web
- A de-facto standard in astrophysics of CRs,
diffuse gamma rays, and indirect DM searches
(Blasi 2007 ICRC rapporteur talk) - The diffuse emission model has become a baseline
for the Fermi/LAT - Freely distributed easy access for the forefront
science - GALPROP cosmic gives gt5000 hits on GOOGLE
- Forum topic are visited gt60000 times since it was
opened in 2006 - The principle authors have written 50 reviews
and refereed papers based on its results (highly
cited) - Consistent funding from NASA since 2001 through
2012 _at_ 160k/year (ATP and APRA programs, PI
IM) funding of the proposals which use GALPROP
3All Particle CR Spectrum
- This is an astonishing observation!
- All particle CR spectrum is almost featureless
- the knee
- the ankle
- GZK cutoff
- These are the only features in gt12 decades in
energy and gt32 decades in intensity! - However, there is a lot of information hidden in
the spectra and abundances of individual CR
species nuclear isotopes, antiprotons,
electrons, positrons (diffuse gamma rays) - CRs are the only probes of the interstellar
material available to us. - The whole physics is involved various branches
of Astrophysics, MHD, shock waves, plasma
physics, atomic, nuclear, particle physics,
exotic physics SUSY
Galactic
Galacticextragalactic
GZK cutoff?
extragalactic
4CR Propagation Milky Way Galaxy
1 pc 3x1018 cm
Optical image Cheng et al. 1992, Brinkman et al.
1993 Radio contours Condon et al. 1998 AJ 115,
1693
100 pc
NGC891
40 kpc
0.1-0.01/ccm
1-100/ccm
Sun
4-12 kpc
Intergalactic space
R Band image of NGC891 1.4 GHz continuum (NVSS),
1,2,64 mJy/ beam
Flat halo model (Ginzburg Ptuskin 1976)
5CR Interactions in the Interstellar Medium
ISM
X,?
synchrotron
IC
ISRF
P He CNO
diffusion energy losses
reacceleration convection
production of secondaries
bremss
p
Flux
LiBeB
He CNO
20 GeV/n
BESS
- CR species
- Only 1 location
- modulation
PAMELA
ACE
helio-modulation
6Elemental Abundances CR vs. Solar System
CR abundances ACE
output
O
Si
C
Na
Fe
S
CNO
Al
CrMn
Cl
LiBeB
F
ScTiV
Solar system abundances
input
7Cosmic Rays vs Diffuse Gamma Rays
Even an unrealistic model (e.g. Leaky-Box) can be
fitted to the CR data, but diffuse emission
requires the CR spectra in the whole Galaxy
8A Model of CR Propagation in the Galaxy
- Gas distribution (energy losses, p0, brems)
- Interstellar radiation field (IC, e energy
losses) - Energy losses ionization, Coulomb, brems, IC,
synch - Nuclear particle production cross sections
- Solve transport equations for all CR species
- Fix propagation parameters
- Gamma-ray production brems, IC, p0
- Precise Astrophysics
9Components of the ISM Views from the Inside
synchrotron 21 cm H I 2.6 mm CO (H2) dust stars
star forming optical n-stars, BHs CRs x gas
10CO maps
- Extend CO surveys to high latitudes
- newly-found small molecular clouds will otherwise
be interpreted as unidentified sources, and
clearly limit dark matter studies - C18O observations (optically thin tracer) of
special directions (e.g. Galactic center, arm
tangents) - assess whether velocity crowding is affecting
calculations of molecular column density, and for
carefully pinning down the diffuse emission
Dame, Hartmann, Thaddeus (2001) Dame Thaddeus
(2004)
11Gas distribution in the Milky Way
Molecular hydrogen H2 is traced using J1-0
transition of 12CO, concentrated mostly in the
plane (z70 pc, Rlt10 kpc) Atomic
hydrogen H I has a wider distribution
(z1 kpc, R30 kpc) Ionized hydrogen H II
small proportion, but exists even in halo (z1
kpc)
Sun
12Distribution of interstellar gas
- Neutral interstellar medium most of the
interstellar gas mass - 21-cm H I 2.6-mm CO (surrogate for H2)
- Differential rotation of the Milky Way plus
random motions, streaming, and internal velocity
dispersions is largely responsible for the
spectrum - Rotation curveV(R) ? unique line-of-sight
velocity-Galactocentric distance relationship
CO
Rotation Curve
Dame et al. (2001)
H I
Clemens (1985)
Kalberla et al. (2005)
W. Keel
13Problems to deal with
- H2 gas
- The conversion factor XN(H2)/WCO is unknown,
most probably variable, and is determined from
the diffuse ?-ray emission itself - HI gas
- Spin temperature is unknown, it can significantly
vary usually used the same temperature 125K for
HI gas for the whole Galaxy - Self absorption (cold gas cloud in front of the
emitting cloud) the optical depth is very large
14More on gas in the Milky Way
Surface mass density of the H2 in M sun pc-2
25
G.C.
- Problems
- Near-far ambiguity
- No velocity information in the Center-Anticenter
direction
Pohl08
15Dealing with neutral gas galactocentric rings
25
G.C.
sum
- Examples of the Galactocentric rings
sum
16Central Molecular Zone
Feriere07
17Inner 3 kpc
18GALPROPunder the hoodPart II
Igor V. Moskalenko (stanford/kipac)
IAS, Princeton/Apr.23, 2009
19ISRF Large Scale Distribution
- Requires extensive modeling
- Distribution of stars of different stellar
classes in the Galaxy - Dust emission
- Radiative transfer
- The z scale height is large, takes 10s of kpc at
R 0 kpc to get to level of CMB
Total
Optical
IR
CMB
R 0,4,8,12,16 kpc
Energy Density
Optical IR (no CMB)
20ISRF Angular distribution
- Calculation
- Calculate intensity maps throughout Galaxy
- stellar component
- dust component
- radiative transfer
- These vary with position AND wavelength
- Different from local intensity distribution
- Isotropic approximation for anything but CMB is
clearly incorrect
Porter Strong
21HESS Observations of Composite SNR G0.90.1
leptonic scenario
SNR at the GC age a few kyr
K-N cutoff
electron index d2.8
photon index 2.400.110.20
Porter,IVM,Strong06 IVM,Porter,Strong06
- Interstellar radiation field in the inner Galaxy
is dominated by the dust (IR) emission and
starlight
22Nuclear Reaction NetworkCross Sections
- Many different isotopes are produced via
spallations of CR nuclei A(p,He)?BX
stable isotopes
Secondary, radioactive 1 Myr K-capture isotopes
Co57
Fe55
Mn54
Cr51
V49
Ca41
Ar37
Cl36
ß-, n
Al26
p,EC,ß
Be7 Be10
p
Plus some dozens of more complicated
reactions But many cross sections are not well
known
n
23Effect of Cross Sections Radioactive Secondaries
Different size from different ratios
T1/2?
Zhalo,kpc
- Errors in CR measurements (HE LE)
- Errors in production cross sections
- Errors in the lifetime estimates
24Nuclear Reaction Network
I
II
III
IV
V
nuc_package.cc
25nuc_package.cc Long-Lived Isotopes
Half-life for bare nucleus Half-life for H-like
atoms
26Total Nuclear Cross Sections
Ekin, MeV/nucleon
Wellisch Axen 1996
27Fe p data vs phenomenology
Large deviations everywhere!
28Fe p data vs codes
- Very large deviations at low energies!
Villagrasa05
29LiBeB Major Production Channels
- Well defined (65)
- C12, O16 -gtLiBeB
- N14 -gt Be7
- (IVM Mashnik 2003)
- Few measurements
- C13,N -gt LiBeB
- B -gt BeB
- Unknown
- LiBeB,C13,N -gt LiBeB
- Tertiary reactions also important! -35
Propagated Abundance Cross-section
O
Li6
C
16
12
Li
B
N
Be
7
13
A
10
9
15
14
11
30Dealing with isotopic production Xsections
- Use data where possible
- Not much data use data in conjunction with the
results of a proper nuclear code - Not much data, no code use Webber or ST
semi-empirical systematics renormalized to the
data - Nothing available use Webber or ST
semi-empirical systematics
31Isotopic Production Cross Sections of LiBeB
- Semi-empirical systematics are not always
correct. - Results obtained by different groups are often
inconsistent and hard to test. - Very limited number of nuclear measurements
- Evaluating the cross section is very laborious
and cant be done without modern nuclear codes. - Use LANL nuclear database and modern computer
codes.
IVM Mashnik 2003
32isotope_cs_eval.dat
33Fitting the CR Abundances with GALPROP
Initial abundances, e.g., solar isotopic
abundances
Fine adjustment of the source abundances NSAOSA
d(ACE-propagated)
Propagation (GALPROP) 64Ni 1H
NSAnew source abundance OSAold source
abundance d0.01-0.001
Propagation parameters
Comparison with ACE data
Solar modulation (force-field)
Diffuse gammas
34CR source isotopic abundances
- Two K-capture isotopes are present in the
sources! -- 41Ca, 53Mn - Could tell us about the origin of CRs -- supports
volatility hypothesis, but needs more analysis
Solar system Reacceleration Plain diffusion
IVM,Strong,Porter07
- The first time that a realistic propagation model
has been used to derive isotopic source
abundances !
35SNR distribution
Milky Way
36Distribution of CR Sources Gradient in the CO/H2
Pulsar distribution Lorimer 2004
CR distribution from diffuse gammas (Strong
Mattox 1996) SNR distribution (Case
Bhattacharya 1998)
sun
XCON(H2)/WCO Histo This work, Strong et
al.04 ----- -Sodroski et al.95,97 1.9x1020
-Strong Mattox96 Z-1 Boselli et
al.02 Z-2.5 -Israel97,00, O/H0.04,0.07
dex/kpc
37Galactic magnetic field
Regular B-field large-scale structure
in the plane
in the halo
Han08
synchrotron emission
- Plane bisymmetrical field with reversals on
arm-interarm boundaries - Halo azimuth B-fields with reversed directions
below and above the plane - Random field Regular field
- Consistent with observations of the synchrotron
emission
38Transport Equation
39Transport Equations 90 (no. of CR species)
sources (SNR, nuclear reactions)
diffusion
convection (Galactic wind)
diffusive reacceleration (diffusion in the
momentum space)
E-loss
radioactive decay
fragmentation
- ?(r,p,t) density per total momentum
boundary conditions
40Simplified equation VHE electrons
Simplified eq.
Cylindrically symmetric solution
2p
Bessel fns
hypergeometric fn
dhalo size aradius
zeros of J0
41Electron propagation solutions
Galactic disk with sources
Galactic halo boundary
N
no energy losses
large energy losses
N0
The Galaxy
42Analytical vs. Numerical solution
- Analytical solutions for simple cases give
insight into the relations between the quantities
involved and are useful for rough estimates, but
in real cases the analytical formulae may become
too complicated so no insight is gained - Electrons and positrons are beyond the analytical
methods because their energy losses are spatially
dependent and different processes are important
in different energy ranges - Numerical solutions are fast and more realistic
- It is unclear whether one would wish to go much
beyond the generalizations discussed above for an
analytically soluble diffusion model. The added
insight from any analytic solution over a purely
numerical approach is quickly cancelled by the
growing complexity of the formulae. With rapidly
developing computational capabilities, one could
profitably employ numerical solutions
J.M.Wallace 1981
43GALPROPunder the hoodPart III
Igor V. Moskalenko (stanford/kipac)
IAS, Princeton/Apr.23, 2009
44Analytical vs. Numerical solution
- Analytical solutions for simple cases give
insight into the relations between the quantities
involved and are useful for rough estimates, but
in real cases the analytical formulae may become
too complicated so no insight is gained - Electrons and positrons are beyond the analytical
methods because their energy losses are spatially
dependent and different processes are important
in different energy ranges - Numerical solutions are fast and more realistic
- It is unclear whether one would wish to go much
beyond the generalizations discussed above for an
analytically soluble diffusion model. The added
insight from any analytic solution over a purely
numerical approach is quickly cancelled by the
growing complexity of the formulae. With rapidly
developing computational capabilities, one could
profitably employ numerical solutions
J.M.Wallace 1981
45Finite Differencing
46Finite Differencing Example
47Tri-Diagonal Matrix
48Coefficients for the Crank-Nicholson Method
49Spatial Grids
- Typical grid steps (can be arbitrary!)
- ?z 0.1 kpc, ??z 0.01 kpc (gas averaging)
- ?R 1 kpc
- ?E x1.2 (log-grid)
50GALPROP Output/FITS files
- Provides literally everything
- All nuclei and particle spectra in every grid
point (x,y,R,z,E) -FITS files - Separately for p0-decay, IC, bremsstrahlung
- Emissivities in every grid point
(x,y,R,z,E,process) - Skymaps with a given resolution (l,b,E,process)
- Output of maps separated into HI, H2, and rings
to allow fitting X-factors, CR source
distribution etc.
51How It Is Done Secondary Particles
- Positrons/electrons pp-gtp,K-gte (MS 1998)
- Dermer 1986 method LE Stecker ?-isobar model
(isotropic decay), HE scaling (inv. x-section
Stephens Badhwar 1981), plus interpolation in
between - Pion decay includes polarization of muons
- Kaon decay scaling (Stephens Badhwar 1981)
- Kamae et al. formalism
- Antiprotons (IM02)
- pp inclusive production x-section (Tan Ng 1983)
- pA, AA-gt pbar scaling using Gaisser Schaeffer
1992 or Simon et al. 1998 results similar (Mori
coeff. will be implemented) - Total inelastic x-section (p TN83, A Moiseev
Ormes 1997) - ppbar annihilation x-section (ppbar)tot
(pp)tot (LE TN83, HE PDG00 Regge
parameterization)
52How It Is Done Nucleons
- Calculated for pA reactions and scaled for aA
(Ferrando et al. 1988) - Calculation of total nuclear cross sections
- Letaw et al. 1983
- Wellisch Axen 1996 (corrected), Zgt5
- Barashenkov Polanski 2001
- Calculation of isotopic production cross sections
- Webber et al. 1993 (non-renormalized,
renormalized) Egt200 MeV/nucleon, essentially
flat - Silberberg Tsao 2000 (non-renormalized,
renormalized) claim that it works at all
energies, but is problematic sometimes - Fits to the available data (LANL, Webber et al.,
etc.) in the form of a function or a table (see
.dat files), but data may not be always available - Use the best of all three, but very time
consuming work - Nuclear reaction network
- Nuclear Data Tables (includes several decay
channels branching) - Standard ß -decay, emission of p, n
- K-capture isotopes are treated separately
53How It Is Done Gammas
- Bremsstrahlung (Koch Motz 1959, SMR2000) many
different regimes - LE 0.01 lt Ekin lt 0.07 MeV nonrelativistic
non-screened brems - Intermediate 0.07 lt Ekin lt2 MeV
- HE Ekin gt 2 MeV arbitrary screening unshielded
charge, 1-, 2-electron atoms (form factors,
Hylleraas, Hertree-Fock wave functions) - Fano-Sauter limit k-gtEkin
- Anisotropic IC (MS2000)
- Takes into account the anisotropic angular
distribution of background photons - Troy is working on new faster implementation
- Neutral pion decay (see secondary
positrons/electrons) - Synchrotron radiation (Ginzburg 1979, Ghisellini
et al. 1988) - Averaging over pitch angle
- Uses total magnetic field (regular random)
- Emissivities uses real H2, H I gas column
densities (rings) - Skymap calculations integration over the line of
sight
54- Propagation in the interstellar medium
55Secondary/Primary Ratio
- Leaky-box model
- fitting path-length distribution -gt free
function
- Diffusion models
- Diffusive reacceleration
- Convection
- Damping of interstellar turbulence
- Etc.
B/C
Ek, MeV/nucleon
Accurate measurements in a wide energy range may
help to distinguish between the models
56 if the mirror is moving towards the incident
particle, the particle gains energy upon
reflection, just as does a tennis ball pushed by
a racket
57Distributed Stochastic Reacceleration
Scattering on magnetic turbulences
Simon et al. 1986 Seo Ptuskin 1994
Dpp p2Va2/D D vR1/3 - Kolmogorov spectrum
Icr
1/3
?E
strong reacceleration
- Fermi 2-nd order mechanism
weak reacceleration
Dxx 5.2x1028 (R/3 GV)1/3cm-2 s-1 Va 36 km
s-1 ? R-d, d1.8/2.4 below/above 4 GV
E
58Diffusion in Galactic magnetic fields
Magnetic turbulences random field
inhomogeneities random walk
59Convection
Galactic wind
Jones 1979
DR0.6
Xe
v
R-0.6
wind or turbulent diffusion
resonant diffusion
E
problem too broad sec/prim peak
Dxx 2.5x1028 (R/4 GV)0.6cm-2 s-1 dV/dz 10 km
s-1 kpc-1 ? R-d, d2.46/2.16 below/above 20 GV
60Damping of Interstellar Turbulence
Kolmogorov cascade
Iroshnikov-Kraichnan cascade
nonlinear cascade
W(k)
- Simplified case
- 1-D diffusion
- No energy losses
dissipation
k
Mean free path
1/1020cm
1/1012cm
Ptuskin2003, 2005
61Dxx Diffusion Coefficient
R0.6
Reacceleration with damping
Plain diffusion
ß-3
Diffusive reacceleration
62Fixing Propagation Parameters
Radioactive isotopes Galactic halo size Zh
E2 Flux
B/C
Carbon
Be10/Be9
Interstellar
Ek, GeV/nucleon
Ek, MeV/nucleon
Zh increase
- Using secondary/primary nuclei ratio flux
- Diffusion coefficient and its index
- Propagation mode and its parameters (e.g.,
reacceleration VA, convection Vz) - Propagation params are model-dependent
- Make sure that the spectrum is fitted as well
Ek, MeV/nucleon
Parameters (model dependent) D 1028 (?/1 GV)a
cm2/s a 0.3-0.6 Zh 4-6 kpc VA 30 km/s
63Secondary/Primary Nuclei Ratio
B/C
B/C
CREAM Ahn08
Jones01
Sub-Fe/Fe
- Secondary/primary nuclei ratio in CR is declining
gt1 GeV/n, not rising!
64B/C ratio _at_ HE
- Different propagation models are tuned to fit the
low energy part of sec./prim. ratio where the
accurate data exist
ACE Ulysses Voyagers
Jones01
Reacceleration
CREAM Ahn08
Standard diffusion
- However, they differ at high energies which will
allow to discriminate between them when more
accurate data will be available
65Effect of Cross Sections Radioactive Secondaries
Different size from different ratios
T1/2?
Zhalo,kpc
- Errors in CR measurements (HE LE)
- Errors in production cross sections
- Errors in the lifetime estimates
66Light isotopes
- Light isotopes are produced mostly by the
CNO-nuclei (but tertiary component from LiBeB is
also important) relatively easy to calculate,
but there is no much data above 1 GeV/n - Also provide constraints on the propagation models
de Nolfo06
67Fitting the CR Abundances with GALPROP
Initial abundances, e.g., solar isotopic
abundances
Fine adjustment of the source abundances NSAOSA
d(ACE-propagated)
Propagation (GALPROP) 64Ni 1H
NSAnew source abundance OSAold source
abundance d0.01-0.001
Propagation parameters
Comparison with ACE data
Solar modulation (force-field)
Diffuse gammas
68CR source isotopic abundances
- Two K-capture isotopes are present in the
sources! -- 41Ca, 53Mn - Could tell us about the origin of CRs -- supports
volatility hypothesis, but needs more analysis
Solar system Reacceleration Plain diffusion
IVM,Strong,Porter07
- The first time that a realistic propagation model
has been used to derive isotopic source
abundances !
69Pamela positron fraction
- Excess in positron fraction is confirmed and
extended to higher energies
sec. production (GALPROP)
Solar modulation
Adriani08 (arXiv0810.4995)
70Fermi LAT ee- spectrum
- no prominent spectral features between 20 GeV
and 1 TeV - significantly harder spectrum than inferred from
previous measurements
- events for e e- analysis required to fail ACD
vetoes for selecting g events resulting g
contamination lt 1 - further cuts distinguish EM and hadron events
rejection 1103 up to 200 GeV 1104 at 1 TeV - energy reconstruction aided by shower imaging
capability of calorimeter - more than 4x106 e-e events in selected sample
- Simple normalized E-3.04 power law fits the Fermi
data reasonably well (c2 9.7, d.o.f. 24,
including estimated systematic uncertainties)
71HESS electron (photon) spectrum
HESS
- Sharp cutoff at 1 TeV
- No strong local sources?
- Features at higher energies?
72Pamela pbar/p ratio
- Pbar/p ratio is consistent with secondary origin
GALPROP
Adriani08 (arXiv0810.4994)
73Pamela pbars
- Absolute pbar flux is consistent with secondary
origin
From M.Boezio talk at SLAC
74Fermi/LAT First 3 months rate skymap
75Diffuse emission at mid-latitudes
sin 100.17 ZH2 50 pc sin 200.34 vs ZHI
400 pc
- Conventional GALPROP model is in agreement with
the LAT data at mid-latitudes (mostly local
emission)
76Morphology of the Diffuse Emission _at_ 150 GeV
IC
p0
RegisUllio09
IC
p0
77Milagro Skymap
- 34 Fermi BSL Galactic sources above declination
of -5o - 14 detected by Milagro above 3s
- FDR Miller 2001 estimates 1 false positive rate
- 5 new TeV sources
- Geminga 6.3s as extended source (2.6o fwhm)
Credit G.Sinnis
78Milagro TeV Observations of Fermi Sources
unID (new TeV source) Fermi Pulsar MGRO
190806 HESS 1908063
IC433 SNR MAGIC VERITAS
Geminga Pulsar Milagro C3
Radio pulsar J063110 (new TeV source)
unID (new TeV source)
Pulsar (AGILE/Fermi) MGRO 201937
Fermi Pulsar g Cygni SNR Fermi Pulsar HESS
203241 MGRO 203141 MAGIC 20324130
Fermi Pulsar Milagro (C4) 3EG 22276122 Boomerang
PWN
G65.10.6 (SNR) Fermi Pulsar (J1958) New TeV
sources
W51 HESS J1923141 SNR
Credit G.Sinnis
79Anisotropic Inverse Compton Scattering
- Electrons in the halo see anisotropic radiation
- Observer sees mostly head-on collisions
Energy density
e-
R4 kpc
small boost less collisions
e-
head-on large boost more collisions
Z, kpc
?
?
?
Important _at_ high latitudes !
sun
80Effect of anisotropic ICS
Ratio anisoIC/isoIC
pole
anti-GC
Intermediate latitudes
GC
- The anisotropic IC scattering plays important
role in modeling the Galactic diffuse emission - Affects estimates of isotropic extragalactic
background
Galactic latitude, degrees
81Extragalactic Gamma-Ray Background
EGRB in different directions
Predicted vs. observed
E2xF
Sreekumar et al. 1998
Elsaesser Mannheim05
Strong et al. 2004
E, MeV
- Blazars
- Cosmological neutralinos
82(No Transcript)
83Contributions to the extragalactic background
?gt100!
Dermer07
84First 3 months of LAT data
PSR J18365925
3C 454.3
PKS 1502106
85(Some) Important Questions to Answer
- How large is the positron fraction at HE (PAMELA)
- Identifies the nature of sources of primary
positrons - If SNRs are the sources of primary positrons,
this should also affect secondary nuclei - Measure the secondary nuclei (PAMELA, CREAM)
- How typical for the local Galactic environment is
the observed Fermi/LAT spectrum - If this is the typical spectrum than the sources
of primary positrons are distributed in the
Galaxy (could be pulsars, SNRs, or DM) - If this spectrum is peculiar than there is a
local source or sources of primary positrons - The answer is in the diffuse gamma-ray emission
(Fermi/LAT) - Dark matter vs Astrophysical source
- Distribution of the IC emission at HE (Fermi/LAT)
- WE HAVE ALL NECESSARY INSTRUMENTS IN PLACE (in
the orbit) TO ANSWER THESE QUESTIONS
86More data expected!
- ATIC and CREAM
- Elemental abundances up to 1015 eV
- PAMELA
- Absolute positron flux
- More on absolute antiproton flux
- Electrons
- Light nuclei
- Fermi Large Area Telescope
- Electrons gt1 TeV
- Diffuse emission (Galactic and extragalactic)
- Keep tuned
- A probe of the CR electron spectrum from the
solar surface to Saturns orbit at 10 AU - A probe of the local interstellar CR proton
spectrum - AMS will it fly in 2010?