# Determine spectral slope of science source

(Difference between revisions)

These notes are from a project to determine the spectral slope of DG Tau, which is a protostellar object.

We assume your data has been properly flagged.

There are 16 continuum windows (500 MHz BW) in this example.

## Determine Spectral Index of Passband Calibrator

set flux=URANUS
set pass=3C84
set vis=uvdata
set bands = ( 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 )
#Brightness temperature of URANUS at our observing frequencies for flux calibration
#Derived using function provided by A. Isella. He derived it from WMAP data (and references therein).
```   #Tb=b1 + b2*log(x) + b3*log(x)^2 + b4*log(x)^3
#b1=154.813 b2=106.598 b3=-91.0474 b4=15.0796
```

set uranustb = ( 125.02683 125.20088 125.37588 125.55175 125.72856 125.90630 126.08495 126.26458 121.98721 121.82910 121.67174 121.51508 121.35952 121.20437 121.04997 120.89625 )

#C8 at center of this D-array

set refant=8

#Prepare data

rm -Rf \$vis uvcat vis=c0800.4D_102HH211.2.miriad out=\$vis select="-auto" options=nopass,nocal,nopol

#Flag visibility data

foreach band ( \$bands )

``` #Select individual window
/bin/rm -Rf \$vis.\$band
uvcat vis=\$vis out=\$vis.\$band select="'win(\$band)'" options=nocal,nopol,nopass
```
``` #Run mfcal on this window on very short timescale to minimize decorrelation
mfcal vis=\$vis.\$band select="source("\$pass")" interval=0.1,10000. refant=\$refant edge=2 line=chan,39,1,1,1
```
``` #Check mfcal solution
#gpplt vis=\$vis.\$band yaxis=amp device=1/xs nxy=4,4
#set text = "Press ENTER .."
#if (\$#argv > 0) set text = "\$* >>"
#printf "\n\$text"
#set junk = (\$<)
```
``` #Apply mfcal solution to correct for variation across the bandpass
/bin/rm -Rf \$vis.\$band.pb
uvcat vis=\$vis.\$band out=\$vis.\$band.pb options=nocal,nopol

#STARTING STEPS FOR BOOTFLUX
```
```   #BANDPASS SOURCE PHASE CALIBRATION
#Run selfcal on bandpass source do derive gain solutions for the bandpass source
#Interval can be short since just applying it to the same source (?)
mselfcal vis=\$vis.\$band.pb select="source("\$pass")" interval=1.01 options=apriori,noscale,phase refant=\$refant line=chan,1,1,39,39

#Check selfcal solution for bandpass source

#Apply selfcal phase solution from the bandpass calibrator to the bandpass calibrator
/bin/rm -Rf \$vis.\$band.\$pass.scp
uvcat vis=\$vis.\$band.pb out=\$vis.\$band.\$pass.scp select="source("\$pass")" options=nopass,nopol

#FLUX SOURCE PHASE CALIBRATION
#Run selfcal on flux source to derive gain solutions for the flux source
#Interval can be short since just applying it to the same source (?)
mselfcal vis=\$vis.\$band.pb select="source("\$flux")" interval=1.01 options=apriori,noscale,phase refant=\$refant line=chan,1,1,39,39

#Check selfcal solution for flux source
```
```   #Apply selfcal phase solution from the flux calibrator to the flux calibrator
/bin/rm -Rf \$vis.\$band.\$flux.scp
uvcat vis=\$vis.\$band.pb out=\$vis.\$band.\$flux.scp select="source("\$flux")" options=nopass,nopol
```
```   #UVCAT bandpass and phase corrected visiblity data for the BANDPASS and FLUX calibrator together
#This is the visibility file to use in bootflux
/bin/rm -Rf \$vis.\$band.btrdy
uvcat vis=\$vis.\$band.\$pass.scp,\$vis.\$band.\$flux.scp out=\$vis.\$band.btrdy options=nopass,nopol,nocal
```
```   #BOOTFLUX STEP
#print out Tb of URANUS at this window
echo \$uranustb[\$band]
#run bootflux and dump output into log file
bootflux vis=\$vis.\$band.btrdy select="source("\$pass,\$flux")" line=chan,1,1,39,39 taver=5.01 primary=\$flux,\$uranustb[\$band] badres=70 log=bootfluxJ.\$pass.log2
#grep that log file and extract the bootstrapped Average Flux and Flux Error of 3C84 at this window
grep Average bootfluxJ.\$pass.log2 | grep Flux | cut -c 17-31 > script1_bootflux.\$pass.\$band
#print the result to the screen
set tempfluxpass = `awk 'NR==1' script1_bootflux.\$pass.\$band`
echo \$tempfluxpass

#clean up
/bin/rm -Rf \$vis.\$band \$vis.\$band.pb \$vis.\$band.\$pass.scp \$vis.\$band.\$flux.scp \$vis.\$band.btrdy bootflux.\$pass.log2
```

end

## Determine Spectral Index of Phase Calibrator

set onsource=IRAS4

set flux=URANUS

set pass1=3C84

set pass2=3C84

set cal=3C84

set vis=uvdata

set bands = ( 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 )

set uranustb = ( 124.89169 125.05010 125.20929 125.36919 125.52988 125.69135 125.85355 126.01654 121.49650 121.35417 121.21234 121.07124 120.93065 120.79076 120.65140 120.51263 )

set refant=8

mfcal vis=\$vis select="source("\$pass1")" interval=0.1,10000. refant=\$refant edge=2 line=chan,624,1,1,1 flux=1.868,103.0,-0.695 /bin/rm -Rf \$vis.fullpb uvcat vis=\$vis out=\$vis.fullpb options=nopol,nocal

foreach band (\$bands)

``` #Select individual window
/bin/rm -Rf \$vis.\$band
uvcat vis=\$vis.fullpb out=\$vis.\$band select="'win(\$band)'" options=nocal,nopol,nopass
```
``` mselfcal vis=\$vis.\$band select="source("\$cal")" interval=1.01 options=noscale,phase,apriori refant=\$refant line=chan,1,1,39,39
#smagpplt vis=\$vis.fullpb yaxis=phase device=3/xs nxy=4,4
/bin/rm -Rf \$vis.\$cal.\$band.scp
uvcat vis=\$vis.\$band out=\$vis.\$cal.\$band.scp select="source("\$cal")" options=nopass,nopol
```
``` mselfcal vis=\$vis.\$band select="source("\$flux")" interval=1.01 options=noscale,phase,apriori refant=\$refant line=chan,1,1,39,39
#smagpplt vis=\$vis.fullpb yaxis=phase device=3/xs nxy=4,4
/bin/rm -Rf \$vis.\$flux.\$band.scp
uvcat vis=\$vis.\$band out=\$vis.\$flux.\$band.scp select="source("\$flux")" options=nopass,nopol
```
``` /bin/rm -Rf \$vis.\$band.btrdy
uvcat vis=\$vis.\$cal.\$band.scp,\$vis.\$flux.\$band.scp out=\$vis.\$band.btrdy options=nopass,nopol,nocal
```
``` bootflux vis=\$vis.\$band.btrdy select="source("\$cal,\$flux")" line=chan,1,1,39,39 taver=5.01 badres=70 primary=\$flux,\$uranustb[\$band] log=bootfluxJ.\$cal.log2
grep Average bootfluxJ.\$cal.log2 | grep Flux | cut -c 17-31 > script2_bootflux.\$cal.\$band
set tempfluxcal = `awk 'NR==1' script2_bootflux.\$cal.\$band`
echo \$tempfluxcal

/bin/rm -Rf \$vis.\$band \$vis.\$cal.\$band.scp \$vis.\$flux.\$band.scp \$vis.\$band.btrdy
/bin/rm -Rf bootfluxJ.\$cal.log2
```

end

/bin/rm -Rf \$vis.fullpb

## Determine Spectral Index of Science Source

set onsource=IRAS03 set flux=URANUS set cal=3C84 set vis=uvdata

set bands = ( 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 )

set uranustb = ( 125.027 125.201 125.376 125.552 125.729 125.906 126.085 126.265 121.987 121.829 121.672 121.515 121.360 121.204 121.050 120.896 )

set sources = ( IRAS03 B1c HH211 IRAS2A )

set refant=8

set fg=11.868

foreach onsource ( \$sources )

```   mfcal vis=\$vis select="source("\$cal")" interval=0.1,10000. refant=\$refant edge=2 line=chan,624,1,1,1 flux=11.868,103.0,-0.695
/bin/rm -Rf \$vis.fullpb
uvcat vis=\$vis out=\$vis.fullpb options=nopol,nocal
```
```   mselfcal vis=\$vis.fullpb select="source("\$cal")" interval=18.01 options=noscale,phase,apriori refant=\$refant line=chan,1,1,624,624
/bin/rm -Rf \$vis.scp
uvcat vis=\$vis.fullpb select="source("\$cal,\$onsource")" out=\$vis.scp options=nopass,nopol
```
```   mselfcal vis=\$vis.scp select="source("\$cal")" interval=40.01 options=noscale,amplitude refant=\$refant line=chan,1,1,624,624 flux=\$fg
/bin/rm -Rf \$vis.\$onsource.scpa
uvcat vis=\$vis.scp select="source("\$cal,\$onsource")" out=\$vis.\$onsource.scpa options=nopass,nopol
```
```   #Remove DGTAUB with uvmodel
#/bin/rm -Rf \$vis.scpam \$vis.scpam2 \$vis.calibrators \$vis.scpamf
#uvmodel vis=\$vis.scpa out=\$vis.scpam select="source("\$onsource")" options=subtract flux=0.015 offset=-28.592,-45.934
#uvmodel vis=\$vis.scpam out=\$vis.scpam2 select="source("\$onsource")" options=subtract flux=0.0010 offset=-28.887,-43.999
#uvcat vis=\$vis.scpa out=\$vis.calibrators select="source("\$cal,\$pass1,\$pass2")" options=nopass,nocal,nopol
#uvcat vis=\$vis.calibrators,\$vis.scpam2 out=\$vis.scpamf options=nopass,nocal,nopol

#Make image of IRAS2A to see peak flux is at expected phase center
/bin/rm -Rf \$vis.fullmap \$vis.fullbm \$vis.fullcln \$vis.\$onsource.scpamres
invert vis=\$vis.\$onsource.scpa select="source("\$onsource")" map=\$vis.fullmap beam=\$vis.fullbm imsize=257,257 cell=1,1 options=systemp robust=0 line=chan,1,1,624,624
clean map=\$vis.fullmap beam=\$vis.fullbm out=\$vis.fullcln gain=0.1 niters=2000
restor model=\$vis.fullcln beam=\$vis.fullbm map=\$vis.fullmap out=\$vis.\$onsource.scpamres
```
```   #Is phase center of DGTau where expected?
#imfit in=\$vis.scpamres object=gaussian region="rel,box(-20,-20,21,21)" > script3_imfit.tbl

#extract ra and dec position of center position from imfit
#grep Declination script3_imfit.tbl | cut -c 35-36 > script3_deccenter1
#grep Declination script3_imfit.tbl | cut -c 38-39 > script3_deccenter2
#grep Declination script3_imfit.tbl | cut -c 41-46 > script3_deccenter3
#set deccenterval1 = `awk 'NR==1' script3_deccenter1`
#set deccenterval2 = `awk 'NR==1' script3_deccenter2`
#set deccenterval3 = `awk 'NR==1' script3_deccenter3`
#grep Ascension script3_imfit.tbl | cut -c 30-35 > script3_racenter1
#grep Ascension script3_imfit.tbl | cut -c 37-38 > script3_racenter2
#grep Ascension script3_imfit.tbl | cut -c 40-45 > script3_racenter3
#set racenterval1 = `awk 'NR==1' script3_racenter1`
#set racenterval2 = `awk 'NR==1' script3_racenter2`
#set racenterval3 = `awk 'NR==1' script3_racenter3`
#Correct phase center
#/bin/rm -Rf temp.scpamf
#uvedit vis=\$vis.scpamf out=temp.scpamf source=\$onsource ra=\$racenterval1,\$racenterval2,\$racenterval3 dec=\$deccenterval1,\$deccenterval2,\$deccenterval3 options=nouv
#/bin/rm -Rf \$vis.scpamf
#uvcat vis=temp.scpamf out=\$vis.scpamf options=nopass,nocal,nopol

foreach band (\$bands)
```

uvamp vis=\$vis.\$onsource.scpa select="source("\$onsource"),win("\$band")" bin=10,5,klam log=script3_uvamp.\$onsource.\$band uvamp vis=\$vis.\$onsource.scpa select="source("\$cal"),win("\$band")" bin=1,50,klam log=script3_uvamp.\$cal.\$band

```   end

/bin/rm -Rf \$vis.fullpb \$vis.scp
/bin/rm -Rf \$vis.scpam \$vis.scpam2 \$vis.calibrators
/bin/rm -Rf \$vis.fullmap \$vis.fullbm \$vis.fullcln \$vis.scpamres
/bin/rm -Rf temp.scpamf
```

end