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:

  1. Installation of IRIS SolarSoft branch, including IRIS database:

    0.1 Set enviroment variables to start SolarSoft including IRIS SolarSoft branch and IRIS database.

  2. Getting the IRIS Mg II h&k data properly calibrated (or any other IRIS data).

  3. Inspection and selection of points of interest to be inverted:

    2.1 Cropping the profiles in wavelength to avoid border effects.

  4. 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.

  5. 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.

  6. Displaying the inversion:

    5.1 Plot observed profile (dashed line) and inverted RP (thick) 5.2 Plot the corresponding Representative Model Atmosphere (RMA)

  7. 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_idl

In 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, $
           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. Because the \(IRIS^{2}\) database is normalized to the continuum, we will normalize the observed data similarly:

IDL> cont_level = 3.454e-6
IDL> obs_mgii = obs_mgii/cont_level
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]', $
           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], $
         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!