How to Make Quick Contour Map from OTF data

Introduction

So, you've run otftool to create a CLASS file and you are thinking "What does my data actually look like?" This cookbook will guide you through looking at spectra and making an integrated intensity contour map of your data. I try to assume little previous knowledge of CLASS or GILDAS, but I don't spend a lot of time explaining the commands or syntax, so you may want to have a copy of the CLASS Manual handy for reference.

I already know how to use CLASS, just show me how to make a map!


To start CLASS type class at the unix shell prompt:

> class
I-SIC_GTLGTR, No user defined logical name table
I-SIC, X-Window mode active
Default macro extension is .class
LAS>

At this point I find it useful to load in a short macro file defining some abbreviations which will save typing later and set up the display as I like it. You can borrow abb14.class, my abbreviation file for the FCRAO 14m telescope if you like. Download it and put it in the directory in which you run CLASS.
So, to load the abbreviations file type:

LAS> @abb14

In general macros can be run by preceeding the macro name with an @ sign.

Next you need a graphics display window:

LAS> dev xauto
I-X, X-Window default version $Revision: 4.49 $ $Date: 1997/10/22 13:23:46 $
I-X, Device DISPLAY :0
I-X, Screen characteristics 16 planes, TrueColor
I-X, Backing Store ALWAYS supported
I-X, Vendor The XFree86 Project, Inc; release 3360
I-X, Protocol X11 ;revision: 0
I-X, Screen size (mm) : height: 347, width: 433
I-X, Screen pixels : height: 1024, width: 1280
W-X, B&W Window, dithering used for images
Graphics device /dev/tty successfully defined as a XAUTO WHITE
LAS>

And a blank "GreG" display window should open.

You can use the help command to get more information on any of the CLASS commands. E.g.:

LAS> help dev

To load in your data, first set the filename containing your map:
This will be the name of the file that you created with otftool.

LAS> file in s171_13co.fcrao
I-CONVERT, File is [Native]
I-INPUT, s171_13co.fcrao successfully opened
LAS>

Next "index" the spectra:

LAS> find
I-FIND, 2752 observations found
LAS>

List them if you like (not a great idea if you have a lot ;)

LAS> list
2748; 1 S171 13CO_1-0 FCRAO-COR -425.000 +1.03E+03 Eq 3915
2749; 1 S171 13CO_1-0 FCRAO-COR -450.000 +1.03E+03 Eq 3915
2750; 1 S171 13CO_1-0 FCRAO-COR -475.000 +1.03E+03 Eq 3915
2751; 1 S171 13CO_1-0 FCRAO-COR -500.000 +1.03E+03 Eq 3915
2752; 1 S171 13CO_1-0 FCRAO-COR -525.000 +1.03E+03 Eq 3915
[snip]

The number on the far left gives the "observation number." You can use this to select individual spectra with the "get" command. E.g.

LAS> get 2751
I-GET, Entry 2751 Observation 2751; 1 Scan 3915

You can look at the header:

LAS> header
2751; 1 S171 13CO_1-0 FCRAO-COR O: 07-FEB-2002 R: 07-FEB-2002
RA: 00:01:23.000 DEC: 68:17:59.00 (1950.0) Offs: -500.000 +1.03E+03 Eq
Unknown Tau: 8.8505E-02 Tsys: 234.5 Time: 0.1163 El: 59.00
N: 1024 I0: 512.5 V0: -10.00 Dv: 6.6415E-02 LSR F0: 110.201370 Df: -2.4414E-05 Fi: 0.00000000
B_ef: 0.000 F_ef: 0.000 G_im: 0.000
H2O : 4.000 Pamb: 988.7 Tamb: 280.3 Tchop: 0.0 Tcold: 0.0
Tatm: 0.0 Tau: 0.000 Tatm_i: 0.0 Tau_i: 0.000 3465- 3524, 3886- 3915,

and plot the spectrum:

LAS> plot
S-CHAR, Fonts loaded

and hopefully you now see your spectrum and some header information displayed in the GreG window...

However more useful is to select spectra based on position. For this use the "find /offset" command:

LAS> find /offset 0 0
I-FIND, 1 observation found

get the spectrum, and plot it:

LAS> get first
I-GET, Entry 968 Observation 968; 1 Scan 3915
LAS> plot

and if you're really lucky your plot should look something like this:

example spectrum

In order to make an integrated intensity map you need to decide on the integration limits. This is why we needed to know how to look at a spectrum before we can make a map. It is also quite useful to have a quick glance at all the spectra in the file before moving on. A good way to do this is with a "stamp" map.

Before you can do this you have to specify the part of the spectrum that you want to look at. In this case the range -20 -> 0 km/s seems a good choice. You also need to specify the y-axis range:

LAS> set mode x -20 0
LAS> set mode y -0.5 4.5

Make sure you have all spectra indexed:

LAS> find

then clear the plotting surface and make the "map"

LAS> clear
LAS> map match /grid

You can use the command "pop" to pop-up individual spectra from the "map."

LAS> pop

then click the middle mouse button in the centre of a spectrum in the map and a second GreG window should open containing a blown-up version of that spectrum. To exit "pop-up" mode type "e" in the GreG window containing the stamp-map and you should get a prompt back in the text console.


Making an Integrated Intensity Contour Map

If you're new to CLASS/GILDAS and have ploughed through all the stuff above then this bit should seem pretty easy! Below I document two methods for producing a contour map from a CLASS file.

Here I'll just list the commands I use in order. Anything after a "!" on a line is a comment.

Method 1

> class

LAS> @abb14 ! or your macro file of choice
LAS> dev xauto ! open a plotting device
LAS> file in fname.fcrao ! open the input file
LAS> find ! index spectra in that file
LAS> find /offset 0 0 ! find the spectrum at the central position
LAS> get first ! get the spectrum
LAS> plot ! plot it

Examine the spectrum to determine the velocity range within which you want to integrate, then:

LAS> find ! get all the spectra back into the index
LAS> print area v1 v2 /output fname.tab

makes an output ascii table containing 4 columns: (1) observation number, (2,3) the position offsets in current angle units, (4) the integrated area between velocities v1 and v2.

At this point you can use the table to make a map in any of your favourite plotting packages. However, I find GreG makes nice contour plots, and is significantly quicker than say, PGPLOT, if you really just want a quick look at your data.

So, exit CLASS, and start up GreG:

> greg
GreG> dev xauto ! open a graphics device
GreG> clear plot ! clear the plotting device
GreG> set plot 20 25 ! set up the plotting window
GreG> set font duplex ! set duplex fonts
GreG> set expand 0.8 ! set the expand factor

The above are the parameters I use and seem to work pretty well...
GreG> column x 2 y 3 z 4 /file "fname.tab" ! read in the data from the table file, taking the x values from column 2, y values from column 3 etc
GreG> limits /reverse x ! calculate the x and y limits and reverse the x-axis (usually the case for celestial co-ordinate systems).
GreG> set box match ! make the x and y aspect ratios equal
GreG> box ! draw a box with labelled axes
GreG> levels 0.2 to 2 by 0.2 ! set the contour levels
GreG> rgdata /blanking -1000 ! grid the data, assuming that the z column contains the values for eaxh x,y position and the data is on a regularly sampled grid. Set any "blank" pixels to -1000.
GreG> rgmap /grey 1 2 ! draw a greyscale map
GreG> rgmap ! overlay contours
GreG> label "RA Offset (arcmin)" /x ! label the x-axis
GreG> label "Dec Offset (arcmin)" /y ! label the y-axis
GreG> hard fname.ps /dev eps grey ! make a postscript copy.

You can put all the above greg commands into a macro, and save yourself a lot of typing ;) The macro name should have the extension .greg and then you can call it with:
GreG> @macroname ! you do not need the extension if it is .greg

As an example, here is my map.greg. Feel free to download it and adapt it as you wish.


Method 2

In class, 'find' the index you want (see method 1 if you don't know how to do this), then:

LAS> grid fname.gdf new tdv(v1,v2) /image beam 49

This produces a file fname.lmv which contains an integrated intensity image, integrated between v1 and v2.

There are a few handy mathematical functions like "tdv" which can be applied as options to the grid command. Use the help command:

LAS> help cube

or refer to the manual for a description.

Now exit CLASS and start up graphic:

> graphic
GRAPHIC> image fname.lmv! loads in your image file
GRAPHIC> extrema /compute ! compute the extrema of your image
GRAPHIC> limits /rgdata ! set the plot limits
GRAPHIC> set box match ! make the x and y aspect ratios equal
GRAPHIC> box /u m ! draw a labelled box
GRAPHIC> plot ! plots half-tone plot
GRAPHIC> wedge ! draws a scale wedge on the right of plot
GRAPHIC> levels 30 to 100 by 10 ! set the contour levels
GRAPHIC> rgmap /per ! draws contours, assuming the levels are percentages of the maximum value

Again you could create a short macro (this time with the extension .graphic) to save yourself some typing. I'll leave this as an exercise for the reader ;)


Copyright Naomi Ridge
Last modified: Thu Apr 4 15:37:04 EST 2002
Valid HTML 4.01!