# Glish script to flag and calibrate a 334 MHz observation of W3
# After calibration, the data are gridded into a spectral line cube.
# HISTORY
# 03Jul10 GIL Filler changed to produce 1024 SP channels
# 03Apr01 GIL Initial version
#
# This script assumes that the AIPS++ dish package has been loaded.
# In GB this is done by typing:
# % gbdish
#             at the linux command prompt

# Import 17 scans into a measurement set for the AIPS++ dish package
gbtDataName := '/home/gbtmirror/gbtdata/T_PF1_03MAR13'    # set data import parameters
msName := 'T_PF1_03MAR13_SP'                    # output measurement set name
firstScan   := 60                               # start of scans to load
lastScan    := 75                               # end of scans to load
d.import( gbtDataName,,,firstScan, lastScan)
# Examine the overall properties of the Measurement set
d.gms()
# define the numbers scans to examine and image
scans := [firstScan:lastScan]

# flag, then examine the first cal On/Off state each scan of the raw data
for (i in scans) { foo := d.flag( i,,,[1:56],T); d.plotr(i,1);}

# Include a special script to flag all the RFI channels for this obs.
# The script below was generated by plotting
include '/users/glangsto/flareGlish/w3Flag.g'
#Type:
 flagChannels( rfiCenters, rfiWidths)
#and:
 flagChannels( bandCenters, bandWidths)
#to flag W3 with center frequency of 334 MHz channels

# An example of flagging is to remove the Band pass edges and noise spikes
# every 256 Channels, which occur in every spectral processor observation
# foo is a dummy variable to catch the excessive output
# In this observation, channels run 1 to 1024

for (i in scans) { foo := d.flag( i,,,[255:257],T);}
for (i in scans) { foo := d.flag( i,,,[511:513],T);}
for (i in scans) { foo := d.flag( i,,,[767:768],T);}
for (i in scans) { foo := d.flag( i,,,[976:1024],T); }

# Include the AIPS++ viewer to flag the bad time ranges.
include 'viewer.g'
# start up the viewer and flag the bad RFI times.
# This step is done with "clicks" and is difficult to document
# in the Gui click on "T_PF1_03MAR13_SP" then "Raster Image"
dv.gui()
# you should see the horizontal strips where the GBT scanned W3
# and verticle purple strips where the RFI was flagged.

# also note a few oddly scintillating horizontal strips where broad
# band RFI was present.  The next step is to flag these data.
# Flagging is easy, but fairly tricky to do the first time.

# In the image GUI, click on "Adjust" and another GUI appears
# In this new GUI, click on "Flagging Options"
# Click on "Channels", "Baselines", "Spectral Windows" and "Polarizations"
# leave un-selected "Times"
# Now the tricky part:  Back in the GUI with the image, make boxes.
# around the RFI  By default this is done with the Right Most Mouse
# button.  DURING FLAGGING DO NOT USE any other button in the image window!
# After making the box, the double click with the same mouse button.
# The RFI should turn light purple.
#
# Click around the other RFI time range.  It will turn light purple.
# NOW REMEMBER!  Go back to the Flagging Options GUI and click the
# save edits to Disk!  You are now done flagging.

# now cal_On - cal_Off calibrate all scans
#d.calib( firstScan)
# create the reference scan from the average of the first and last scans

d.sclear();                          # clear average
d.plotc( firstScan);                 # load average data 
d.accum();                           # add to average

d.plotc( lastScan);                  # load average data 
d.accum();                           # add to average
d.ave();
d.show();
mapRefScan := d.gsget();                # get the reference scan

# now for all scans, do (sig - ref)/ref calibration 
for (i in scans) {d.SRcal( i, mapRefScan); d.plotc(i)};

# now examine every scan
#for (i in scans) {d.plotc(i)};

# prepare to start making an image
include 'regionmanager.g'

# produce an image.
# the parameters below use channels 60 to 950, averaging 5 channels
# to produce an image cube 90 by 90 by 178 pixels.
d.imagems(msName,60,950,5,imname='w3-4.im',cellx='4arcmin',celly='4arcmin',nx=90,ny=90,gridfn='PB')

#next use the default viewer and click the refresh button.
#next click on "Raster Image" and after a few seconds an image will appear
#adjust to your hearts content, then print!