IDL procedures


There are a large number of IDL utility procedures for the analysis of CUTLASS and other data which are available for use and are stored in subdirectories in /people/cutlass/idl. Many of these routines are used by the IDL data visualisation Widgets (Plot_Radar, Rawful, etc.) and also by command line driven analysis packages (go, plot_image, plot_merge etc.). A subset of these routines (an expanding subset, hopefully) are documented here. Documentation on the main analysis programs themselves and the specialist routines therein, are documented separately, in the appropriate links above

This hypertext is version based on the APL original. Constructive criticism to yxo


CONTENTS


STARTUP.PRO
Device definition, compilation of useful routines
SD_DEVICE_INIT.PRO
Device functions and utilities
RADAR DATA PROCEDURES
IDL Functions and procedures for use with radar data
DATA STRUCTURES
Radar data structures defined within IDL
GENERAL FUNCTIONS
Functions implemented by genlib.pro
Raw Data Functions
Functions implemented by rawlib.pro
Fit Data Functions
Functions implemented by fitlib.pro
ACF Functions
Functions implemented by acflib.pro
Plotting ACFs
plotting of autocorrelation functions etc. in acflib.pro
Other Procedures
This software has never been implemented at Leicester

STARTUP.PRO

Device definition, compilation of useful routines

The procedure STARTUP contains some features that are designed to make IDL easier to use, particularly when switching from one device to another. Most of the utilities designed to make device switching easier are described below in "SD_DEVICE_INIT.PRO." The procedure also compiles library routines needed for e.g. plot_radar use. The procedure file STARTUP provides the following capabilities:


SD_DEVICE_INIT.PRO

Device functions and utilities

The routines contained in this IDL procedure are used to switch from one device to another. It contains the following routines:

init_plot
This prints out a little bit of documentation.
init_default
Reinitializes the device to the default device defined by the new system variable !default.
set_default
Changes the default device to have the characteristics of the current device.
init_x
Changes the device to X-Windows, and loads the color table #12 using the library routine loadct.
init_4100
Changes the device to a textronix 4100 series color terminal (this is useful with some terminal emulators). NOTE: there is no procedure to initialize a plain taktronix 4012 series terminal, since all you have to do for that is issue the command set_plot,tek
init_ps [,/aspect][,/landscape]
Changes the device to a standard postscript printer. the printer is set up to use the Helvetica font at 14 points. the name of the output file is IDL.PS. normally, the printer is set up for portrait mode. If the keyword, /landscape, is included the printer will be set to landscape mode. If the keyword, /aspect, is included the aspect ratio is be set to 1:1, so that circles plot as circles. a local environment variable is also created containing the extension PS of the file name.
init_epsf [,/aspect]
Changes the device to standard postscript, but uses the encapsulated postscript format. This is useful if the output file is to be included in another document (like a word processor). the name of the plot file is IDL.EPSF. a local environment variable is also created containing the extension epsf of the file name.
output_ps [,quename] [,/color]
Closes the file plot file and sends the plot to the queue specified by the environment variable IDL_CPRINT or IDL_BWPRINT, Optionally, a single argument can be included to specify a different queue. For example, output_ps, "colorscript" will send the plot to the queue for the color postscript printer. If the keyword /color is used, the plot will be sent to the queue identified by the variable IDL_CPRINT (which at Leicester points to the HP DeskJet). NOTE: The device will be reset to the default device as defined by the system variable, !default, after the plot file has been sent to the printer.

Radar Data Procedures

IDL Functions and procedures for use with radar data

A suite of routines for reading and manipulating radar data exist in the libraries genlib.pro, rawlib.pro, acflib.pro and fitlib.pro . These routines contain several functions and procedures for reading radar data from both raw data files and .FIT files. The procedures often call external programs using the CALL_EXTERNAL mechanism. The following external programs are required:

Current external calls from IDL
Executable image Environment Variable Entry Point Purpose
libacf.so SD_LIB_ACF ckrng_idl Remove lags with cross range noise
libfit.so SD_LIB_FIT fitropen_idl Open .FIT and.INX files
libfit.so SD_LIB_FIT fit_close_idl Close .FIT and .INX files
libfit.so SD_LIB_FIT find_fit_rec_idl Find a record according the time
libfit.so SD_LIB_FIT read_fit_idl Read one logical data record from a .FIT file
libfit.so SD_LIB_FIT block_read Reads a block of data from FIT file
libfit.so SD_LIB_FIT beam_read Reads a beam of data from FIT file
libpgm.so SD_LIB_PGM cnvcoord_idl Convert from ACGM to geog. coords and vice-versa
libpgm.so SD_LIB_PGM mlt_idl Convert from UT time in seconds from the beginning of the year to MLT
librbpos.so SD_LIB_FIT rbpos_idl Calculate position given radar beam & range
librbpos.so SD_LIB_RBPOS geochrt_idl Calculate position and radar k vector
libmag.so SD_LIB_MAG magcmp_idl Returns the bx, by, and bz (local coords) at the point lat, lon, elev
libgen.so SD_LIB_GEN median_filter_idl Smooth a velocity map with a median filter
libraw.so SD_LIB_RAW rawropen_idl Open a raw (.DAT or .CMP) data file
libraw.so SD_LIB_RAW raw_close_idl Close a raw data file
libraw.so SD_LIB_RAW read_raw_idl Read a raw data record
libraw.so SD_LIB_RAW find_raw_rec_idl Find a raw data record according to time
libraw.so SD_LIB_RAW raw_file_info_idl Print raw data record information

Data Structures

Radar data structures defined within IDL

There are two data structures that are defined by rawlib.pro and fitlib.pro and each of these structures contains the substructure that defines the radar parameters. The structure for the raw data is defined in rawdef.pro as follows:

  r = {rawdata, $
        p: {parms, rev: {rev_no, major: byte(0), minor: byte(0)}, $
                nparm: 0, $
                st_id: 0, $
                year: 0, $
                month: 0, $
                day: 0, $
                hour: 0, $
                minut: 0, $
                sec: 0, $
                txpow: 0, $
                nave: 0, $
                atten: 0, $
                lagfr: 0, $
                smsep: 0, $
                ercod: 0, $
                agc_stat: 0, $
                lopwr_stat: 0, $
                radops_sys_ress: 0, $
                noise: long(0), $
                radops_sys_resl: lonarr(2), $
                intt: 0, $
                txpl: 0, $
                mpinc: 0, $
                mppul: 0, $
                mplgs: 0, $
                nrang: 0, $
                frang: 0, $
                rsep: 0, $
                bmnum: 0, $
                xcf: 0, $
                tfreq: 0, $
                scan: 0, $
                mxpwr: long(0), $
                lvmax: long(0), $
                usr_resL1: long(0), $
                usr_resL2: long(0), $
                cp: 0, $
                usr_resS1: 0, $
                usr_resS2: 0, $
                usr_resS3: 0}, $
        pulse_pattern: intarr(PULSE_PAT_LEN), $
        lag_table: intarr(2, LAG_TAB_LEN), $
        combf: bytarr(COMBF_SIZE), $
        pwr0: lonarr(MAX_RANGE), $
        acfd: lonarr(2, LAG_TAB_LEN, MAX_RANGE), $
        xcfd: lonarr(2, LAG_TAB_LEN, MAX_RANGE) $
        }

The structure for the fit data is defined in fitdef.pro as follows:

f = {fit, rec_time: long(0), $ p: {parms, rev: {rev_no, major: byte(0), minor: byte(0)}, $ nparm: 0, $ st_id: 0, $ year: 0, $ month: 0, $ day: 0, $ hour: 0, $ minut: 0, $ sec: 0, $ txpow: 0, $ nave: 0, $ atten: 0, $ lagfr: 0, $ smsep: 0, $ ercod: 0, $ agc_stat: 0, $ lopwr_stat: 0, $ radops_sys_ress: 0, $ noise: long(0), $ radops_sys_resl: lonarr(2), $ intt: 0, $ txpl: 0, $ mpinc: 0, $ mppul: 0, $ mplgs: 0, $ nrang: 0, $ frang: 0, $ rsep: 0, $ bmnum: 0, $ xcf: 0, $ tfreq: 0, $ scan: 0, $ mxpwr: long(0), $ lvmax: long(0), $ usr_resL1: long(0), $ usr_resL2: long(0), $ cp: 0, $ usr_resS1: 0, $ usr_resS2: 0, $ usr_resS3: 0}, $ pulse_pattern: intarr(16), $ lag_table: intarr(2,48), $ combf: bytarr(80), $ pad1:long(0), $ noise_lev: double(0.0), $ noise_lag0: double(0.0), $ noise_vel: double(0.0), $ pwr_lag0: dblarr(75), $ slist: intarr(75), $ nsel: 0, $ qflag: lonarr(75), $ pad2:long(0), $ pwr_l: dblarr(75), $ pwr_l_err: dblarr(75), $ pwr_s: dblarr(75), $ pwr_s_err: dblarr(75), $ vel: dblarr(75), $ vel_err: dblarr(75), $ width_l: dblarr(75), $ width_l_err: dblarr(75), $ width_s: dblarr(75), $ width_s_err: dblarr(75), $ stnd_dev_l: dblarr(75), $ stnd_dev_s: dblarr(75), $ stnd_dev_phi: dblarr(75), $ gscat: intarr(75), $ pad3: 0, $ x_qflag: lonarr(75), $ pad4: long(0), $ x_pwr_l: dblarr(75), $ x_pwr_l_err: dblarr(75), $ x_pwr_s: dblarr(75), $ x_pwr_s_err: dblarr(75), $ x_vel: dblarr(75), $ x_vel_err: dblarr(75), $ x_width_l: dblarr(75), $ x_width_l_err: dblarr(75), $ x_width_s: dblarr(75), $ x_width_s_err: dblarr(75), $ phi0: dblarr(75), $ phi0_err: dblarr(75), $ elev: dblarr(75), $ elev_low: dblarr(75), $ elev_high: dblarr(75), $ x_stnd_dev_l: dblarr(75), $ x_stnd_dev_s: dblarr(75), $ x_stnd_dev_phi: dblarr(75), $ num_lags: intarr(75)} The following give an example of using these structures.

rawdef,r fitdef,f status=readraw(r) print, "The version number is",r.p.rev.major,r.p.rev.minor t1=r.p.year t2=r.p.month t3=r.p day t4=r.p.hour t5=r.p.minut t6=r.p.sec status=readfit(f,inx_unit,fit_unit,t1,t2,t3,t4,t5,t6) print,"velocity at last range=",f.vel(f.p.nrang-1);


General Functions

Functions implemented by genlib.pro

cnvtime

	calling sequence:  t = cnvtime(yr,mo,day,hour,min,sec)

This function converts the time given as year, month, day, etc. to the time in seconds from the beginning of the year. It returns a 32-bit integer value.

cnvcoord

	calling sequence: outpos = cnvcoord(inpos, [/GEO]

OR


	calling sequence:  outpos = cnvcoord(inlat, inlong, 
				height, [/GEO]

This function converts the input position given in latitude, longitude and height, to an output position in latitude, longitude and height. The default conversion is from geographic to ACGM geomagnetic, but this can be changed to a conversion from magnetic to geographic co-ordinates by specifying the keyword /GEO.

mlt

	calling sequence:  mt = mlt(year, t, mlong)

This function converts from the time given in seconds from the beginning of the year (t) to a floating point value giving the magnetic local time. The magnetic longitude must be specified in the third argument.

rbpos

	calling sequence:  pos = rbpos(range, 
					[/CENTER], [/GEO]

OR


	calling sequence:  pos = rbpos(range, data=data_struct, 
						[/CENTER], [/GEO]

OR

	calling sequence:  pos = rbpos(range, height=h, station=st, 
		beam-=bm, lagfr=lag, smsep=sample_sep, [/CENTER],[/GEO]

This routine determines the co-ordinates of a beam and range cell. the position is normally returned in ACGM geomagnetic co-ordinates, unless the /GEO keyword has been specified (in which case the co-ordinates are given in geographic). the returned position is normally a 3 x 2 x 2 table, with the first index giving (0) latitude, (1) longitude, (2) radial distance from earth's centre. The 2nd and 3rd indices indicate the 4 corners of the range/beam cell:

range cell quadrilateral

If the keyword /CENTER is present, then the co-ordinates for the center of the beam/range cell are given and the returned value is only a 3 element array.

In the first form of the calling sequence, the only argument is the range number (first range is1). All the other data are taken from the structure Fit_data. This additional information includes the station id, the beam number, the lag to the first range and the range separation.

In the second from of the calling sequence, the additional data needed are taken from the data structure specified by the keyword "data." For example, to take the data from the structure used when reading raw data you would use: pos = rbpos(range,data=raw_data).

In the third form of the calling sequence any (or all) of the additional parameters can be specified by keywords. Using this form of the call does not require that any actual radar data be read into the data structure fit_data (or any other structure).

NOTE: The only required argument is the range. This argument can either be a single value, or it can be a vector containing a set of range gates. For example, to calculate all the positions of the centers for 50 ranges you could use:

range = indgen(50) + 1
pos = rbpos(range,/CENTER)

This would return an array pos(3,50). If the /CENTER keyword had not been used, the returned array would have been pos(3,2,2,50).

Raw Data Functions

Functions implemented by rawlib.pro

rawropen

	calling sequence: status = rawropen(filename)
					

This function attempts to open the raw data file specified by filename and the associated index file. If both files are successfully opened, status returns 1. If the index file cannot be opened, status returns 0, but the data file will still be open. The file will be searched for in the list of directories specified by environment variable SD_RAWROPEN_PATH. If the extension is not given, the extension .DAT is assumed for the data file. You may, however, specify an extension of .CMP in order to read Halley radar compressed raw data files.

If no filename is specified IDL will attempt to use the X-windows widget pickfile to select the file from a list of files. This will only work on an X-windows device.

raw_close

	calling sequence:  status = raw_close(fileptr)

The file specified by fileptr will be closed. If the close was successful, status will return with a value of 1. If no argument is specified, the last raw file opened will be closed.

findrawrec

	calling sequence:  status=find rawrec(t_req,[t_found], [fileptr];

This routine will locate the time specified by the first argument, t_req, (or the first record after the specified time) in the raw data file. The time t must be given in seconds of the year. If the argument t_found is included, the actual time found is returned in that argument. NOTE: you must define the argument t_found before you call the routine: For example:

t_found = 0L
status = findrawrec(t_req, t_found)

If the argument fileptr is included, the files specified by that file pointer will be read. If the fileptr argument is not used, then the most recently opened raw data file is assumed (this will be the files specified by the pointer in the common block, rawdata_com.

readraw

	calling sequence:  status = readraw([fp], [fdata]) 
			

fp is the file pointer

fdata is the raw data structure

This function reads a record of raw data. If the read was successful the status will be 1 on return.

If no arguments are specified, then the file specified by rawfileptr in the common block rawdata_com will be used.

rawdef

	calling sequence: rawdef,r
This procedure creates a new variable with the name specified by the argument and gives that variable the structure associated with raw data.

Fit Data Functions

Functions implemented by fitlib.pro

fitropen

	calling sequence:  status = fitropen(fname)
The function opens a .INX and .FIT file with the name specified by the first argument. The directory does not need to be specified. The file can exist in the directories specified by the environment variable SD_FITROPEN_PATH . You should specify only the file name and not the extension.

If no file name is specified, IDL will attempt to use the x-windows widget pickfile to select the file to be opened. This will only work on an x-windows device.

fit_close

	calling sequence:  status = fit_close(fileptr)
The fit and index files specified will be closed. If the optional argument is not included, the file pointer in the common block fitdata_com will be closed. If the argument is included the files specified by the argument are closed.

findfitrec

	calling sequence:  status - findfitrec(treq, [t_found], [fileptr])
This function locates a .FIT record by using the index file. The time t_req must be in seconds from the beginning of the year. If the optional argument t_found is specified, the actual time found is returned in that argument. If the optional argument fileptr is specified, the .INX file specified by this argument is used. Otherwise, the .INX file specified in common block fitdata_com is used.

readfit

	calling sequence: status = readfit([fp], [fdata])
			
This function reads a logical record from a .FIT file. If no arguments are specified, the next record in the file is read into the structure fit_data in the common block fitdata_com. If the argument f is specified, the data is read into the structure specified by the argument.

fitdef

	calling sequence:  fitdef,f
This procedure creates a new variable with the name specified by the argument and gives that variable the structure associated with .FIT data.

ACF Functions

Functions implemented by acflib.pro

pwrspec

	calling sequence: spec = pwrspec(acf)
This function is not normally called by the user. Instead this function is called by another function, acf_pwrspec.

This function determines the power spectrum associated with an autocorrelation function. The acf must be a complex vector. The return power spectrum will be a real vector. The length of the returned power spectrum will be a power of 2, but the acf can be of any size. The length of spec will be at least 4 times that of the original acf.

spectitle

	calling sequence: spectitle,range [,bm]
This procedure creates the plot titles to be used in a plot of a Doppler power spectrum. the range number is required. If the beam number is not specified, the value is taken from the structure raw_data in the common block rawdata_com. The routine sets the plot title in the form: "rng: xxx bm: xx time: yyyy/mm/dd hhmm:ss". the x title is "Doppler Velocity" and the y title is "Power".

vvect

	calling sequence: v = vvect(spec, [f, tau])
This function returns a vector that can be used as the x-array to be plotted when plotting a Doppler power spectrum. The argument spec is a spectrum created by the function acf_pwrspec. It is used to determine the size of the vector which is to be returned. If the value of f (the radar frequency) and tau (the fundamental lag separation) are not specified, they are taken from the appropriate values in the structure raw_data in the common block rawdata_com.

acf_pwrspec

	calling sequence: spec = acf_pwrspec(r1, r2 [,acfarray[)
This function uses the data from the structure raw_data in the common block rawdata_com to determine the Doppler power spectra for a sequence of ranges from r1 to r2. In the process of determining the spectrum, missing and bad lags have to be identified and filled in (by linear interpolation). If the variable acfarray is specified, the fixed up ACFs are returned in that variable. This routine calls the function pwrspec to actually calculate the power spectrum for a given (fixed up) ACF.

plot_spec

	 calling sequence: plot_spec, v, spec, r1 [,bm]
This procedure plots a sequence of spectra produced by the function acf_pwrspec. The vector, v, is determined by the function vvect. The argument r1 specifies the range number of the first range that is present (it should correspond to the argument r1 in the function acf_pwrspec). The beam number is optional. If it is not specified, the value from the data structure raw_data is used. The titles for the plots are created by the procedure spectitle.

ident_badlags

	calling sequence: ident_badlags,range
This procedure is called to identify all bad lags, including lags that are bad because, 1) a sample is lying on a transmitter pulse, 2) a sample is contaiminated by power from another range, 3) the lag is bad because its power is anomalous.

Plotting ACFs

routines for plotting autocorrelation functions etc. in acflib.pro

These routines in acflib.pro provide for the plotting of autocorrelation functions, the power and phase at each lag of an acf, and the lag0 power at all ranges. The following procedures are available.

plt_acf
plots the ACF
plt_pwr
plots the lag0 power of ACF
plt_acf_pwr
plots the ACF and the power of ACF at each lag
plt_both
plots both: the plt_acf_pwr and plt_phase
plt_phase
plots the phase

plt_acf, plt_pwr, plt_acf_pwr, plt_both and plt_phase always require the range as an argument and as presently implemented, always plot the data from the common block rawdata_com. For example,

		IDL>plt_both, 5
plots the power and phase for each lag of range 5.


Other Procedures

This software has never been implimented at Leicester, some source can be found in radarlib.pro, the rest in the APL ftp area

fplot

	calling sequence: fplot, x, y, valid

This procedure plots a set of data specified by the vectors x and y. the vector valid is used to indicate where gaps occur in the data. Valid ,must have the same dimensions as x and y. Data points corresponding to non-zero values of valid will be plotted, while data values corresponding to zeros in valid will be skipped. If the plot symbol is 0 (i.e. you are plotting lines), there will be gaps in the line plot.

foplot

	calling sequence: foplot, x, y, valid
This is the same fplot but it does not erase the current plot before plotting the new data.

acf_fit

calling sequence: A = acf_fit(range, cl, cs, lambda, sigma, omega)
This routine takes as input, a specified range, and produces as output: 1) the lag0 power for the lambda fit as determined from the parameters pwr_l and the noise_lev parameter. 2) the lag0 power for the sigma fit as determined from the parameters width_s and noise_lev. , 3) the decay parameter lambda determined from the width_lparameter, 4) the decay parameter sigma determined from the width_s parameter, 5) the omega parameter (angular frequency) associated with the Doppler velocity, and 6) the complex ACFs derived from these parameters:

acf equations

The return value of the function is A = cmplxarr(nlags, 2) where nlags is the number of lags in the acf and the second index refers to the lambda fit (0) or the sigma fit (1).

The functions and procedures below are used to read IMP-8 magnetometer data and plot these data. There are a number of functions which are used to convert from EBCDIC to ASCII and from IBM floating point format to VAX floating point format. These routines will not be defined further, since the user will normally not need to call them - they are used internally by the other procedures in this set.

The following procedures and functions are provided:

imp8_mag_read

calling sequence:  imp8_mag_read,filename,data
The file specified by the first argument will be opened and the data will be read into a structure with the name specified by the second argument. The file must be one of the IMP8 files produced by the NDADS: archive at the Goddard Space Flight Center. The data are read into a large structure with the following definition:

	{mag,
		year: 0L,
		day:  0L
		milli_sec: 0L,
		data_q: 0L,
		orbit_number: 0L,
		bit_rate_flag: 0L,
		Pseq_count: 0L,
		Fill: 0L
		Housekeeping: 0L,
		Field: flatarr(2)	;field strength
		F_latlon: fltarr(2),	;field lat & long in S/C coords
		variance: fltarr(6),
		num seq: 0L,
		num_detail: 0L,
		trajectory: lonarr(2),
		geomag_latlon: fltarr(2)	;this is the geomag. lat/long of S/C
		geo_se_pos: fltarr(3),	;this is the GSE position in Re
		radial_dist: 0.0,
		gsm pos: fltarr(2)	;Y & Z GSM position in Re
		sun_latlon: fltarr(2)	;geomag.lat & long of the sun
		moon_se_pos: fltarr(3)	;moonÕs position in GSE
		rot_se_sm: fltarr(3,3),	;rotation matrix for GSM to GSE
		rot_ci_se: fltarr(3,3),	;rot matrix for Celestial to GSM
		fill_2: fltarr(2),
		spin_angle: fltarr(2),
		f_se_sm_latlon: fltarr(2,2)	;field angles in GSE (0) & GSM (1)
		B_se_sm: fltarr(3.2) }	;field components in GSE (0) & GSM (1)

This structure is replicated enough times to accommodate the entire file, and the entire data set is returned in the second argument to the procedure.

read_imf


	calling sequence: data = read_imf(filename)
This function is used to read the files produced by the program IMPREAD, which reads the data from magnetic tape. This data file contains records which have the following information:
1) Seconds (time in seconds from beginning of day)
a) Bx, By, Bz (in GSM)
c) R, Y, Z (S/C position in GSM)
This data is read into a floating point array dimensioned to (7,nrec), which is returned as the value of the function.

imf_comp

	calling sequence: imf_comp, data,hr,bx,by,bz,r,y,z
This procedure decomposes the data in the first argument into the 7 vectors specified by the remaining arguments. The resulting vectors are: the time in hours, the magnetic field Bx, By, and Bz components in GSM coordinates, and the spacecraft position in GSM coordinates (radial distance, Y & Z)

plot_imf


	calling sequence:  plot_imf,hr,bx,by,bz [,mtitle]
This procedure plots the three components of the IMF on a single plot. The arguments are: the time in hours, the three IMF vectors, and an option title to be used on the plot.

plot_position

 calling sequence:  plot_position,hr,r,y,z [,mtitle]
This procedure plots the satellite position in GSM coordinates with an optional title.

plot_impdata

	calling sequence:  plot_impdata,filename,data, [mtitle] 
					[,/archive] [,/gse]
This procedure provides the quck and easy way of reading the file and plotting the data. The file specified by the first argument is opened and the data are read into the variable data. If the keyword /archive is specified the data file is assumed to be one of the IBM type files that comes from NDADS: otherwise the file is presumed to have been created with the IMPREAD program. If the data comes from the archive file, then the keyword /gse may be used to specify that you want to plot the data in GSE coordinates instead of the default GSM coordinates. The optional variable mtitle specifies a title for the plots.
Back to contents list
SuperDARN documentation