There is an example function runningmjdmedian which computes a running median as a function of time instead of as a function of integration number. Invoke it as

foo=runningmjdmedian(mjd,z,dt=1.0/24.0)

You can bin the data as follows. Suppose you have an x array dmjd, and a y array tant, and want to average it down to 300 bins:

bindata,dmjd,tant,dmjdbinned,tantbinned,etant,300

The 3rd and 4th arguements are return values; the 5th argument (etant) contains an estimate of the errorbar on the binned y-values based on the scatter. The last argument is the number of bins you want. This number doesn't need to divide evenly into the original number of elements (though the last binned data point may end up with a much larger error bar in effect).

Data Editing Examples

You can subtract a running median of length nchunk integrations from an array as follows

newdata=runningmedian(myarray,n=nchunk)

You can "clip" the outlier values (those greater, say, than badvalue) from an array as follows:

goodindices=nonoutlierindices(mydata,clip=badvalue)
newdata=mydata[goodindices]

what this actually has done in practice is to identify elements of the original array which have good values (values < badvalue), and return a list of indices to that original array which satisfy that condition. These indices are then used to access the original array, returning an array which does not include the bad values.

NOTE: anything you do to the data array, you must also do to corresponding arrays like right ascencion and declination, or else the data will not be properly registered with their coordinates!!

If you have Radio Frequency Interference in your data, you will want to do something like this:

• to subtract a very narrow-buffer running median over the data to get rid of all structure on scales of the beam or larger-- that's the real sky structure-- leaving only things which vary too fast to be astronomical: presumably the interference.
• then identify elements of the median subtracted data which correspond to RFI events. You can use the good-elements array indices to access the original data structure.

For the astronomical analysis you will want to run a longer-timescale running median over the data.

Elements of the Data Structure

The main elements of the gbt data structure of interest are

• dmjd -- modified julian date , a number with units of days
• calibdata[Nif, Ninteg] -- the data
• gbtra, gbtdec -- J2000 right ascencion and declination, both in units of degrees, where the telescope was pointed at a given integration.

mydata=(*1).calibdata[1,*]

mydata=(*q).calibdata[1,*]

Note that the GBT is at latitude 38:25:59.236 degrees north.

Groups are as follows

Group 1 observes monday night; group 2 tuesday night; group 3 wednesday night. all sesssions are 7-10pm.

Advisors are Larry Morgan (lmorgan # nrao.edu ); Brian Mason (bmason # nrao.edu); Chris Salter (csalter # naic.edu); and Jim Condon (jcondon # nrao.edu)

GBT MOTION CONSTRAINTS: Telescope will be at 192.441 Az and can be moved between 30 and 80 degrees elevation. We must restrict the number of elevation changes to 3 within a single observing session.

plot,(*q).calibdata[0,*]-median((*q).calibdata[0,*]),/ynozero

; First put ra and the first IF's data into standard IDL arrays myra=(*q).gbtra mydata=(*q).calibdata[0,*] ; then plot them plot,myra,mydata-median(mydata),/ynozero ; load a color table and overplot the other IF's data

oplot,(*q).calibdata[1,*]-median((*q).calibdata[1,*]),color=55

mydata=(*1).calibdata[1,*] oplot,myra,mydata-median(mydata),color=55

The calibdata array will be the basis of your analysis. You can apply standard IDL routines to it. See, for example

The calibdata array will be the basis of your analysis. You can apply standard IDL routines to it. DOcumentation for standard IDL routines, and for the routines in the NASA's very useful IDL astronomy library, can be found at:

Data analysis will be a little easier if you use a scan duration that's long enough that more than 5, ideally more like 10, beams will have drifted through in one scan. This is the last (third) argument of Track().

It will be best if your individual scans have many beams-- how long will this be?. This is the last (third) argument of Track(), given in seconds.

Reminder to the advisors: when observing, we need to have the antenna in "monitor only" mode. From Frank Ghigo: in the "Drives" screen, select "monitor only", "doScan", and select the correct optics (prime focus or Gregorian), and correct receiver. in the "pointing" screen, enable the pointing model. Note that when using the L-band receiver, the subreflector should be moved to X,Y,Z = 0 position, i.e, not in the stowed position used for prime focus.

The calibdata array will be the basis of your analysis. You can apply standard IDL routines to it. See, for example

http://www.physics.nyu.edu/grierlab/idl_html_help/idl_alph.html and http://idlastro.gsfc.nasa.gov/homepage.html

%META:TOPICINFO{author="BrianMason" date="1183752960" format="1.0" version="1.1"}% %META:TOPICPARENT{name="HandsonProjects"}%

Continuum Hands-On Project

1. Background

The objectives of the project, and the conceptual background, are described here

2. Observing

You will use the accounts cont1...contN. Use this and the password given to log in to the unix command line.

Start astrid, the observing application. Select project SDCONT in the drop-down window.

Data analysis will be a little easier if you use a scan duration that's long enough that more than 5, ideally more like 10, beams will have drifted through in one scan. This is the last (third) argument of Track().

3. Data Reduction

Copy into your home directory the file /users/bmason/sdcontidl (or just invoke it from the command line as ~bmason/sdcontidl)

This is a wrapper that starts up regular IDL (not GBTIDL) and points the path variables at some scripts that are handy for analyzing continuum data.

Here is an example session with commentary

; point IDL at the project & session name you got your data
setpipeopts,gbtproj='SDCONT_02'
; make an index of scan number to filename. You use this index
;  to access your data
indexscans,si
; print out a text summary of your observations
summarizeproject,myproj

; read scan #17 into a variable called q
readdcrdrift,si[17],q
; calibrate the data in q assuming an on-off and the Tcal value in the receiver calibration files
calibdcr,q
; show what IFs are present in the data
;  note that the first column in this printout is 
showifs,q

; plot the data... slightly awkward syntax is needed since q is a pointer to a structure
;  it's all standard IDL though.
plot,(*q).calibdata[0,*]-median((*q).calibdata[0,*]),/ynozero
loadct,12
oplot,(*q).calibdata[1,*]-median((*q).calibdata[1,*]),color=55

-- BrianMason - 06 Jul 2007