4. Hands-on tutorial of IRIS² in IDL.¶
This tutorial is intended to help you to understand the concepts behing the \(IRIS^{2}\) code step-by-step in a hands-on approach.
A detailed description including limitations, errors, function description, video tutorials, and how to use \(IRIS^{2}\) in IDL is given in |IRIS2| in SolarSoft.
A formal description of the concepts, errors, and limitations of \(IRIS^{2}\) can be found in Sainz Dalda et al. 2019.
The fundamental steps to make an inversion of IRIS Mg II h&k lines based on Representative profiles Inverted by STiC are:
Installation of IRIS SolarSoft branch, including IRIS database:
0.1 Set enviroment variables to start SolarSoft including IRIS SolarSoft branch and IRIS database.
Getting the IRIS Mg II h&k data properly calibrated (or any other IRIS data).
Inspection and selection of points of interest to be inverted:
2.1 Cropping the profiles in wavelength to avoid border effects.
Loading \(IRIS^{2}\) database:
3.1 Understanding the elements in \(IRIS^{2}\) database. 3.2 Interpolation of the RPs in the \(IRIS^{2}\) database to the observed wavelengths.
Looking for the closest inverted RP to the selected profiles:
4.1 Calculate Euclidean distance between the observed profiles and all the RPs in the \(IRIS^{2}\) database.
4.2 Look for the closest (mininum Euclidean distance) inverted RP to the observed profile.
4.3 Look for the closest (mininum Euclidean distance) PCA coefficients of the inverted RP to the PCA coefficients of the observed profile.
Displaying the inversion:
5.1 Plot observed profile (dashed line) and inverted RP (thick) 5.2 Plot the corresponding Representative Model Atmosphere (RMA)
Analazing, interpreting, and eventually enjoying the results.
4.1. Installation of minimal IRIS SolarSoft branch (including \(IRIS^{2}\))¶
To follow this tutorial you will need to have installed IRIS branch in your SolarSoft installation. A minimal installation of IDL SolarSoft including IRIS branch and the last version of IRIS2 (temporary until it is officially released) can be found here.
After you have downloaded the tar files and uncompressed them in your_target_directory, you should have the following directories there:
your_target_directory> ls
|
|__ data
| |__ 20141201_154438_3800008053 20160114_230409_3630008076
| |__ 20141201_172438_3800008053 20170605_170000_3662506342
| |__ 20151226_111400_3620008076 20180423_050928_3620106077
|
|__ iris2ssw
| |__ shared
|
|__ results
| |__ iris2model_20141201_154438_3800008053_raster_t000_r00000.sav
| |__ iris2model_20141201_172438_3800008053_raster_t000_r00000.sav
| |__ iris2model_20151226_111400_3620008076_raster_t000_r00000.sav
| |__ iris2model_20160114_230409_3630008076_raster_t000_r00000.sav
| |__ iris2model_20170605_170000_3662506342_raster_t000_multi.sav
| |__ iris2model_20180423_050928_3620106077_raster_t000_multi.sav
| |__ iris2model_20180423_050928_3620106077_raster_t000_r00003.sav
|
|__ ssw
| |__ gen iris offline packages site vobs
|
|__ sswdb
|__ iris
|__ iris2
where your_target_directory can be, e.g., ‘/Users/asainz/ciencia/IRIS-01’
The following step is to set the enviroment variables correctly to run this minimal IRIS SolarSoft installation. For that, you first need to know in what kind of shell session you are working by typing:
your_target_directory> echo $SHELL
If you are working under a csh or tcsh session, then you should do:
your_target_directory> setenv SSW your_target_directory/ssw
your_target_directory> setenv SSWDB your_target_directory/sswdb
your_target_directory> source SSW/gen/setup/setup.ssw
your_target_directory> setenv SSWDB your_target_directory/sswdb
If you are working under a bash session, then you should do:
your_target_directory> export SSW='your_target_directory/ssw'
your_target_directory> export SSWDB='your_target_directory/sswdb'
your_target_directory> $SSW/gen/setup/setup.ssw
your_target_directory> export SSWDB='your_target_directory/sswdb'
Note you may need to set the SSWDB variable twice. Now, you can start an IDL SolarSoft session by typing:
your_target_directory> sswidl
If that command (in fact, it is an alias) is not found, then you may try:
In csh or tcsh session:
your_target_directory> source $SSW/gen/setup/ssw_idlIn csh or tcsh session:
your_target_directory> $SSW/gen/setup/ssw_idl
If for any reason, sswidl does not work, we suggest you to include the definition of SSWDB and SSWDB in your .cshrc, .tcshrc, or .bashrc file.
These script shell files are executed any time a new shell session is started.
4.2. Verifying the minimal \(IRIS^{2}\) SolarSoft and the temporary \(IRIS^{2}\) code installation¶
For the sake of clarity, let’s create a variable called my_dir with your_target_directory, for instance:
IDL> my_dir = '/Users/asainz/IRIS-01/'
Let’s verify that the \(IRIS^{2}\) database is properly addressed:
IDL> dir_sswdb = getenv('SSWDB')
IDL> print, dir_sswdb
The output should be something similar to:
/Users/asainz/ciencia/IRIS-01/sswdb
Let’s include the last version of \(IRIS^{2}\) (temporary until \(IRIS^{2}\) is officially released) before the IRIS SolarSoft branhc:
IDL> !path = expand_path('+'+ my_dir)+':'+!path
Let’s verify the code will be properly found:
IDL> findpro, 'iris2'
Procedure iris2 found in directory /Users/asainz/ciencia/IRIS-01/iris2ssw/
Procedure iris2 found in directory /Users/asainz/ciencia/IRIS-01/ssw/iris/idl/lmsal/iris2/
IDL> findpro, 'heliocentric'
Procedure heliocentric found in directory /Users/asainz/ciencia/IRIS-01/iris2ssw/shared/
As you can see, iris2 function is now found in the iris2ssw directory first, then in the IRIS SolarSoft branch.
In addition, heliocentric, which is originally located in Yohkoh SolarSoft branch, it is found in a subdirectory called shared.
4.3. First \(IRIS^{2}\) inversion¶
Let’s look for the data:
IDL> data_dir = my_dir+'/data/'
IDL> files = find_files('iris_l2*raster*fits', data_dir+'*/')
IDL> dim_files = size(files)
IDL> for j = 0, dim_files[1]-1 do print, j,' - ', strmid(files[j], strpos(files[j],'iris_l2'),100)
You should expect something like this:
0 - iris_l2_20141201_154438_3800008053_raster_t000_r00000.fits
1 - iris_l2_20141201_172438_3800008053_raster_t000_r00000.fits
2 - iris_l2_20151226_111400_3620008076_raster_t000_r00000.fits
3 - iris_l2_20160114_230409_3630008076_raster_t000_r00000.fits
4 - iris_l2_20170605_170000_3662506342_raster_t000_r00093.fits
5 - iris_l2_20170605_170000_3662506342_raster_t000_r00094.fits
6 - iris_l2_20170605_170000_3662506342_raster_t000_r00095.fits
7 - iris_l2_20170605_170000_3662506342_raster_t000_r00096.fits
8 - iris_l2_20170605_170000_3662506342_raster_t000_r00097.fits
9 - iris_l2_20170605_170000_3662506342_raster_t000_r00098.fits
10 - iris_l2_20170605_170000_3662506342_raster_t000_r00099.fits
11 - iris_l2_20170605_170000_3662506342_raster_t000_r00100.fits
12 - iris_l2_20170605_170000_3662506342_raster_t000_r00101.fits
13 - iris_l2_20170605_170000_3662506342_raster_t000_r00102.fits
14 - iris_l2_20170605_170000_3662506342_raster_t000_r00103.fits
15 - iris_l2_20170605_170000_3662506342_raster_t000_r00104.fits
16 - iris_l2_20170605_170000_3662506342_raster_t000_r00105.fits
17 - iris_l2_20170605_170000_3662506342_raster_t000_r00106.fits
18 - iris_l2_20170605_170000_3662506342_raster_t000_r00107.fits
19 - iris_l2_20170605_170000_3662506342_raster_t000_r00108.fits
20 - iris_l2_20170605_170000_3662506342_raster_t000_r00109.fits
21 - iris_l2_20170605_170000_3662506342_raster_t000_r00110.fits
22 - iris_l2_20170605_170000_3662506342_raster_t000_r00111.fits
23 - iris_l2_20170605_170000_3662506342_raster_t000_r00112.fits
24 - iris_l2_20180423_050928_3620106077_raster_t000_r00003.fits
4.3.1. Calibration of IRIS Level 2 data¶
We need to calibrate the IRIS data to physical units. The function get_info_irisl2 allows to recover information and calibrated data from IRIS Level 2.
This code is an user-friendly wrapper of other functions developed by J.P. Wuelser and S.L. Freeland (iris_get_response), G. Vissers (iris_calib_spectrum), and A. Sainz Dalda (iris_flux_radcal).
We can get some information from an IRIS Level 2 files with:
IDL> data_iris = get_info_irisl2(files[10])
IDL> help, data_iris
** Structure <30a3408>, 22 tags, length=35408, data length=35308, refs=1:
FILE_NAME STRING '/Users/asainz/ciencia/IRIS-01/data/20170605_170000_3662506342/iris_l2_201706'...
INDEX STRUCT -> <Anonymous> Array[1]
LAMBDA_FILE DOUBLE Array[626]
DATE STRING '05-Jun-17 19:17:43'
POS_XY DOUBLE Array[2]
LONLAT DOUBLE Array[2]
MU DOUBLE 0.94271849
SPCSCL DOUBLE 0.025459999
SPXSCL DOUBLE 0.16635000
SPYSCL DOUBLE 1.9954402
DATE_IN_FILENAME
STRING '20170605_170000'
IRIS_OBS_CODE STRING '3662506342'
RASTER_INFO STRING '00099'
EXTENT_ARCSEC_ARCSEC
DOUBLE Array[4]
EXTENT_PX_PX FLOAT Array[4]
EXTENT_PX_ARCSEC
DOUBLE Array[4]
EXTENT_TIME_ARCSEC
DOUBLE Array[4]
EXTENT_TIME_PX DOUBLE Array[4]
EXTENT_OPT DOUBLE 129.42030
EXTENT_OPT_COORDS
STRING Array[2]
LIST_EXTENT_COORDS
OBJREF <ObjHeapVar4125(LIST)>
FACTOR STRUCT -> <Anonymous> Array[1]
Now, let’s inspect what is in factor structure:
IDL> help, data_iris.factor, /str
** Structure <30c6408>, 41 tags, length=24984, data length=24980, refs=2:
FILENAME STRING Array[1]
VER_RF2 STRING 'L12-2017-04-23'
BLD_VERS STRING 'V9R0X'
FLUX_UNITS STRING 'Hz (cgs): erg cm^-2 s^-1 sr^-1 Hz^-1'
INFO_SOURCE_FILE
STRING Array[1]
DN2PHOTON_1343 FLOAT 4.00000
DN2PHOTON_2814 FLOAT 18.0000
DN2PHOTON_2826 FLOAT 18.0000
DN2PHOTON_2832 FLOAT 18.0000
DN2PHOTON_C_II_1336
FLOAT 4.00000
DN2PHOTON_FE_XII_1349
FLOAT 4.00000
DN2PHOTON_MG_II_K_2796
FLOAT 18.0000
DN2PHOTON_O_I_1356
FLOAT 4.00000
DN2PHOTON_SI_IV_1403
FLOAT 4.00000
FACTOR_1343 DOUBLE Array[52]
FACTOR_2814 DOUBLE Array[146]
FACTOR_2826 DOUBLE Array[145]
FACTOR_2832 DOUBLE Array[113]
FACTOR_C_II_1336
DOUBLE Array[84]
FACTOR_FE_XII_1349
DOUBLE Array[52]
FACTOR_MG_II_K_2796
DOUBLE Array[626]
FACTOR_O_I_1356 DOUBLE Array[100]
This structure contains all the information needed to calibrate the data contained in that IRIS Level 2 file.
If you are interested in to get the data calibrated, and we definitively are, then you can use the keyword factor_and_data equals 1:
IDL> data_iris = get_info_irisl2(files[10], factor_and_data=1)
IDL> help, data_iris.factor, /str
** Structure <30dfc08>, 50 tags, length=73816728, data length=73816724, refs=2:
FILENAME STRING Array[1]
VER_RF2 STRING 'L12-2017-04-23'
BLD_VERS STRING 'V9R0X'
FLUX_UNITS STRING 'Hz (cgs): erg cm^-2 s^-1 sr^-1 Hz^-1'
INFO_SOURCE_FILE
STRING Array[1
---> DATA_1343 FLOAT Array[52, 778, 16]
---> DATA_2814 FLOAT Array[146, 778, 16]
---> DATA_2826 FLOAT Array[145, 778, 16]
---> DATA_2832 FLOAT Array[113, 778, 16]
---> DATA_C_II_1336 FLOAT Array[84, 778, 16]
---> DATA_FE_XII_1349
FLOAT Array[52, 778, 16]
---> DATA_MG_II_K_2796
FLOAT Array[626, 778, 16]
---> DATA_O_I_1356 FLOAT Array[100, 778, 16]
---> DATA_SI_IV_1403 FLOAT Array[164, 778, 16]
DN2PHOTON_1343 FLOAT 4.00000
DN2PHOTON_2814 FLOAT 18.0000
DN2PHOTON_2826 FLOAT 18.0000
DN2PHOTON_2832 FLOAT 18.0000
DN2PHOTON_C_II_1336
FLOAT 4.00000
DN2PHOTON_FE_XII_1349
FLOAT 4.00000
DN2PHOTON_MG_II_K_2796
FLOAT 18.0000
DN2PHOTON_O_I_1356
FLOAT 4.00000
DN2PHOTON_SI_IV_1403
FLOAT 4.00000
FACTOR_1343 DOUBLE Array[52]
FACTOR_2814 DOUBLE Array[146]
FACTOR_2826 DOUBLE Array[145]
FACTOR_2832 DOUBLE Array[113]
FACTOR_C_II_1336
DOUBLE Array[84]
FACTOR_FE_XII_1349
DOUBLE Array[52]
FACTOR_MG_II_K_2796
DOUBLE Array[626]
FACTOR_O_I_1356 DOUBLE Array[100]
FACTOR_SI_IV_1403
DOUBLE Array[164]
TEXP_ARR_1343 DOUBLE Array[16]
TEXP_ARR_2814 DOUBLE Array[16]
TEXP_ARR_2826 DOUBLE Array[16]
TEXP_ARR_2832 DOUBLE Array[16]
TEXP_ARR_C_II_1336
DOUBLE Array[16]
TEXP_ARR_FE_XII_1349
DOUBLE Array[16]
TEXP_ARR_MG_II_K_2796
DOUBLE Array[16]
TEXP_ARR_O_I_1356
DOUBLE Array[16]
TEXP_ARR_SI_IV_1403
DOUBLE Array[16]
---> WL_1343 DOUBLE Array[52]
---> WL_2814 DOUBLE Array[146]
---> WL_2826 DOUBLE Array[145]
---> WL_2832 DOUBLE Array[113]
---> WL_C_II_1336 DOUBLE Array[84]
---> WL_FE_XII_1349 DOUBLE Array[52]
---> WL_MG_II_K_2796 DOUBLE Array[626]
---> WL_O_I_1356 DOUBLE Array[100]
---> WL_SI_IV_1403 DOUBLE Array[164]
As you can see, now there are a set of variables called data_spectralrange, where spectralrange is the label name of the spectral range
observed by IRIS in that observation.
Because this tutorial and \(IRIS^{2}\) are focus in Mg II h&k lines, we can use the keyword mgii_only equals 1:
IDL> data_iris = get_info_irisl2(files[10], factor_and_data=1, mgii_only=1)
IDL> help, data_iris
** Structure <1ad7c08>, 22 tags, length=31190448, data length=31190348, refs=1:
FILE_NAME STRING '/Users/asainz/ciencia/IRIS-01/data/2017060'...
INDEX STRUCT -> <Anonymous> Array[1]
LAMBDA_FILE DOUBLE Array[626]
DATE STRING '05-Jun-17 19:17:43'
POS_XY DOUBLE Array[2]
LONLAT DOUBLE Array[2]
MU DOUBLE 0.94271849
SPCSCL DOUBLE 0.025459999
SPXSCL DOUBLE 0.16635000
SPYSCL DOUBLE 1.9954402
DATE_IN_FILENAME
STRING '20170605_170000'
IRIS_OBS_CODE STRING '3662506342'
RASTER_INFO STRING '00099'
EXTENT_ARCSEC_ARCSEC
DOUBLE Array[4]
EXTENT_PX_PX FLOAT Array[4]
EXTENT_PX_ARCSEC
DOUBLE Array[4]
EXTENT_TIME_ARCSEC
DOUBLE Array[4]
EXTENT_TIME_PX DOUBLE Array[4]
EXTENT_OPT DOUBLE 129.42030
EXTENT_OPT_COORDS
STRING Array[2]
LIST_EXTENT_COORDS
OBJREF <ObjHeapVar19914(LIST)>
FACTOR STRUCT -> <Anonymous> Array[1]
IDL> help, data_iris.factor, /str
** Structure <1ad6408>, 10 tags, length=31180024, data length=31180020, refs=2:
FILENAME STRING Array[1]
VER_RF2 STRING 'L12-2017-04-23'
BLD_VERS STRING 'V9R0X'
FLUX_UNITS STRING 'Hz (cgs): erg cm^-2 s^-1 sr^-1 Hz^-1'
INFO_SOURCE_FILE
STRING Array[1]
---> DATA_MG_II_K_2796
FLOAT Array[626, 778, 16]
DN2PHOTON_MG_II_K_2796
FLOAT 18.0000
FACTOR_MG_II_K_2796
DOUBLE Array[626]
TEXP_ARR_MG_II_K_2796
DOUBLE Array[16]
---> WL_MG_II_K_2796 DOUBLE Array[626]
4.3.2. Inspecting the Mg II h&k data cube¶
Let’s take look to a medium-size raster. You can select points in the data cube By pressing S while tvg function is running. Both the selecteped points location and the data corresponding to these points are returned by tvg, in this example in sel variable, as a structure or an array of structures if more than 1 point is selected:
IDL> data_iris = get_info_irisl2(files[2], factor_and_data=1, mgii_only=1)
IDL> iris_mgii = transpose(data_iris.factor.DATA_MG_II_K_2796, [2,1,0])
IDL> sel = tvg(iris_mgii)
IDL> help, sel
SEL STRUCT = -> <Anonymous> Array[3]
IDL> help, sel, /str
** Structure <75f8468>, 3 tags, length=2128, data length=2128, refs=1:
X_DISPLAY LONG 205
Y_DISPLAY LONG 200
DATA_XY FLOAT Array[530]
IDL> dim_sel = size(sel)
IDL> .r
IDL> for j=0, dim_sel[1]-1 do begin
IDL> plot, wl_mgii, sel[j].data_xy
IDL> wait, 5
IDL> endfor
IDL> end
Warning
The function tvg is in a very early stage. It will be updated and realeased soon.
Let’s crop the profiles in wavelentgh to avoid border effects, and transform the wavelentgh values to the air. Let’s see how they look with respect to the positions of the Mg II h&k lines
IDL> obs_mgii = sel.data_xy[50:450,*]
IDL> wl_obs = vac2air(wl_mgii[50:450])
IDL> mgII_k_pos_k3 = 2795.528
IDL> mgII_h_pos_h3 = 2802.704
IDL> mgII_uv_triplet_pos23 = 2797.930 + (2797.998 - 2797.930)
IDL> .r
IDL> for j=0, dim_sel[1]-1 do begin
IDL> plot, wl_obs, obs_mgii[*,j]
IDL> miver, [mgII_k_pos_k3, mgII_h_pos_h3, mgII_uv_triplet_pos23], line=2
IDL> wait, 5
IDL> endfor
IDL> end
4.3.3. Loading the \(IRIS^{2}\) database¶
Now, we are going to load the \(IRIS^{2}\) database:
IDL> wl_db = readfits(dir_sswdb+'/iris/iris2/latest/wl_iris2.v2.fits')
IDL> prof_mu_db = readfits(dir_sswdb+'/iris/iris2/latest/inv_mu_iris2.v2.0.fits')
IDL> model_db = readfits(dir_sswdb+'/iris/iris2/latest/mod_iris2.v2.0.fits')
IDL> help, wl_db, prof_mu_db, model_db
WL_DB DOUBLE = Array[473]
PROF_MU_DB DOUBLE = Array[474, 50599]
MODEL_DB DOUBLE = Array[4, 39, 50599]
IDL> prof_db = prof_mu_db[0:472,*]
IDL> mu_db = prof_mu_db[473,*]
Note
The dimension of wl_db correspond to spectral sampling in \(AA\).The 1st dimenstion of prof_mu_db corresponds to 473 spectral samples and the last value (474) to the \(\mu\) (cosine of the heliocentric angle) of the observation of that profile. The 2nd dimension corresponds to the number of samples in the \(IRIS^{2}\) database. The 1st dimension of mode_db corresponds to the physical variables \(T,\ v_{los},\ v_{turb},\ and\ electron\ density\ (n_{e})\). The 3rd corresponds to the number of samples in the \(IRIS^{2}\) database.
Let’s take a look at a profile and its corresponding model atmosphere:
IDL> !p.multi =[0,2,3,0]
IDL> j_db = 500
IDL> window, xs=700, ys=700
IDL> plot, wl_db, prof_db[*,j_db], xtit='Wavelength [AA]', charsize=2, $
IDL> pos = [.1,.7,.95,.95], /norm, /dev
IDL> miver, [mgII_k_pos_k3, mgII_h_pos_h3, mgII_uv_triplet_pos23]
IDL> plot, indgen(10), xstyle=4, ystyle=4, color=0 ; plots nothing
IDL> scales = [1e3, 1e5, 1e5, 1]
IDL> titles = ['T [kK]', 'vlos [km/s]', 'vturb [km/s]', 'ne [cm^-3]']
IDL> ltau = findgen(39)*.2 - 7.6
IDL> for i = 0, 3 do plot, ltau, model_db[i,*,j_db]/scales[i], title=titles[i], xtitle='log(tau)', charsize=2
IDL> !p.multi=0
4.3.4. Lookind for the closest RP in the \(IRIS^{2}\) database¶
Let’s interepolate the \(IRIS^{2}\) database to the wavelengths observed:
IDL> dim_prof_db= size(prof_db)
IDL> dim_obs_mgii=size(obs_mgii)
IDL> prof_db2obs = dblarr(dim_obs_mgii[1], dim_prof_db[2])
IDL> for i = 0, dim_prof_db[2]-1 do prof_db2obs[*, i] = interpol(prof_db[*,i], wl_db, wl_obs)
Now, let’s look for the closest RP in the \(IRIS^{2}\) database interpolated to the observation spectral samples:
IDL> euc_dist = calc_dist(obs_mgii, prof_db2obs, mio=1)
IDL> w = min(transpose(euc_dist), closest, dim=1)
IDL> closest = closest mod dim_prof_db[2]
IDL> print, closest
IDL> .r
IDL> for j = 0, dim_sel[1]-1 do begin
IDL> w = min(euc_dist[j,*], closest)
IDL> print, j, closest
IDL> endfor
IDL> end
We can accelarate the look-up process working in the Principal Component Analysis (PCA) space. Thus, we calculate: the PCA space defined by the
\(IRIS^{2}\) database interepolated to the observed profiles; the coefficients of the interepolated \(IRIS^{2}\) database in the that PCA space;
the coefficients of the observed profiles in the interpolated \(IRIS^{2}\) PCA space; the Euclidean distance between the PCA coefficients of the observed
profiles and the PCA coefficients of the interpolated profiles of the \(IRIS^{2}\) database; finally, we look for
the closest match. The number of PCA coefficients n_pca is determined by the user. Usually, a number larger than 20 is enough to represent accuaretly
all the profiles in a data set. A larger number of coefficients will reproduce the profile with a larger accueracy, including the noise of the profile. We should try to avoid that:
IDL> n_pca = 50
IDL> print, 'Using PCA coef: ', n_pca
IDL> mypca, transpose(prof_db2obs), eval, evec
IDL> project_pca, transpose(prof_db2obs), eval, evec, n_pca, prof_db2obs_pca_pca, prof_db2obs_pca_coef
IDL> prof_db2obs_pca_coef = transpose(prof_db2obs_pca_coef)
IDL> project_pca, transpose(obs_mgii), eval, evec, n_pca, obs_mgii_pca, obs_mgii_pca_coef
IDL> obs_mgii_pca_coef = transpose(obs_mgii_pca_coef)
IDL> euc_dist = calc_dist(obs_mgii_pca_coef, prof_db2obs_pca_coef, mio=1)
IDL> .r
IDL> for j =0, dim_sel[1]-1 do begin
IDL> w = min(closest[j,*], aux)
IDL> print, j, aux
IDL> endfor
IDL> end
IDL> ; Vectorized approach (faster)
IDL> w = min(transpose(euc_dist), closest_pca, dim=1)
IDL> closest_pca = closest_pca mod dim_prof_db[2]
IDL> print, closest_pca
IDL> print, closest
Several tests have demostrated that for n_pca\(> 50\), the difference between the closest RPs found by using the whole profile, given in the
example by closest, and the closest RPs by using the PCA coefficients, given in the example by closest_pca, is less than 0.1%.
This approach was introduced by Rees, D.E. et al. 2000.
4.3.5. Inversion of Mg II h&k lines based in RPs inverted by STiC¶
Let’s plot the observed data, the closest RP in the \(IRIS^{2}\) databas, and its corresponding model atmosphere:
IDL> .r
IDL> window, xs=700, ys=700
IDL> !p.multi =[0,2,3,0]
IDL> scales = [1e3, 1e5, 1e5, 1]
IDL> titles = ['T [kK]', 'vlos [km/s]', 'vturb [km/s]', 'ne [cm^-3]']
IDL> ltau = findgen(39)*.2 - 7.6
IDL> for j = 0, dim_sel[1]-1 do begin
IDL> plot, wl_obs, obs_mgii[*,j]*cont_level, xtit='Wavelength [AA]', $
IDL> charsize=2, pos = [.1,.7,.95,.95], /norm, /dev, line=2, xstyle=1
IDL> oplot, wl_obs, prof_db2obs[*,closest[j]]*cont_level
IDL> miver, [mgII_k_pos_k3, mgII_h_pos_h3, mgII_uv_triplet_pos23], line=2
IDL> plot, indgen(100), xstyle=4, ystyle=4, color=0 ; plots nothing
IDL> for i = 0, 3 do plot, ltau, model_db[i,*,closest[j]]/scales[i], $
IDL> title=titles[i], xtitle='log(tau)', charsize=2
IDL> wait, 15
IDL> endfor
IDL> !p.multi=0
IDL> end
Congratulations! You have done your first inversion of IRIS Mg II h&k lines based on the Representative Profiles inverted by STiC. Now, it is when the fun starts!