Difference between revisions of "Flare Imaging"

From EOVSA Wiki
Jump to: navigation, search
(Self-calibration)
 
(7 intermediate revisions by the same user not shown)
Line 1: Line 1:
 
== Preparation ==
 
== Preparation ==
First, identify the time when the flare in question happened. There are many ways to do this, but a nice way is to use the [http://sprg.ssl.berkeley.edu/~tohban/browser/ RHESSI Browser]  and check "EOVSA Radio Data" on the upper left. Use the time you identified to find the corresponding IDB file(s) under /data1/eovsa/fits/IDB/yyyymmdd/. Then, on pipeline, it is advised to go to your working directory and <b>copy the IDB data into your directory</b>. Note, <b>never, NEVER, work directly on the IDB data in the original data directory!</b> Here I use /data1/eovsa/fits/IDB/20170711/IDB20170711201620 (2017 July 11 C flare) as an example.  
+
First, identify the time when the flare in question happened. There are many ways to do this, but a nice way is to use the [http://sprg.ssl.berkeley.edu/~tohban/browser/ RHESSI Browser]  and check "EOVSA Radio Data" on the upper left. Use the time you identified to find the corresponding IDB file(s) under /data1/eovsa/fits/IDB/yyyymmdd/. Then, on pipeline, it is advised to go to your working directory and <b>copy the IDB data into your directory</b>. Note, <b>never, NEVER, work directly on the IDB data in the original data directory!</b> Here I use /data1/eovsa/fits/IDB/20170821/IDB20170821202020 (2017 Aug 21 C flare) as an example.
  
 
== Calibration ==
 
== Calibration ==
Line 7: Line 7:
 
casa
 
casa
 
</pre>
 
</pre>
1. Perform polarization rotation and attenuation calibration.
 
<pre>
 
idbfile='IDB20170711201620'
 
from eovsapy import pipeline_cal as pc
 
pc.udb_corr(idbfile,calibrate=True)
 
</pre>
 
This will create a new IDB dataset named IDB20170711201620_1
 
  
2. Import Miriad file into CASA.
+
* Import Miriad file into CASA.
 
<pre>
 
<pre>
importeovsa(idbfiles='IDB20170711201620_1', doscaling=False)
+
importeovsa(idbfiles='IDB20170821202020')
 
</pre>
 
</pre>
Remember it is important to set doscaling=False to disable scaling cross-correlation visibilities to auto-correlation. The output is automatically named "IDB20170711201620_1.ms".
+
idbfiles could be a string of the IDB file name, or a list of IDB files. The output is automatically named "IDB20170821202020.ms".
  
3. Perform reference phase calibration and multiband delay calibration (daily phase calibration). First, go to [http://ovsa.njit.edu/EOVSA/cal_status.php this page] and check if reference calibrations and daily phase calibrations have already been in the database. If not, create them following the instructions on [http://www.ovsa.njit.edu/wiki/index.php/Reference_Gain_Calibration reference calibration], [http://www.ovsa.njit.edu/wiki/index.php/Daily_Gain_Calibration daily phase calibration], and [http://www.ovsa.njit.edu/wiki/index.php/Total_Power_Calibration total power calibration].
+
* Perform reference phase calibration and multiband delay calibration (daily phase calibration). First, go to [http://ovsa.njit.edu/EOVSA/cal_status.php this page] and check if reference calibrations and daily phase calibrations have already been in the database. If not, create them following the instructions on [http://www.ovsa.njit.edu/wiki/index.php/Reference_Gain_Calibration reference calibration], [http://www.ovsa.njit.edu/wiki/index.php/Daily_Gain_Calibration daily phase calibration], and [http://www.ovsa.njit.edu/wiki/index.php/Total_Power_Calibration total power calibration].
  
 
Note for now there is an annoying issue in CASA 5.1.0 that clean does not run after any parallelized tasks. So it would be safe to <b>exit CASA before running calibeovsa if you set doimage=True</b>.
 
Note for now there is an annoying issue in CASA 5.1.0 that clean does not run after any parallelized tasks. So it would be safe to <b>exit CASA before running calibeovsa if you set doimage=True</b>.
 
<pre>
 
<pre>
calibeovsa(vis='IDB20170711201620_1.ms', caltype=['refpha','phacal'], doimage=True)
+
calibeovsa(vis='IDB20170821202020.ms', caltype=['refpha','phacal'], doimage=True)
 
</pre>
 
</pre>
This will calibrate the input visibility, write out calibration tables, and apply the calibration. If doimage=True, a quicklook image will be produced (by integrating over the entire time).
+
This will calibrate the input visibility, write out calibration tables under /data1/eovsa/caltable/, and apply the calibration. If doimage=True, a quicklook image will be produced (by integrating over the entire time).
  
Note both importeovsa and calibeovsa are now a CASA task. If you want to use more functionality of them and/or see more information of each parameter, you can follow the normal procedure of [https://casa.nrao.edu/casadocs/casa-5.0.0/usingcasa/casa-tasks using a CASA task].
+
Note both importeovsa and calibeovsa are (customized) CASA tasks. If you want to use more functionality of them and/or see more information of each parameter, you can follow the normal procedure of [https://casa.nrao.edu/casadocs/casa-5.0.0/usingcasa/casa-tasks using a CASA task].
 +
 
 +
==Self-calibration==
 +
As EOVSA does not calibrate as often as other general-purpose radio interferometers (e.g., VLA and ALMA), self-calibration is often needed in bringing out details of the flare at multiple frequencies. An example script, run in SunCASA, for doing self-calibration of the 2017 Aug 21 flare at ~20:20 UT can be found [https://github.com/binchensun/casa-eovsa/blob/master/slfcal_example.py at this Github link]. Here are some explanations on steps taken in this script.
 +
 
 +
===Step 1: Preparation===
 +
EOVSA observes the full solar disk, so we can in principle image the entire solar disk and perform self-calibration. However, as flares usually happen in an active region within a limited field of view and outshines everything else on the disk, self-calibrating the entire field of view is normally not needed (unless there are multiple bright sources on the Sun simultaneously, particularly if the flare of interest is relatively weak). The first step of selfcal is to make a clean image for a short duration of the flare, when the source structure is relatively simple and the signal-to-noise ratio is high across as many frequencies as possible. We use the full Sun view to find out the better phase center and appropriate field of view for later steps.
 +
 
 +
<pre>
 +
# ============ Prior definitions for spectral windows, antennas, pixel numbers =========
 +
spws=[str(s+1) for s in range(30)] # Use all 30 spectral windows available for imaging (the case for pre-2019 data)
 +
antennas='0~12'
 +
npix=512
 +
# parameters specific to the event (found from step 1)
 +
phasecenter='J2000 10h02m59 11d58m14'
 +
xran=[280,480]
 +
yran=[-50,150]
 +
 
 +
# =========== Step 1, doing a full-Sun image to find out phasecenter and appropriate field of view =========
 +
if dofullsun:
 +
    #initial mfs clean to find out the image phase center
 +
    im_init='fullsun_init'
 +
    os.system('rm -rf '+im_init+'*')
 +
    clean(vis=slfcalms,
 +
            antenna='0~12',
 +
            imagename=im_init,
 +
            spw='1~15',
 +
            mode='mfs',
 +
            timerange=trange,
 +
            imagermode='csclean',
 +
            psfmode='clark',
 +
            imsize=[npix],
 +
            cell=['5arcsec'],
 +
            niter=1000,
 +
            gain=0.05,
 +
            stokes='I',
 +
            restoringbeam=['30arcsec'],
 +
            interactive=False,
 +
            pbcor=True,
 +
            usescratch=True)
 +
 
 +
    hf.imreg(vis=slfcalms,imagefile=im_init+'.image',fitsfile=im_init+'.fits',
 +
            timerange=trange,usephacenter=False,verbose=True)
 +
    clnjunks = ['.flux', '.mask', '.model', '.psf', '.residual']
 +
    for clnjunk in clnjunks:
 +
        if os.path.exists(im_init + clnjunk):
 +
            shutil.rmtree(im_init + clnjunk)
 +
 
 +
    from sunpy import map as smap
 +
    from matplotlib import pyplot as plt
 +
    eomap=smap.Map(im_init+'.fits')
 +
    eomap.data=eomap.data.reshape((npix,npix))
 +
    eomap.plot_settings['cmap'] = plt.get_cmap('jet')
 +
    eomap.plot()
 +
    eomap.draw_limb()
 +
    eomap.draw_grid()
 +
    plt.show()
 +
    viewer(im_init+'.image')
 +
</pre>
  
==selfcalibration==
 
This section needs to be filled in.
 
 
 
 
=[[DG Notes on 20170904 Event]]=
 
=[[DG Notes on 20170904 Event]]=

Latest revision as of 14:01, 15 May 2019

Preparation

First, identify the time when the flare in question happened. There are many ways to do this, but a nice way is to use the RHESSI Browser and check "EOVSA Radio Data" on the upper left. Use the time you identified to find the corresponding IDB file(s) under /data1/eovsa/fits/IDB/yyyymmdd/. Then, on pipeline, it is advised to go to your working directory and copy the IDB data into your directory. Note, never, NEVER, work directly on the IDB data in the original data directory! Here I use /data1/eovsa/fits/IDB/20170821/IDB20170821202020 (2017 Aug 21 C flare) as an example.

Calibration

Start CASA in your working directory (it takes a minute or so to load)

casa
  • Import Miriad file into CASA.
importeovsa(idbfiles='IDB20170821202020')

idbfiles could be a string of the IDB file name, or a list of IDB files. The output is automatically named "IDB20170821202020.ms".

Note for now there is an annoying issue in CASA 5.1.0 that clean does not run after any parallelized tasks. So it would be safe to exit CASA before running calibeovsa if you set doimage=True.

calibeovsa(vis='IDB20170821202020.ms', caltype=['refpha','phacal'], doimage=True)

This will calibrate the input visibility, write out calibration tables under /data1/eovsa/caltable/, and apply the calibration. If doimage=True, a quicklook image will be produced (by integrating over the entire time).

Note both importeovsa and calibeovsa are (customized) CASA tasks. If you want to use more functionality of them and/or see more information of each parameter, you can follow the normal procedure of using a CASA task.

Self-calibration

As EOVSA does not calibrate as often as other general-purpose radio interferometers (e.g., VLA and ALMA), self-calibration is often needed in bringing out details of the flare at multiple frequencies. An example script, run in SunCASA, for doing self-calibration of the 2017 Aug 21 flare at ~20:20 UT can be found at this Github link. Here are some explanations on steps taken in this script.

Step 1: Preparation

EOVSA observes the full solar disk, so we can in principle image the entire solar disk and perform self-calibration. However, as flares usually happen in an active region within a limited field of view and outshines everything else on the disk, self-calibrating the entire field of view is normally not needed (unless there are multiple bright sources on the Sun simultaneously, particularly if the flare of interest is relatively weak). The first step of selfcal is to make a clean image for a short duration of the flare, when the source structure is relatively simple and the signal-to-noise ratio is high across as many frequencies as possible. We use the full Sun view to find out the better phase center and appropriate field of view for later steps.

# ============ Prior definitions for spectral windows, antennas, pixel numbers =========
spws=[str(s+1) for s in range(30)] # Use all 30 spectral windows available for imaging (the case for pre-2019 data)
antennas='0~12' 
npix=512
# parameters specific to the event (found from step 1)
phasecenter='J2000 10h02m59 11d58m14'
xran=[280,480]
yran=[-50,150]

# =========== Step 1, doing a full-Sun image to find out phasecenter and appropriate field of view =========
if dofullsun:
    #initial mfs clean to find out the image phase center
    im_init='fullsun_init'
    os.system('rm -rf '+im_init+'*')
    clean(vis=slfcalms,
            antenna='0~12',
            imagename=im_init,
            spw='1~15',
            mode='mfs',
            timerange=trange,
            imagermode='csclean',
            psfmode='clark',
            imsize=[npix],
            cell=['5arcsec'],
            niter=1000,
            gain=0.05,
            stokes='I',
            restoringbeam=['30arcsec'],
            interactive=False,
            pbcor=True,
            usescratch=True)

    hf.imreg(vis=slfcalms,imagefile=im_init+'.image',fitsfile=im_init+'.fits',
             timerange=trange,usephacenter=False,verbose=True)
    clnjunks = ['.flux', '.mask', '.model', '.psf', '.residual']
    for clnjunk in clnjunks:
        if os.path.exists(im_init + clnjunk):
            shutil.rmtree(im_init + clnjunk)

    from sunpy import map as smap
    from matplotlib import pyplot as plt 
    eomap=smap.Map(im_init+'.fits')
    eomap.data=eomap.data.reshape((npix,npix))
    eomap.plot_settings['cmap'] = plt.get_cmap('jet')
    eomap.plot()
    eomap.draw_limb()
    eomap.draw_grid()
    plt.show()
    viewer(im_init+'.image')

DG Notes on 20170904 Event