AIA DEM Workshop @ NAOJ
Mark Cheung, LMSAL x NAOJ!
cheung@lmsal.com!

Schedule••

Monday March 28th•

10 am - 12 noon: Lecture (interactive)•

1pm onward: Data analysis / lab exercises

Tuesday March 29th Schedule is flexible•

BUT 1pm to 2:30pm Shimojo-san ALMA Solar Town Hall•

Data analysis / lab exercises•

Show and tell•“Office hours”

Outline
1. Aims
2. What do we see in the EUV channels of AIA?
3. Introduction to the problem of Differential Emission
Measure (DEM) Inversions
4. Description of sparse solver for AIA data
5. Tutorial on how to use the AIA_SPARSE_EM package
6. Practice exercises
7. Example science problems to tackle3Outline
1. Aims
2. What do we see in the EUV channels of AIA?
3. Introduction to the problem of Differential Emission
Measure (DEM) Inversions
4. Description of sparse solver for AIA data
5. Tutorial on how to use the AIA_SPARSE_EM package
6. Practice exercises
7. Example science problems to tackle
4

Aims
By the end of this workshop, I hope you will have someanswers to the following:•

How to find information about AIA data analysis?

•

What do SDO/AIA EUV images show?

•

What is the differential emission measure (DEM)
problem?

•

How do we perform DEM inversions on SDO/AIA
data?

•

How to perform quantitative investigations of thephysical evolution of the solar corona?
5

A taste of what we can do with AIA: Below is adifferential emission measure inversion.

6

Useful links / documentation
SDO User Guide!
http://www.lmsal.com/sdouserguide.html(Does not yet include how to do DEMs)!

AIA DEM package (in beta)!
http://tinyurl.com/aiadem!

Documentation about AIA!
https://www.lmsal.com/sdodocs/• AIA Instrument Paper (Lemen et al. 2012)• Initial Calibration of AIA (Boerner et al. 2012)!

CHIANTI User Guide!
$SSW/packages/chianti/doc/cug.pdf

SDO User Guide!

1 Introduction
5.2 Querying the HER
1.1 Synopsis of the SDO Mission
5.2.1 Basic Query Syntax
1.2 About this Guide
5.2.2 Adding Filters to Queries
1.3 Data Products from AIA
5.3 Retrieving Data from the JSOC
1.4 Data Products from EVE
5.3.1 Getting Dataseries Names
5.3.2 Selecting Records
1.5 Data Products from HMI
5.3.3 Exporting Data
1.6 SDO Mission Data Policy
5.3.4 Exporting Specific Segments
2 How to Browse SDO Data
5.4 Requesting a Cutout Using SSW
2.1 “The Sun Now”
5.5 Using the SSW Interface to the VSO
2.2 “Sun In Time”
5.5.1 Querying the VSO
2.3 SolarMonitor
5.5.2 Downloading Data from the VSO
2.4 The Helioviewer Project
5.6 Uncompressing the Image Data
3 How to Find SDO Data
5.6.1 Using read_sdo.pro
3.1 Heliophysics Events Knowledgebase
5.6.2 Using read_sdo.pro with a Shared Library
3.2 iSolSearch
6 How to Process AIA Data
6.1 Using aia_prep.pro to Align and Derotate AIA Data
4 How to Get AIA and HMI Data
6.1.1 Aligning Full-Disk Images
4.1 How Data are Organized at the JSOC
6.1.2 Aligning Cropped Images and Cutouts
4.2 The lookdata Tool
6.2 Co-Aligning HMI Data with AIA Data
4.2.1 Getting Dataseries Names
6.3 Despiking and/or Respiking AIA Data
4.2.2 Selecting Records
6.4 Point Spread Function Deconvolution
4.2.3 Exporting Data
6.5 Getting AIA Filter Response Data
4.2.4 JSOC Image Timing Details
6.6 Miscellaneous AIA Tasks
4.3 The Cutout Service
6.6.1 Using AIA Standard Scaling and Colors
4.4 The Virtual Solar Observatory
6.6.2 Creating AIA Tri-Color Images
4.5 Synoptic Data
6.6.3 Plotting an AIA Light Curve
6.6.4 Interfacing with plot_map.pro
4.6 Uncompressing the Data
(Dom Zarro’s Mapping Software)
4.6.1 Compiling and Installing imcopy
6.6.5 Checking the QUALITY Keyword
5 How to Use SolarSoft with SDO Data
7 Frequently Asked Questions
5.1 Installing or Upgrading SolarSoft
8 Lists of Useful Links
5.1.1 Doing a Clean Install
9 Version History of this Guide
5.1.2 Updating Packages for an Existing Installation
10 Contributors to this Guide
5.1.3 Using a Proxy Server with SSW

http://www.lmsal.com/sdouserguide.html

Outline
1. Aims
2. What do we see in the EUV channels of AIA?
3. Introduction to the problem of Differential Emission
Measure (DEM) Inversions
4. Description of sparse solver for AIA data
5. Tutorial on how to use the AIA_SPARSE_EM package
6. Practice exercises
7. Example science problems to tackle
9What do we see in AIA/EUV images?

10What do EUV images from SDO/AIA show?
;
;!
;
;!
;
;!
;
;
;
;
;!
;
;
;
;
;
;
;
;

Name: aia_get_response!
Purpose: return SDO/AIA instrument response data structures!
Input Parameters:!
Output:!
function returns AIA wavelength response, spectral model, or temperature response!
Descriptions of the structures returned by this routine are available via: !
http://sohowww.nascom.nasa.gov/solarsoft/sdo/aia/response/README.txt!
Calling Examples:!
IDL> effarea=aia_get_response(/area, /dn)
; return per wavelength effective area !
IDL> effarea=aia_get_response(/area,/full,/dn) ; same with component details (filter/ccd...)!
IDL> emiss=aia_get_response(/emiss, /full); return default CHIANTI model with line list!
IDL> tresp=aia_get_response(/temp,/dn,/evenorm) ; return temperature response including!
; constant scale factors to give good!
; overall agreement with SDO/EVE!

See section 6.5 of SDO User Guide for extensive explanation.
11

From aia_get_response.pro in SSW
12

Fe X
Fe IX
Fe IX
Cont.

190.04189.94188.50

6.05
5.85
5.85

0.06
0.06
0.07
-

0.04
-

0.05

0.04

211 Å

Cr IX
Ca XVIFe XVIIFe XIVFe XIIIFe XIIIFe XIII
Fe XI
Fe X
Ni XI
Cont.

210.61
208.60204.67211.32202.04203.83209.62209.78207.45207.92

5.95
6.7
6.6
6.3
6.25
6.25
6.25
6.15
6.05
6.1

0.07
0.11
0.05
0.03
0.08

0.13
0.05
0.05
0.12
0.03
0.04

0.39
0.07
0.05
0.07

0.09
0.07
0.12
0.41

304 Å

He II
He II
Ca XVIII
Si XI
Cont.

303.786
303.781302.19303.33

4.7
4.7
6.85
6.2

0.33
0.66
-

0.32
0.65
-

0.27
0.54
0.11
-

0.29
0.58
0.05
-

335 Å

Al X
Mg VIII
Mg VIII
Si IX
Si VIIIFe XVIFe XIV
Fe X
Cont.

332.79335.23338.98341.95319.84335.41
334.18184.54

6.1
5.9
5.9
6.05
5.95
6.45
6.3
6.05

0.05
0.11
0.11
0.03
0.04
0.13
0.08

0.11
0.06
0.06
0.03
0.04
0.15
0.05

0.86
0.04
-

0.81
0.06

coronal hole, quiet Sun, active region and flare plasma

communication), with the exception of the 171 Å and 335 Åchannels for which updated versions were used obtained from
Table 1. Predicted AIA count rates.
Solarsoft(12 July, 2010).
Ion

3. Results

λ
Å

T pa
K

Fraction of total emission
CH
QS
AR
FL

Mg VIIIspectral
94.07
5.9 contribute
0.03 - more
- than 3%
Table941Ålists thoselines which
Fe XX
7.0 for- CH, QS,- AR- and flare
0.10
to the total emissionin 93.78
each channel
XVIII
6.85 -contribution
0.85
plasma. AlsoFeincluded
is93.93
the fractional
of0.74the conFe X for any94.01
6.05 the0.63
0.72 contributes
0.05 tinuum emission
case wherecontinuum
93.47
0.04 - Synthetic
more than 3%FetoVIIIthe total
emission5.6in a channel.
specFe VIII
93.62
5.6
0.05 tra for each ofCont.
the channels are displayed0.11
in Fig.0.12
2 - 8.0.17
For every
channel each spectrum has been divided by the peak intensityof theWeakerspectrabeen -scaled-by
131strongest
Å O VIspectrum.
129.875.45
0.04have
0.05
Fe XXIII
0.07
factors indicatedin each132.91
figure. 7.15 XXI- 93.930.83
The 94 ÅFechannel
is128.75
expected17.05
to observethe- Fe XVIII
Å
Fe
VIII
130.94
5.6
0.30
0.25
0.09
∼ 6.85] in flaring regions. For both AR and flaren line [log T[K]Fe VIII
131.24
5.6
0.39 0.33 0.13 plasma (see Fig.
2)
thedominant
contribution
comes from theCont.
0.11 0.20 0.54 0.04
Fe XVIII 93.93 Å line. However, for the CH and QS spectrum (see
Å dominant
Ni XIV contribution
171.37
6.35
Fig. 1712) thecomes- from -the Fe 0.04
X 94.01- Å, line [log T ∼ Fe
X
174.536.05 0.03 6.05].
Fe IX
171.075.85 0.95 0.92 1 0.80 0.54
In flaringCont.
regions the 131 Å channel- is expectedtoobserve
0.23
e
the Fe XX 132.84 Å and Fe XXIII 132.91 Å lines. However, forn
193 Åspectrum
O V (see 192.90
0.03 contribution- comesthe flareFig. 3) the5.35
dominant
M
Ca 128.75
XVII
192.85
6.75- combined
0.08
from
the
Fe
XXI
Å
line
[log
T
∼
7.05].
Theconm
Ca XIV
193.87
6.55 0.04 tribution
of
the
Fe
XX
132.84
Å
and
Fe
XXIII
132.91Ålines
is less
Fe XXIV
192.03
7.25 0.81
thanten
percent
of
the
total
emission.
ForCHobservationsr
Fe XII
195.12
6.2
0.08 0.18 0.17 -the1
to be 6.2
dominated
Fe VIII0.17
lines [log
Fe is
XIIexpected
193.510.09 by0.19
c 131 Å channel
the dominant
contribution
Fe XIIour simulations,
192.396.2
0.04 0.09
0.08 -forl T ∼ 5.6]. FromFe XI (see188.23
6.15come
0.09from0.10
Fig. 3) does
Fe VIII0.04
lines, -but
- CH and QS plasma
Fe XI contribution
192.83 from
6.15
0.05 0.06- For -thecontinuum
emission.
. with a significant
Fe XI
188.306.15 0.04 0.04 e
Fe X
190.046.05 0.06 0.04 Å2
Fe IX
189.945.85 0.06 m
Fe IX
188.505.85 0.07 Cont.

-

-

0.05

From O’Dwyer et al. 2010

The count rates are normalised for each channel. Coronal hole(CH), quiet Sun (QS), active region (AR) and flare (FL) plasma.
a
T p corresponds to the log of the temperature of maximum abundance.

0.04 13

Outline
1. Aims
2. What do we see in the EUV channels of AIA?
3. Introduction to the problem of Differential Emission
Measure (DEM) Inversions
4. Description of sparse solver for AIA data
5. Tutorial on how to use the AIA_SPARSE_EM package
6. Practice exercises
7. Example science problems to tackle
14

Statement of the Problem

The Atmospheric Imaging Assembly (AIA, Lemen et al. 2012;
Boerner et al. 2012) instrument onboard NASA’s Solar
Dynamics Observatory (SDO, Pesnell et al. 2012) is a suite offour normal-incidence reflecting telescopes that image the Sun
in Coronalseven EUV
channels,channels and Arcades
one visible
Mass
Sourcetwo
for UV
Post-Eruption
wavelength channel.
!
4. EMISSION MEASURE DETERMINATION USING AIA
The aim of this and many other studies is to extract thermal
4.1. theStatementof theinformation about
Sun’s opticallythinproblem
corona using the EUV
Thephysicalcalibratedproperties(i.e. dark-subtracted,
flat-fielded
can be observations.
related to theof opticallythin coronal pand exposure time normalized) count rate yi in the i-th EUV
channel is related to the thermal distribution of coronal plasma
Z
1
by:

yi =

Ki (T ) DEM(T )dT,0where Ki(T) is the temperature response function (see nexttime-normalized
pixel
value
in
the
i-th
AIA
channel
(in
units
of
D
slide):
R1

function (in units of DN cm5 s 151 pixel 1 ) and DEM(T ) =

n

1
ed
to
the
physicalproperties
ofopticallythincoronal plasma as an integ
4.1.
Statementoftheproblem(T ) DEM(T )dT,yi =Ki (T ) DEM(T )dT,
ons can be related to the physical propertiesof optically thin coronal plas0Z 11i-th AIA 1
theexposuretime-normalized
pixel
value
in
the
channel(in)uihe i-thyiAIA
channel
(in
units
of
DN
s
pixel
),
K
(T
i
=
Ki (T ) DEM(TZ 1)dT, R 1 5 1
1
ature response
function
(in
units
of
DN
cm
s
pixel
)
and
DEM(T
)
1
2
0
pixel and) and DEM(T= ) DEM(T
n
(T
)dz,
where
DEM(
y
K
(T
)dT,
e
5i = 1)dTi0easure
(in
units
of
cm
K
)
of
the
plasmaalong
aline-of-sight.
n
(T
)
e
1
1
0
1pixel value in the i-th AIA channel (in units of DN s
ized
pixel
),
K
K
) whereoftemperature
theDEM(T)plasmaalong
aline-of-sight.
n
(T
)
is
the
elec
R
a certain
T
.
Let
the
temperature
range
bedividedintonne
is1 the differential
emission measuree2 units of
1 (in
5
1
nure
units
of
DN
s
pixel
)
and
DEM(T
)
=
n
(T
)dz
is
the
di↵
-5
-1 cm
time-normalized
pixel
value
in
the
i-th
AIA
channel
(in
units
of
DN
e
cm
K
)
of
plasmaalong
theline-of-sight.
n
(T)
is
the
electron0e
et1 the temperature range 5ben 1divided
intonneighboring
R1 2
Z
1
T
+
T
X
) of
the
plasmaalong
aDN
line-of-sight.
nje (TT.)j)istheelectronnumber
onse
function
(in
units
cm
s pixelandDEM(T= 0 nede(Tnumberdensityofof
plasmaat temperatureThe
challengeis) toy=K(T)DEM(T
)dT,
5
1
i
i
et
the
temperature
range
bedividedintonneighboring
bins,
sothe
that
givenplasma
a set ofalong
EUV measurements
y. ne (T ) isunits ofsolve
cm forKDEM(T)) of thealine-of-sight.
elect
Tj
j=1
mperatureT. Let the temperature range be divided into n neighboring b
ZnTj
j + TjLet
range be divided into n neighboring bins,
XtheTtemperature
j-thyi temperature
bin K
has
range
T
2
[T
,
T
+
T
).
The
aim
is
to
use
A
Z
j
j
j
n
=
(T
)DEM(T
)dT,
T
+
T
so
that:
i
j
j
X
K
(T
)DEM(T)dT,ieasure !distribution(eitherindi↵erential
or
integralform).
Tjj=1
yi =
Ki (T )DEM(T )dT,
that Ki!(T ) = Kij is piecewiseTj constant in each j-th temperature bin, so
j=1
ange T where
2 [Tj ,theTj +to nuseAIAT measurements
to retrj ). The aim
j-th Ttemperaturebinis has
range
∈
[T
,T
+∆T
).
j
j
j
X
in
di↵erential
or
integralform).
raturehas range
2to[Tpiecewise, Tj +AIA
T
).
The
aim
is
to
use
AIA
measurTj ).binTheaimmeasurements
toretrievj use
jAssumingthat KTiis(T)isconstantin
each
temperatureyi =Kij EMj , where
ise constant
in
each
j-th
temperaturebin,soribution(eitherindi↵erential
or
integralform).
tegralform).
bin, we have:j=1= Kij is piecewise
constant
in
each
j-th
temperature
bin,
so
Z
n
Tj + Tjh j-th temperaturebin, soX
n EM =yi =Kij EMj , where
DEM(T )dT.
Xjyi =Kij EMj ,Twhere
j
j=1
16
Z

Statement of the Problem

that Ki (T ) = Kij is piecewise constant in each j-th temperature bin, so
range T 2 [Tj , Tj + Tj ). The aim is to nuse AIA measurements to retXin di↵erential or integral form).
yi =Kij EMj , where
wise constant in each j-th temperature bin, so
j=1
Z Tj + TjnXyi =Kij EMj , where EMj =DEM(T )dT.

Statement of the Problem: y = Kx
j=1

Tj

Za Tmatrix
The
above
is
j + Tjequation of the form y = Kx, whereori knowledge about
the form of DEM(T ) over the range [Tj , Tj + Tj ), it
EM
=
DEM(T
)dT.
• KK
is(T
an
m
x piecewisen response
matrix*,
with
each
row
corresponding to temperathejn that)beconstant
anddefinea
DEM-weightedi
Tj response function of one AIA channel
temperature
this assumption
since such information is not available.
•writteny is an m-tuple
correspondingof AIA count (rates), and
n
be
in
the
followingform
rm of DEM(T ) over the range [Tj , Tj + Tj ), it is in principle possible to• x is an n-tuple with components EMj.
constant and define a DEM-weighted temperature response function. H!
y = Kx,information is not available.
The He II line in the 304 Å channel is not well-modeled by CHIANTI
ng forms an (Warren,m ⇥ n matrix
with157,
components, y is not
an usedm-tuplecorrespondingt
2005, ApJ
147) so it K
is ij
usuallyfor DEM
analysis.
m =So6mwhen
using
94, we
131,
171,
193,
2116 temperatureand 335 Å channels)= 6 for
AIA. theUsually
want
more
than
bins. For m andy = Kx,EMj .< n, the matrix equation y = Kx represents an underdetermined system.
<n !(i.e. moretemperaturebins), Eq.to(5)is underdetermined.
mponents
Kij , ythan
is an6 m-tuple
corresponding
measurements
by the A
measureinversions
and
thus
far,
researchers
haveattempted
totacklet*Matrix
elements
dependson
basis
functionsusedforcomputingthe131, 171, 193, 211 and 335 Å channels) and x is an n-tuple with com
That
is, the DEM solution is chosen to be one such that
integral

✓
◆
17
m
2
mperature bins), Eq. (5) is underdetermined.
This
is
a
well-knownproXUsual Approach: 𝛘-squared Minimization!
Function to minimize: | y - Kx |2 or |(y - Kx)/𝞂|2!
Basically, minimize difference between observed and predicted
counts. The benefits of a least-squares approach is that it leads toEuler-Lagrange equations that can be used to seek (global or local)
minima.
!
For an overdetermined system, we know no single model will fit all
the data. So 𝛘-squared minimization is ideal. However, forunderdetermined systems such an approach can be subject to theperils of overfitting.
!
Usual way to get around this:Parameterization: e.g. Guennou et al (2012a,b),
xrt_dem_iterative2.pro (M. Weber in SSW, see also Cheng et al
2012)
Regularization: e.g. Hannah & Kontar (2012), Plowman et al. (2013)18

Outline
1. Aims
2. What do we see in the EUV channels of AIA?
3. Introduction to the problem of Differential Emission
Measure (DEM) Inversions
4. Description of sparse solver for AIA data
5. Tutorial on how to use the AIA_SPARSE_EM package
6. Practice exercises
7. Example science problems to tackle
19

K~x max(~y ~⌘ , 0),0.
Cheung et al. 2015: The ~xSparse
Solution

13) andprovidesometoleranceforthesolutiontodeviatef
We(14)
address
the inverse
problem usinganapproachdifferentpose ⌘jthan= 2chi-squared
yj . Other
than the The
problemstated satisfyingby (13) - (14minimization.
set ofas
solutionsthe underdetermined matrix equation y = Kx lies in an affine
are imposed.
n#subspace of R .
such that:minimize

We pick the solution x within this subspacenX

xj subject to K~x = ~y , ~x0.

jThe linear program above finds a solution that minimizes the
L1-norm of the nsolution vector (c.f. Candes 2006). This is notXchi-squared minimization.
minimizex
subject
to
K~
x

~
y
+
~
⌘
,
,
~
x
0,
j
If K is the response function sampled by Dirac Delta functions atjspecific temperatures,
this is equivalent to minimizing the total
EM. If other basis functions are K~
used,
there yis no
simplex max(~~⌘ , 0).
corresponding physical interpretation.
!
!
20

The Sparse Solution
We are unaware of physical principles pertaining tocoronal plasma that motivate the optimization problem
posed above. However, this choice has some important
benefits:
1)It does not overfit (consistent with the principle ofparsimony, i.e. Ockham’s Razor).
2)It ensures positivity of the solution (if solutions
exist).
3)It is an L1-norm minimization problem, so we can use
standard techniques from compressed sensing (c.f.
Candes & Tao 2006).
BTW the L1-norm of a vector x = 𝝨 |xi|
4) Speed: O(104) solutions / sec with single IDL thread.
213) and (14)providesometoleranceforthesolutiontodeviatefromsatpand (14) provide some tolerance for the solution to deviate
itions(13)ose ⌘j = 2 yj . Otherp than the problem as stated by (13) - (14), no ot
we⌘j = 2 yj . Other than the problem as stated by (13) - (14are choose
imposed.
M curve) are imposed. Handling noisenXnXminimizex
subject
to
K~
x
=
~
y
,
~
x
0.
j
In practice,minimize
measurement
uncertainties
imply
that0.thexsubject
to
K~
x
=
~
y
,
~
x
j
j
equality y = Kx may jnot be satisfied. So our method
solves the followed modified linear program:nXnXminimize
xj subject to K~x  ~y + ~⌘ ,minimize
xj subject to K~x  ~y + ~⌘ ,

j

j

~x

0, K~x max(~y ~⌘ , 0).
~x 0, K~x max(~y~⌘ , 0).

The vector η is a measure of the uncertainty in the
count rate and provides tolerance for the predicted
counts (Kx) to deviate from the observed values (y). To
enforce positive counts the lower bound is set to max(yη, 0).

22𝛘-squared
minimization

Sparse inversion

94

94

335

131

335

131

211

171

vs.
211

171
193

193

Crossmarks are observed counts
Subspace of allowed solutions in our method
23

Basis Functions for DEM
APPENDIXQUADRATURE SCHEME

Let i = 1, 2, ..., m denote the index over a set of wavelength band channels and/or line spectra. Let the DEM function
be written in terms of a set of positive semidefinite basis functions {bj (log T ) 0 | k = 1, 2, ..., l}, viz.
DEM(log T ) =lXbk (log T )xk ,

(A1)

k=1

with quadrature coefficients xk

0. Approximating the integrals in equation (1) as sums in log T space, we have
yi =

n X
l
X

Kij Bjk xk

log T,(A2)

j=1 k=1

where j = 1, 2, ..., n is the index over temperature bins, Kij = Ki (log Tj ) and Bjk = bk (log Tj ). The response matrix
K = (Kij ) has dimensions m ⇥ n. The basis matrix B = (Bjk ) has dimensions n ⇥ l, with the k-th column vector
corresponding to the k-th basis function bk (log Tj ). Defining the dictionary matrix D = KB, the set of integral
equations (1) can be written in matrix form:~y = D~x,(A3)

Temperature Bins

In	 practice	 we	 solve	 this	 ➜

where the sought-after solution vector ~x is an l-tuple with components xk log T (k = 1, 2, ..., l). When the numberof basis functions exceeds the number of image channels (i.e. l > m), the linear system Eq. (A3) is underdetermined.
For the results in this paper, we use an equidistant grid in log T with log T = 0.1 ranging from log T = 5.5 to 7.5(i.e. n = 21). Over this temperature grid, the set of Dirac-delta basis functions {bDirac| k = 1, ..., n} iskBasis matrix B

bDirac(log Tj ) = 1, if log Tj = log Tk ,k
= 0, otherwise.

(A4)
(A5)

Recall that the basis matrix B consists of column vectors corresponding to basis functions. So for the set of Dirac-delta
functions BDirac = I (the identity matrix).
In addition to Dirac-delta functions, we also use basis functions consisting of truncated Gaussians. Each Gaussianfunction of width a generates a set of basis functions {bak | k = 1, ..., n}, where(log Tj log Tk )2abk (log Tj ) = exp, if | log Tj log Tk |  1.8a.
(A6)a2= 0, otherwise.
(A7)
The Gaussian basis functions are truncated (i.e. set to zero) for values of log Tj outside the temperature grid usedafor inversions. The corresponding basis matrix for this set
24 is denoted B . Di↵erent sets of basis functions can be
combined by concatenating their associated basis matrices. For the inversions shown here, we use the combined basis

y
94 DN/s
131 DN/s
171 DN/s
193 DN/s
211 DN/s
335 DN/s

Geometrical Interpretation of Problem
D = KB=We seek to build a column vector y by linear
combinations of columns of the matrix D=KB.
xj are the coefficients of the linear
combination. More columns of KB fromwhich to chose than components of y !
=> underdetermined problem.!
!

Do “basis pursuit” by minimizing L1 norm
of x.
25

x

Application to real AIA data

26

Warren, Brooks & Winebarger (2011)

Side benefit: Image Denoisingy (AIA 131 Level 1.5)

Dx# (from inversion)

Outline
1. Aims
2. What do we see in the EUV channels of AIA?
3. Introduction to the problem of Differential Emission
Measure (DEM) Inversions
4. Description of sparse solver for AIA data
5. Tutorial on how to use the AIA_SPARSE_EM package
6. Practice exercises
7. Example science problems to tackle
30

How to use the Sparse DEM code
http://tinyurl.com/aiadem•

Instructions in the readme file•

Download *genx files, which contain sample AIA
6-channel data.

•

Download aia_sparse_em_init.pro & exercise1.pro➡ .compile aia_sparse_em_init➡ .r exercise1 ; Example of DEM inversion➡ .r exercise2 ; Example of DEM sensitivity to noise
31

; Author: Mark Cheung; Purpose: Sample script for running sparse DEM inversion (see Cheung; et al. 2015) for details.
; Revision history: 2015-10-20 First version;
2015-10-23 Added note about lgT axis!
files = file_list('.','*genx')
aiadatafile = files[0]; Restore some prepared AIA level 1.5 data (i.e. already been aia_preped)
restgen, file=aiadatafile, struct=s!

!
; Initialize solver.
; This step builds up the response functions (which are; time-dependent) and the basis functions.
; If running inversions over a set of data spanning over a few hours; or even days, it is not necessary to re-initialize!
;IF (0) THEN BEGIN; As discussed in the Appendix of Cheung et al. (2015), the inversion; can push some EM into temperature bins if the lgT axis goes below; logT ~ 5.7 (see Fig. 18). It is suggested the user check the; dependence of the solution by varying lgTmin.
lgTmin = 5.7 ; minimum for lgT axis for inversion
dlgT = 0.1 ; width of lgT bin
nlgT = 21 ; number of lgT bins
aia_sparse_em_init, timedepend = s.oindex[0].date_obs, /evenorm, $!
!
use_lgtaxis=findgen(nlgT)*dlgT+lgTmin!
lgtaxis = aia_sparse_em_lgtaxis(); ENDIF!
; We will use the data in s.img to invert. s.img is a stack of level; 1.5, exposure time normalized AIA pixel values.
exptimestr = '['+strjoin(string(s.oindex.exptime,format='(F8.5)'),',')+']'!

!
; This defines the tolerance function for the inversion.
; y denotes the values in DN/s/pixel (i.e. level 1.5, exposure time
; normalized)!
; aia_bp_estimate_error gives the estimated uncertainty for a given; DN / pixel for a given AIA channel. Since y contains DN/pixel/s,; we have to multiply by exposure time first to pass; aia_bp_estimate_error. Then we will divide the output of; aia_bp_estimate_error by the exposure time again.
; If one wants to include uncertainties in atomic data, suggest to add; the /temp keyword to aia_bp_estimate_error()
tolfunc = 'aia_bp_estimate_error(y*'+exptimestr+', [94,131,171,193,211,335], $!
num_images=‘+strtrim(string(s.binning^2,format='(I4)'),2)+')/'+exptimestr!

!
; Do DEM solve.
; Note!!! The solver assumes the third dimension is arranged according to[94,131,171,193,211,335]
aia_sparse_em_solve, s.img, tolfunc=tolfunc, tolfac=1.4, oem=emcube, status=status,
coeff=coeffOutput of exercise1.pro

status[Nx, Ny] - Indicating where a DEM solution was found, 0 is yes
coeff[Nx, Ny, 84] - Coefficients of basis functions used to express
DEM; i.e. the x in equation y = KBx.
emcube[Nx, Ny, 21] - DEM solution. i.e. Bx.
34

Interactively inspect DEM profiles•

aia_sparse_dem_inspect, coeff, emcube, status, /ylog

35

36

Outline
1. Aims
2. What do we see in the EUV channels of AIA?
3. Introduction to the problem of Differential Emission
Measure (DEM) Inversions
4. Description of sparse solver for AIA data
5. Tutorial on how to use the AIA_SPARSE_EM package
6. Practice exercises
7. Example science problems to tackle37

Active Region Thermal Structure

Warren, Winebarger & Brooks (2012) devised a method toseparate the ‘warm’ and ‘hot’ components of emission fromthe AIA 94 Å channel. They use this empirical fit:“for this value of f [=0.3], the coefficients to the polynomial fit
are −7.31×10−2, 9.75 × 10−1, 9.90×10−2, and −2.84×10−3.”3839Workshop Practice Exercise 1:
AR 11190

From: Warren, Winebarger & Brooks (2012)

Go to http://www.lmsal.com/~cheung/AIA/ar11190/
Download a genx file and plot the DEM distribution inregions marked by the green boxes.
!

How do the DEMs in these boxes evolve?
What about the DEM averaged over the AR? (Take carewith pixels that have no DEM solutions).
40

Workshop Practice Exercise 2:
AR 11190

From: Warren, Winebarger & Brooks (2012)

Go to http://www.lmsal.com/~cheung/AIA/11190/
Download a genx fileStarting from the DEM solution, generate AIA 94 images
separately for the ‘warm’ and ‘hot’ components.
!

Hint: use
PRO AIA_SPARSE_EM_EM2IMAGE, em, image=image,
status=status
41

Workshop Practice Exercise 3:
Examine impact of choice of basis functions•

By default, aia_sparse_em_init uses 4 sets of basis functions(Dirac + 3 sets of truncated gaussians).

•

Default keyword is bases_sigmas = [0,0.1,0.2,0.6]•

Test how choosing other basis sets changes DEM results. e.g.
•

Only Dirac Delta functions bases_sigmas = [0]•Dirac + Gaussians of width 0.1 bases_sigmas = [0,0.1]•

Dirac + Gaussians of width 0.1, 0.2 and 0.3 bases_sigmas= [0,0.1,0.2,0.3]

Workshop Practice Exercise 4:
AIA + XRT DEM inversions•

Joint AIA-XRT DEM inversions•

Download aiaxrt.pro from http://tinyurl.com/aiadem••.compile aia_sparse_em_init•.r aiaxrt

This script solves uses sample AIA data to solve
for a DEM cube. Then it synthesizes AIA and XRT
images. Then it uses AIA+XRT for DEM inversion.

Outline
1. Aims
2. What do we see in the EUV channels of AIA?
3. Introduction to the problem of Differential Emission
Measure (DEM) Inversions
4. Description of sparse solver for AIA data
5. Tutorial on how to use the AIA_SPARSE_EM package
6. Practice exercises
7. Example science problems to tackle: 13:15pm today
44

Science Cases•

Thermal evolution of active regions fromemergence to decay•

Thermodynamic structure of reconnection
outflows•

From chromospheric evaporation to catastrophic
condensation

Science Case 1) Emerging Flux
DEM movieof the emergenceof AR 11726

Other panels:
EM in various log T bins

Lower right panel only
Greyscale:
Blos from HMI
Green: 6MK EM
Yellow/Red: 10 MK EM

DEM movieof the emergenceof AR 11158

Black pixels arewhere no solutions
were found4 January 1

Science Case 2) Reconnection Outflows
Krucker & Battaglia

Selected scientific studies of this flare:

10000.0

Y (arcsecs)-200

above-the-loop-top• Patsourakos, Vourlidas
thermalfootpoint
1000.0M7.7 flare on
2012-07-19
photons s-1 cm-2 keV-1

-180

-220

& Stenborg, 2013, ApJ, 764,
125: Prior confined eruption produced a pre-existingcoronal flux rope, which then erupted to give the M7.7
100.0
flare.

•10.0Wei Liu, Chen & Petrosian, 2013, ApJ, 767, 168:
Detailed timeline of sequence of events including timing-260
Figure 2. Time proﬁles of the three ﬂares: the top panels show the 1.0
timeof the WL ﬂuxof(summedarbitrary units) and non-thermal HXR
andproﬁles
propagation
EUV overandenhancement;
X-ray sources.
-240

ﬂux at 30–80 keV in red and blue, respectively,
with the linear GOES curve shown schematically in gray in the background for timing reference. Below, the apparent
AIA 193A-280
altitude of the peak of the WL emission
above
limb is shown in•
red.
TheLiu,black dottedlinesMNRAS,
give the FWHM
of a 1309:
Gaussian ﬁtSame
to the HMI
source time
above thelimb; the
Rui
2013,
434,
onset
for
RHESSI
6-8thekeVthin black line is the deconvolved source size derived by a rough deconvolutionof the observed size with the HMI PSF from Yeo et al. (2014). The blue points give
0.1
RHESSI 30-80 keV
HXR andmicrowave
RH data)the HXR
No errorbars are shown(Nobeyama
for the HMI centroids
as they bursts.
are dominated by systematic-300centroid positions with error bars from statistical uncertainties.
10100
uncertainties (∼70
900 km).
920
940
960
980 1000 1020

[keV]
• Krücker energy& Battaglia, 2014, ApJ, 780, 107: Ratio of
copy observations of the 2012 July 19 flare. On the left, the time profile of the non-thermal HXR flux at 30–80 keV
thermal protons (from AIA) to non-thermal electrons
schematically in gray in the background. The time interval used to reconstruct the RHESSI image shown in the
XR peak during the impulsive phase of the flare. The thermal emissions in(fromthe 6–8 keVrange (green
contours
at
RHESSI
HXR)
in loop-top
source is of order 1.

X (arcsecs)

on of the main flare loops also seen in the 193 Å AIA image. The non-thermal HXR emissions come from thet 30%, 50%, 70%, 90%), but also from above the main flare loop as outlinedby theCheng
30–80 keV&
contours.
Relative• Sun,
Ding,
2014, ApJ, 786, 73: Performed ae is shown in contours at 2%, 2.5%, 3%. The plot on the right gives the imaging spectroscopy results: The black
the combined coronal sources; the observed footpoint spectrum is given bycrosses.detailed
The green andred curves
are
very
DEM
(imaging)
analysis of this flare. They
onal sources, while the power-law fit to the footpoints is given in blue.
used xrt_iterative_dem2.pro (code by M. Weber, nonne journal.)

drastically increased
Atmospheric Imaging
rd the Solar Dynamic
best diagnostics ofservations of a GOES
aneously observed byevent provides many., Liu et al. 2013; Liulinear least square inversion with splines).

RHESSI images also reveal theexistence ofsources.
• Krücker
etfootpointal., 2015,
ApJ, 802, 19: HMI white light (WL)
The northern footpoint dominates over the southern footpoint,
enhancement
at footprint
co-incident with HXR
but this asymmetry is likely becausepart of the emission
fromthe flare ribbons is occulted by footpoint
the solar disk.source.
STEREO Interpretation:
EUVI
energy deposition by
images reveal that, if seen from
Earth, the solarelectrons
disk occultsis radiated away and does not
non-thermala larger part of the southern ribbon than the northern ribbon
cause
evaporation.
(Patsourakos et al. 2013; Liu et al.
2013);chromosphericindicating the weakersouthern footpoint could indeed be due to occultation. WhileTotal Emission Measure [cm-5]

Emission Measure-weighted log TTotal Emission Measure [cm-5]

Emission Measure-weighted log TTotal Emission Measure [cm-5]

Emission Measure-weighted log TTotal Emission Measure [cm-5]

Emission Measure-weighted log TShiota et al. (2005, ApJ,
634, 663):• 2.5D MHD simulationof of the eruption of apre-existing flux ropet r i g g e r e d b y fl u x
emergence.
• Similar scenario as
modeled by Chen &Shibata (2000, ApJ,
545, 524) but withfield-aligned thermal
conduction.
• Both temperature anddensity are initially
uniform (dimensionlessvalue of unity).

Shiota et al. (2005, ApJ,
634, 663):• 2.5D MHD simulation ofof the eruption of a preexisting flux rope
triggered by flux
emergence.
• Similar scenario as
modeled by Chen &Shibata (2000, ApJ,
545, 524) but with fieldaligned thermal
conduction.
• Both temperature anddensity are initially
uniform (dimensionlessvalue of unity).

Dashed contours: Total EM =1029 cm-5
Solid contours: Total EM =1030 cm-5
Chromospheric evaporation?

Downward mass pumping from!
reconnection outflow?

Influence of efficiency of thermal conduction / system size
𝜿0 = 0𝜿0 = 0.27
𝜿0 = 2.7
Log10 (T/2 MK)

Log10 (N/109)

Extension of
study byTakasao et al.
(ApJ, 2015,
805, 135)vz [170 km/s][18 s]

β = pg / pm = 0.08

Influence of efficiency of thermal conduction / system size
𝜿0 = 0𝜿0 = 0.27
𝜿0 = 2.7
Log10 (T/2 MK)

1. Higher compression region(“blob”) for stronger thermal
conduction
Log10 (N/109)

vz [170 km/s][18 s]

β = pg / pm = 0.08

Influence of efficiency of thermal conduction / system size
𝜿0 = 0𝜿0 = 0.27
𝜿0 = 2.7
Log10 (T/2 MK)

2. Enhanced evaporation forstronger thermal conduction
Log10 (N/109)

vz [170 km/s][18 s]

β = pg / pm = 0.08

Long Duration
GOES C6.7 flareSome remarks• AIA differential emission measure inversions (e.g. using the method of

Cheung et al. 2015; http://tinyurl.com/aiadem) can be used toquantitatively study the thermal evolution of flare plasma.
• The efficiency of thermal conduction / system size impacts:

1.The density enhancement in deceleration region of the reconnection
outflow, and2.The amount of evaporated material.
• Observations can possibly help us test whether classical Spitzer thermal

conduction is appropriate or should be modified (e.g. Imada, Murakami &Watanabe, PoP, 2015, 22, 10; Wang et al., 2015, ApJL, 811, L13).
• Future work:• Should measure 1. and 2. in a systematic study of flares to test

sensitivity to system size / efficiency of thermal conduction.
• Use 3D MHD simulations of flares for forward modeling of observables.

ographLab)

Science Case 3) ChromosphericEvaporation & Condensation
ObservationsodSpectra covering temperatures from 4,500 K to 10 MK
Images covering temperatures from 4,500 K to 65,000 K
Baseline: 5s for slit-jaw images, 1s for 6 spectral windows, rapid rasteringdPredicted Count Rates

at

IRIS is designed to probe
the chromosphere andtransition region, but it
also sees flare plasma(Fe XXI @10 MK).
!

What about thetemperature range inbetween?

ph!

Joint observations
between IRIS andAIA tofill the temperature gap.

Observing
FromModes
the IRIS poster @ iris.lmsal.com

At http://www.lmsal.com/data.html, !
go to ‘List of Flares Observed with IRIS’ !
compiled by Kathy Reeves.

20140917 19:26 C7.5 flare - behindthe limb, nice big loop of Fe
XXI visible, red shift in
Fe XXIIRIS SJI 1330: Diffuse, long-lived loops reaching into thecorona are generally attributed to Fe XXI 1354 Å emission.

IRIS SJI 1330: Likely Fe XXI 1354 Å emission.
C II 1336

O I 1356

Si IV 1403

Mg II k 2796

GOES X-ray

SJI Total DN/image

SJI 1330

IRIS SJI 1330: Likely Fe XXI 1354 Å emission.
C II 1336

O I 1356

Si IV 1403

Mg II k 2796

GOES X-ray

SJI Total DN/image

SJI 1330

Other panels:
EM in various log T bins
• Tell-tale signs ofchromosphericevaporation!
• Loops filled withplasma at 10 MK andabove!
• Loops cool to lower
log T bins!
• At time (~20:29 UT)
when plasma cools
down to log T/K ~ 5.8,
coronalcondensations in SJI
1330 begin to appear.

Lower right panel only
Greyscale:
Blos from HMI
Green: 6MK EM
Yellow/Red: 10 MK EM

Other panels:
EM in various log T bins
• Tell-tale signs ofchromosphericevaporation!
• Loops filled withplasma at 10 MK andabove!
• Loops cool to lower
log T bins!
• At time (~20:29 UT)
when plasma cools
down to log T/K ~ 5.8,
coronalcondensations in SJI
1330 begin to appear.

Lower right panel only
Greyscale:
Blos from HMI
Yellow/Green: 6MK EM
Red: 10 MK EM

Other panels:
EM in various log T bins
• Tell-tale signs ofchromosphericevaporation!
• Loops filled withplasma at 10 MK andabove!
• Loops cool to lower
log T bins!
• At time (~20:29 UT)
when plasma cools
down to log T/K ~ 5.8,
coronalcondensations in SJI
1330 begin to appear.

Lower right panel only
Greyscale:
Blos from HMI
Yellow/Green: 6MK EM
Red: 10 MK EM

IRIS SJI 1330: Coronal condensations appear at about sametime (~20:29 UT) as when AIA sees sub-MK plasma.