Vertical-Current Approximation NonLinear Force-Free Field (VCA-NLFFF) Code Tutorial [Version 2020 Feb 12]

This tutorial is intended for SolarSoftWare (SSW) users who want to run the "Vertical-Current Approximation Nonlinear Force-Free Field" (VCA-NLFFF) code using SDO/HMI and SDO/AIA (or other) data to obtain a nonlinear force-free field solution of the magnetic field in an active region or a flare. Descriptions and performance tests of the VCA-NLFFF code are given in the following publication (and references therein):

Aschwanden M.J. (2016), The Astrophys. J. Suppl. Series Vol. 224, article 25, (32 pp). The Vertical Current Approximation Nonlinear Force-Free Field Code - Description, Performance Tests, and Measurements of Magnetic Energies Dissipated in Solar Flares
../eprints/2016_vca_nlfff.pdf"

Aschwanden,M.J. 2013, SP 287, 323-344, A nonlinear force-free magnetic field approximation suitable for fast forward-fitting to coronal loops. I. Theory,
../eprints/2013_fff1.pdf"

Aschwanden, M.J. and Malanushenko,A. 2013, SP 287, 345-367, A nonlinear force-free magnetic field approximation suitable for fast forward-fitting to coronal loops. II. Numeric Code and Tests,
../eprints/2013_fff2.pdf"

Aschwanden, M.J. 2013, SP 287, 369-389, A nonlinear force-free magnetic field approximation suitable for fast forward-fitting to coronal loops. III. The Free Energy
../eprints/2013_fff3.pdf"

Aschwanden, M.J., DePontieu, B., and Katrukha, E. 2013, Entropy, 15(8), 3007-3030, Optimization of Curvi-Linear Tracing Applied to Solar Physics and Biophysics
www.mdpi.com/1099-4300/15/8/3007"

Aschwanden, M.J. 2016, ApJSS 224, 25 (32pp), The Vertical Current Approximation Nonlinear Force-Free Field Code - Description, Performance Tests, and Measurements of Magnetic Energies Dissipated in Solar Flares
../eprints/2016_vca_nlfff.pdf"

Warren, H.P., Crump,N.A., Ugarte-Urra,I., Sun,X., Aschwanden,M.J., and Wiegelmann,T. 2018, ApJ 860, 46, (13pp), Towards a quantitative comparison of magnetic field extrapolations and observed coronal loops
../eprints/2018_warren.pdf"

Aschwanden, M.J. and Wang,T.J. 2020, ApJ (in press), Torsional slow-mode oscillations discovered in the magnetic free energy during solar flares
../eprints/2020_tors.pdf"

IDL/SSW Implementation

The entire VCA-NLFFF code is written in IDL and requires outside the IDL software a file SETUP.SSW_ENV that contains the different SSW instrument packages from STEREO/SECCHI, AIA/SDO, and HMI/SDO: SETENV SSW_INSTR "GEN STEREO SECCHI AIA HMI" . The IDL routines of the VCA-NLFFF code are contained in the package $/SSW/PACKAGES/MJASTEREO/IDL/ , which are automatically included in the IDL path with the command SETENV SSW_INSTR "STEREO" . The source codes of the main IDL procedures that run the VCA-NLFFF code have the names NLFFF_*.PRO The acronym MJA stands for the initials of the author. The commands given here in uppercase letters should work also with lowercase lettering.

X-Class flare on 2014 March 29


The following IDL script is an example of a VCA-NLFFF run for a single time step. The code is modular and each procedur "nlfff_*.pro" can be run separately, but they have to be run in the correct time order as given in the following example. Once all tasks are completed, the results are stored in the "savefile" and can be displayed later (using the short cut GOTO,DISP) without re-running the tasks.

The following IDL script can be downloaded either from the SSW library or from here:
nlfff_tutorial.pro"
Before you run this script inside IDL, create first a work directory dir='~/work/' and a data directory datadir='~/SDODATA/' An event catalog will automatically be created for all time steps of an event. Choosing a name, RUN='runZ', allows to distinguish runs of the same events with different control parameter settings.


RUN SPECIFICATION_________________________________
iflare = 1 ;1,...,nev = flare ID number (from any catalog)
it = 1 ;1,...,nt, time step (nt=duration/cadence)
cube = 0 ;0=no B-cube 1=B-cube
test = 1 ;0=none, 1=test, 2=detailed
io = 0 ;0=screen, 1=ps, 2=eps, 3=col.ps
run = 'runZ' ;name of run
dir = '~/work/' ;work directory name
datadir = '~/SDODATA/' ;data directory name
catfile = 'event_catalog_'+run+'.txt' ;catalog file name


INPUT DEFINITION _________________________________
tstart = '2014-03-29T18:17:00' ;start time (UT)
cadence = 360. ;time cadence (s)
duration = 4500. ;total time duration (s)
noaa = '12017' ;active region number
helpos = 'N07W30' ;heliographic position
instr = 'AIA' ;instrument (AIA, IRIS)
wave8 =[94,131,171,193,211,335,0,0] ;array with 8 wavelengths (Angstroem)
fov0 = 0.3 ;field-of-view (solar radii)
hmin = 0.003 ;minimum altitude (solar radii)
hmax = 0.200 ;minimum altitude (solar radii)
nsm1 = 1 ;hiphass filter (pixel)
nmag_p = 200 ;number of potential sources
nmag_np = 10 ;number of non-potential sources
amis = 20.0 ;misalignment angle limit (deg)
nitmin = 10 ;min number of iterations
nitmax = 15 ;max number of iterations
nstruc = 1000 ;max number of loops per wavelength
prox_min = 10.0 ;proximity of loop to next magn source
lmin_wid = 5.0 ;minimum loop length/width ratio
rmin_wid = 5.0 ;minimum curvature radius/width ratio
qthresh1 = 0.0 ;threshold flux (in medians)
qthresh2 = 0.0 ;threshold filter flux (in medians)
ngap = 0 ;gap in tracing pixels
nh = 5 ;number of altitude profiles per loop
nseg = 5 ;number of segments per loop
ncurv = 5 ;number of curvature models in altitude profile

input={catfile:catfile,tstart:tstart,cadence:cadence,duration:duration,$
noaa:noaa,helpos:helpos,instr:instr,wave8:wave8,fov0:fov0,$
hmin:hmin,hmax:hmax,nsm1:nsm1,nmag_p:nmag_p,nmag_np:nmag_np,$
amis:amis,nitmin:nitmin,nstruc:nstruc,prox_min:prox_min,$
lmin_wid:lmin_wid,rmin_wid:rmin_wid,qthresh1:qthresh1,$
qthresh2:qthresh2,ngap:ngap,nh:nh,nseg:nseg,ncurv:ncurv}

A template IDL script can be downloaded either from the SSW library or from here:
nlfff_tutorial.pro"
All following examples can be run by modifying this IDL script template, such as by commenting or un-commenting of the desired script lines.


COMPUTATION TASKS _________________________________

nlfff_cat, input,catfile,dir
nlfff_init, input,catfile,it,run,dir,savefile
nlfff_hmi, dir,savefile,test,datadir
nlfff_aia, dir,savefile,test,wave_,datadir
nlfff_fov, dir,savefile,test,datadir
nlfff_magn, dir,savefile,test,datadir
nlfff_tracing, dir,savefile,test,datadir
nlfff_fit, dir,savefile,test
nlfff_field, dir,savefile,test
nlfff_cube, dir,savefile,test,cube
nlfff_merit, dir,savefile,test
nlfff_energy, dir,savefile,test
nlfff_disp, dir,catfile,iflare,it,run,io
nlfff_disp2, dir,catfile,iflare,it,run,overlay,ngrid,bmin


The task NLFFF_CAT generates a catalog "event_tutorial_runZ.txt" with NT=duration/cadence = 4500 s / 360 s = 13 time steps, defining the UT times and heliographic positions for each time step. The task NLFFF_INIT reads the time, heliographic position, and wavelengths of the requested time step IT=1 from the catalog and generates an IDL savefile with the name "20140329_181700_AIA_001_runZ.sav" for this run, in which all results will be saved. The task NLFFF_HMI reads the corresponding HMI/SDO image from the Stanford/JSOC archive, if the magnetogram file does not aready exist in the specified local directory. The HMI Fitsfile is named "20140329_181700_HMI.fits". The procedure NLFFF_AIA downloads the 6 requested AIA full-disk images from the Lockheed Martin or Stanford/JSOC archive, if the datafiles do not already exist in the specified local directory. The procedure NLFFF_IRIS reads the corresponding IRIS slit-jaw images if necessary. The data input is now complete and the work directory should contain the following data files:

If this code is not able to download the datafiles directly, please download them from JSOC and make sure that the files are renamed according to the naming convention used here. For this particular example you can download the data files from the follwing webpage:

20140329_181700_AIA_131.fits"
20140329_181700_AIA_171.fits"
20140329_181700_AIA_193.fits"
20140329_181700_AIA_211.fits"
20140329_181700_AIA_335.fits"
20140329_181700_AIA_94.fits"
20140329_181700_HMI.fits"

The task NLFFF_FOV extracts subimages in the specified field-of-view with size FOV=0.3 solar radii and centered around the heliographic position N07W30. The subimages are saved as separate FITS files and the control display looks like this (if the parameter TEST=1 is set):

The task NLFFF_MAGN extracts the corresponding subimage from the HMI magnetogram and decovolves the line-of-sight magnetogram into the specified NMAG_P=200 potential field magnetic sources. A control display (if TEST=1 is set) is shown with the original magnetogram (top left), the decomposed model of the magnetogram (bottom left), the difference between the observed and the model diagram (top right panel), as well as an East-West scan through both the original (red curve) and model magnetogram (white curve) (bottom right panel).

The task NLFFF_TRACING traces then in each of the 8 AIA wavelength images curvi-linear features and displays the tracings superimposed on one of the AIA images (if TEST=1 is set). The task NLFFF_FIT is then the main computation task that fits a nonlinear force-free solution to the geometry of the traced loops, by minimizing the median 3D misalignment angle. Typical computation times are about 1 minute. The task NLFFF_FIELD computes then a specific data set of modeled field lines that intersect the midpoints of the automatically traced loop segments. The procedure NLFFF_CUBE can be used to calculate an entire 3D cube of B-field vectors, if CUBE=1 is set, otherwise a minimum cube with only nz=3 bottom layers is calculated. The routine NLFFF_MERIT calculates also some figures of merit, and the task NLFFF_ENERGY integrates the non-potential and free energy over the entire computation volume. In a final step, a quick display of the results can be produced with the procedure NLFFF_DISP, which shows the obtained 3D loop trajectories and field lines from two different aspect angles (line-of-sight and 90 deg rotated to North. A postscript file is automatically produced from this display, with the filename "20140329_181700_AIA_001_runA_col.ps":

This standard plot is created with the procedure NLFFF_DISP.PRO. Other overlay plots can be created with the procedure NLFFF_DISP2.PRO, using different overlay options or field line selections. The parameter OVERLAY=0 shows the HMI magnetogram as background. The choice of OVERLAY=171 will show the 171 A EUV image as background on a logarithmic color scale, while OVERLAY=-171 (negative wavelength) which show the highpass filtered 171 A image as background. For the choice of field line selections, either the field lines that intersect the midpoints of the traced loop segments can be chosen (standard option in NLFFF_DISP and NLFFF_DISP2), or a rectangular grid of footpoints, say NGRID=2 will select footpoints in a 20 x 20 grid, with a threshold set by the minimum magnetic field strength in units of Gauss, for instance BMIN=10. The output will be stored in a postscript file with the name "20140329_181700_AIA_001_runZ2.sav'. Note that NLFFF_DISP produces a plot names *runZ_col.ps, while the NLFFF_DISP2 produces a plot with *runZ2_col.ps.


DISPLAY OVERLAYS AND FIELD LINES_______________________________________
overlay = 0 ;HMI magnetogram
ngrid = 0 ;field lines intersecting midpoints of traced loops
bmin = 10 ;(G)
nlfff_cat, input,catefile,dir
nlfff_disp, dir,catfile,iflare,it,run
nlfff_disp2, dir,catfile,iflare,it,run,overlay,ngrid,bmin,datadir


DISPLAY EUV IMAGE AND FIELD LINES WITH FOOTPOINT GRID____________________
overlay = 171 ;171 A AIA/SDO image on log scale
ngrid = 20 ;field lines grid with 20 x 20 footpoints
bmin = 10 ;(G)
nlfff_disp2, dir,catfile,iflare,it,run,overlay,ngrid,bmin


DISPLAY HIGHPASS-FILTERED EUV IMAGE AND FIELD LINES_______________________________________
overlay = -171 ;171 A AIA/SDO image highpass-filtered
ngrid = 20 ;field lines grid with 20 x 20 footpoints
bmin = 10 ;(G)
nlfff_disp2, dir,catfile,iflare,it,run,overlay,ngrid,bmin



(A) Active Region NOAA 11110 observed on 2010 Sept 29


In this example we use a larger field of view (FOV=0.3 solar radii), only one wavelength INPUT.WAVE8=[0,0,171,0,0,0,0,0], a larger highpass filter NSM1=3, more magnetic sources, NMAG_P=60, NMAPG_NP=60, the largest range of acceptable misalignment angles INPUT.AMIS=90 (degrees), and a threshold of QTHRESH2=1.0 standard deviations in the filtered image.

(B) Active Region NOAA 11337 observed on 2011 Nov 8

(C) Active Region NOAA 11459 observed on 2012 Apr 21

(D) Active Region NOAA 11158 observed on 2011 Feb 15

(E) Active Region NOAA 11158 observed on 2011 Feb 13







E-mail: aschwanden@lmsal.com - Markus J.Aschwanden (Lockheed Martin Solar & Astrophysics Lab.)