NAME:
acorexample2 - Looking at calibration data by month.
At the end of each month, all of the calibration scans for the
month are processesed and stored in an idl archive file. This includes
x102, x101,x113,gui calibrate button runs, etc.. . The processing includes
the standard mm_mueller0 processing (which applies the current mueller
matrix to correct the data). It does not run mueller2_5 (these routines
would recompute the mueller matrix elements).
The output of the processing is an array mm[] of {mueller}
structures (one for each pattern done). See mmrestore() for a
description of this structure.
The data is stored in the directory:
/share/megs/phil/x101/x102/runs
The files are:
cyymmdd1_yymmdd2.dat - these are ascii files with a list of the files
and scans within each file for the month.
yymmdd1 is the first date of the file,
yymmdd2 is the last date of the file
cyymmdd1_yymmdd2.sav - these are the idl save files that contain the
processed data.
rcvrN.dat - An ascii file containing a list of all of the scans for
receiver N for all of the data files.
As of 15apr02 the list of save files was:
c010101_010131.sav c010501_010531.sav c010901_010930.sav c020101_020131.sav
c010201_010228.sav c010601_010630.sav c011001_011031.sav c020201_020228.sav
c010301_010331.sav c010701_010731.sav c011101_011130.sav c020301_020331.sav
c010401_010430.sav c010801_010831.sav c011201_011231.sav c020401_020415.sav
Use the routine mmgetarchive to retrieve all or a subset of this data.
eg:
nrecs=mmgetarchive(010501,020501,mm,rcvnum=5)
You can then create a subset of the data with mmget.
You can plot the data with:
mmplotgtsb .. gain,tsys,sefd, beamwidth
mmplotcsme .. coma, sidelobes,mainbeam efficiency
mmplotpnterr .. pointing error
(See /pkg/rsi/local/libao/phil/Cor2/acorexample2.pro)
NAME:
arch_getdata - get cor data using the archive tbl arrays
SYNTAX: n=arch_getdata(slAr,slfileAr,indAr,b,type=type,incompat=incompat,
han=han,missing=missing)
ARGS:
slAr[l] : {sl} or {corsl} returned by arch_gettbl
slFileAr[m] : {slInd} returned by arch_gettbl
indAr[] : long indices into slar to return
KEYWORDS:
type : int 0 first rec (hdr and data) of each scan1
: 1 hdrs from first rec of each scan
: 2 all recs each scan (hdr and data)
: 3 average rec (hdr and data) each scan
: 4 all hdrs from scan
RETURNS:
l : int number of elements in b
b : depending on type it can be {corget} or {hdr}
slind[n]: long the index into slar for each element returned in b.
incompat[p] long indices in indAr that were not returned because the
datatype differs from that of the first record
missing[m] long indices in indAr that were not in the file location
that the database had recorded.
DESCRIPTION:
After using arch_gettbl() and possibly where(), call this routine
to read the header and data from disc. What you get back is determined by
the keyword type and ind. The number of elements in b[] can be greater
than the number of indices in ind[] (eg you asked for all of the records
of the scans, or you are returning just headers). The slind keyword
array has the same number of elements as b. It contains the index into
slAr for each elements of b.
EXAMPLES
;get all data for jan02->apr02 cband
nscans=arch_gettbl(020101,020430,slAr,slFileAr,rcv=9)
; get the start of all the on/off pairs with cal on /off
pattype= 1 ; on/off position with cal on,off
nfound=corfindpat(slar,indar,pattype=pattype)
n=arch_getdata(slar,slfilear,indar,bon,type=3,incompat=incompat)
n=arch_getdata(slar,slfilear,indar+1,boff,type=3,incompat=incompat)
; get just the headers from the on scan of each calon/off pair
pattype= 6 ; on/off position with cal on,off
nfound=corfindpat(slar,indar,pattype=pattype)
n=arch_getdata(slar,slfilear,indar,hdrcon,type=1,incompat=incompat,
slind=slind)
NOTE:
When returning just headers, each header of each board is returned
as a separate entry in b[]. Use slind to figure out which scan each
belongs to.
(See /pkg/rsi/local/libao/phil/Cor2/arch_getdata.pro)
NAME:
arch_getmap - input a map from the archive
SYNTAX: nstrips=arch_getmap(yymmdd1,yymmdd2,projId,srcNm,m,cals,$
calSpc,procNm=procNm,polABrd=polABrd,polBBrd=polBBrd,$
rcvNum=rcvNum,smpPerStrip=smpPerStrip,norev=norev,$
avgsmp=avgsmp,han=han,verb=verb,$
useslar=useslar,slar=slar,slfilear=slfilear
ARGS:
yymmdd1: long starting date in the archive to search:eg 030501 .
yymmdd2: long ending date in the archive to search:eg 030801 .
projId : string project id you used for the map (eg 'a1731').
srcNm : string source name used for the map.
KEYWORDS:
procNm: string procedure name used for taking data. Supported names
are: 'cormap1','cormapdec','cordrift'. The default is
'cormap1'
polABrd: int correlator board number to take polA (1 thru 4). default
is 1.
polBBrd: int correlator board number to take polB (1 thru 4). default
is 1. (This assumes two pols per board).
rcvNum: int restrict search to this receiver number. The default
is all receivers found.
smpPerStrip:long samples per strip in map. Default is to read it from
the first strip found. Use this parameter if the
first strip found is not complete.
norev: When you drive both directions, the routine will normally
reverse the odd strips so that the data,header, and cals of
adjacent strips will line up in RA. Setting norev to true
will override this default. The data will be returned in
the order it is taken. The default is false.
avgsmp: long Number of samples to average together. Handy if you've
oversampled each strip. If this number does not divide
evenly into the samples/strip, then samples at the end
will be dropped. avsmp=0 or 1 is the same as no averaging
The default is no averaging.
han: if set then hanning smooth the data on input. The default
is no hanning smoothing.
useslar: if set the routine will use the slar,slfilear passed in
by the user (rather than the default archive).
slar[]: {corsl} slar passed in by user (if useslar is set) or passed
back by this routine
slfilear[]: {slind} slfilear passed in by user (if useslar is set) or
passed back by this routine.
RETURNS:
nstrips: long The number of strips returned. This may be more than the
requested number of strips in the file if the user
redid some of the strips.
m[2,smpPerStrip,nstrips]:{} array containg the map info. There is 1 entry
per strip taken. See cormapinp() for a description of
this structure.
cals[n*nstrips]:{} cal info for each strip. n is 1 for 1 cal per strip
2 for 2 cals per strip (but note the warning below
for 2 cals per strip).
calSpc[n*nstrips]:{corget} if supplied then return the actual cal spectra
for each strip (calOn followed by calOff). Use this for
recomputing the cal values with a different mask.
DESCRIPTION:
On the fly mapping can be done with cormap1,cormapdec, or the cordrift
routines. Users could use cormapinp() or cormapinplist() to input the
map from one or more files. The routine arch_getmap() provides the
same functionality as cormapinplist() with the added benefit that it
locates all of the strips for the map from the data archive.
To identify the map data in the archive you need to enter:
1.The range of dates in the archive to search (yymmdd1,yymmdd2).
2.The project id used for the map.
3.The source name used for the map.
4.If the map was not taken with cormap1 you need to supply
the keyword procNm='cormapdec' or procNm='cordrift'
5.If the same source was mapped with more than one receiver (eg 430,lb)
then you need to enter the receiver number to use. If not it will
combine all of the data it finds. receiver numbers are:
1=327,2=430,3=610,5=lbw,6=lbn,7=sbw,8=sbh,9=cb,11=xb,12=sbn.
6.By default the routine will take polA, and polB from correlator board
1. You can change this with the keyword polABrd=2,polBBrd=3
arch_getmap searches the archive for the locations of the requested
dataset. It creates a list of files and starting scannumbers that it
then passes to cormapinp() (via cormapinplist). cormapinp() will stop
reading a file if it finds a strip with fewer than the requested number of
strips, or it has input the last strip of the map. arch_getmap() is
careful in searching the complete file looking for any strips that the
user has redone. These strips will be included in the returned data.
This makes it possible for the returned array m[] to have more strips
than the number requested in the map. You can find the stripnumber of
each strip from: stripNum=m[0,0,*].h.proc.iar[4]
EXAMPLE:
the source 'ca' was mapped by project a1731 at 430 using cormap1 (which
drives in ra) and cormapdec (which drives in dec). The data was taken
May through august of 2003. Input both maps.
src='ca'
polABrd=1
polBBrd=2
rcvnum=2 ; the 430 dome receiver
yymmdd1=030501
yymmdd2=030801
projId='a1731'
procNm='cormap1'
nstrips=arch_getmap(yymmdd1,yymmdd2,projId,srcNm,mra,procNm=procNm,$
polABrd=polABrd,polBBrd=polBBrd,rcvnum=rcvnum,$
/verb)
procNm='cormapdec'
nstrips=arch_getmap(yymmdd1,yymmdd2,projId,srcNm,mdec,procNm=procNm,$
polABrd=polABrd,polBBrd=polBBrd,rcvnum=rcvnum,$
/verb)
SEE ALSO:
cormapinp(),cormapinplist().
Note:
The routine will now search backwards for the first calon/off of the
map (if there was one).
You need to have an idea how much memory the maps will take so you
don't ask for more memory that your computer can provide. As
an example: the cormapdec map above used:
(1024chan)*(4bytes/chan)*(2pol)*(257smp/strip)*(65 strips) is 137 megaBytes
Remember that the routine will input all of the strips done (including
any repeated strips).
(See /pkg/rsi/local/libao/phil/Cor2/arch_getmap.pro)
NAME:
arch_getmapinfo - get info on maps from the archive
SYNTAX: nmaps=arch_getmapinfo(yymmdd1,yymmdd2,mapI,projId=projid,srcNm=srcNm,$
procNm=procNm,rcvNum=rcvNum,useslar=useslar,$
slar=slar,slfilear=slfilear
ARGS:
yymmdd1: long starting date in the archive to search:eg 030501 .
yymmdd2: long ending date in the archive to search:eg 030801 .
KEYWORDS:
projId : string if provided then limit the maps to this project id.
srcNm : string if provided then limit the maps to this source name.
procNm: string if provided then limit the maps to this data taking
procedure name. Supported names are:
'cormap1','cormapdec','cordrift'.
rcvNum: int restrict search to this receiver number.
useslar: if set the routine will use the slar,slfilear passed in
by the user (rather than the default archive).
slar[]:{corsl} slar passed in by user (if useslar is set) or passed
back by this routine
slfilear[]:{slind} slfilear passed in by user (if useslar is set) or
passed back by this routine.
RETURNS:
nmaps: long The number of maps found
mI[nmaps]:{} array of mapinfo structures, one foreach map
DESCRIPTION:
Mapping can be done with the datataking routines: cormap1,cormapdec,
and cordrift. The routine cormapinp() can be used to input a single
map if the user knows the filename and scannumber. The routine
arch_getmap() can be used to retrieve all map data given a number of
constraints (date range, projid, srcname, etc..). This data is returned
as an array of strips. For some projects, the amount of mapping data is
too large to recall at once. In this case, arch_getmapinfo() can be
used to get a summary of the maps that are in the archive, without
retrieving the actual data. You can than use the summary info
to decide which maps to input.
arch_getmapinfo() defines a map as a set of scans in a mapping procedure
that are contiguous in time and have increasing strip numbers. It returns
an array of map Info strucutures (one for each map found).
If a map has 120 samples per strip and 36 strips and the user
did one complete map and a second map of strips 1 through 32. This
routine would return a mapInfo structures with two entries.
The mapInfo structure contains:
** Structure <8266244>, 12 tags, length=1052, data length=1050, refs=1:
SRCNM STRING 'MRK335' Name of source
PROCNM STRING 'cordrift' Name of datataking procedure
SMPSTRIP INT 120 requested samples in a strip
REQSTRIPS INT 36 Requested strips in a map
RACENHR FLOAT 0.105417 Ra center in hours
DECCENDEG FLOAT 19.9528 Dec center in degrees
FIRSTSTRIP INT 1 First strip returned(count from 1)
NUMSTRIPS INT 36 number of strips returned
FNAME STRING '/proj/a1552/corfile.15nov01.a1552.1'
SCANNUM LONG 131900007 scan number 1st sample,1st strip
COMPLETEMAP INT 1 map is complete
HDR STRUCT -> HDR Array[1] hdr brd 1, 1st sample,1st strip
EXAMPLE:
the source 'MRK335' was mapped by project a1552 at lband using cordrift.
The data was taken nov01,dec01. To get the mapinfo data:
srcnM='MRK335'
projId='a1552'
procNm='cordrift'
yymmdd1=011101
yymmdd2=031231 ; search over a larger date range
nmaps=arch_getmapinfo(yymmdd1,yymmdd2,mI,projId=projId,srcNm=srcNm,$
procNm=procNm)
help,mI
MI STRUCT = -> Array[29] ; there were 29 maps
It turns out that the maps had two separate centers. You could find
them by looking at racenhr,deccendeg:
SEE ALSO:
arch_getmap(),cormapinp(),cormapinplist().
(See /pkg/rsi/local/libao/phil/Cor2/arch_getmapinfo.pro)
NAME:
arch_getonoff - get on/off -1 data from archive
SYNTAX: nfound=arch_getonoff(slAr,slfileAr,indar,bout,infoAr,$
incompat=incompat,han=han,scljy=scljy,verbose=verbose)
ARGS:
slAr[n] : {sl} or {corsl} returned by arch_gettbl
slFileAr[m] : {slInd} returned by arch_gettbl
indar[] : long indices into slar for scans to return
KEYWORDS: (for corposonoff)
han: if set then hanning smooth
scljy: if set then scale to jy using gain curve for date of scan.
if not set then return in Kelvins.
verbose:if set then call correcinfo for each scan we process
RETURNS:
nfound : long number of on/off-1 returned in bout
infoAr[2,n]: long info on each on/off -1 found. Values are:
info[0,n] status of the returned on/off's
1 - on/off -1 returned as requested.
2 - on/off -1 but only the off was in ind[]. This could
occur if you asked for all the scans within the radius
of some ra/dec and the off fell within the radius
but not the on. You might then look for a negative
going galaxy in the on/off-1.
<0- no cal was found so the units are Tsys.
returnes -1 or -2
info[1,n] this is the index into indar[] where this on/off
was found (could be the on or the off).
bout[n] : on/off-1 scaled to jy, K, or tsys (see info[]). There is
one elemetn for each on/off pair found
incompat[p] long indices in indAr that were not returned because the
datatype differs from that of the first on/off pair.
Since bout[] is an array, each element must have the same
format (eg. number of boards, lags/sbc, etc..).
DESCRIPTION:
After using arch_gettbl() and possibly where(), call this routine
to process all of the on,off position switch records included in the
subset slar[ind]. It will search for all of the on or offs in ind[] and
process them via corposonoff(). If only an ON or an OFF is found in ind[]
it will still try to process the pair. By default it will use the cal
to return the data in kelvins. The /scljy will return the data in janskies.
If no cal scan is found, then the data is returned in units of Tsys.
The infoar[2,*] returns information on the returned data. The first
index tells the format of the data, the second index is the ptr into
slar where the onPos was found.
EXAMPLES
;get all data for jan02->apr02 cband
nscans=arch_gettbl(020101,020430,slAr,slFileAr,rcv=9)
; select all the on/off pairs
pattype= 1 ; on/off position with cal on,off
nfound=corfindpat(slar,indar,pattype=pattype
n=arch_getonoff(slar,slfilear,indar,b,type=3,incompat=incompat)
NOTE:
The routine expects to find the on,off, and cals in the slar. You should
pass the entire slar with indar as the thing that does the subsetting.
Don't pass a subset of slar into this routine: eg..
arch_gettbl(010101,011230,slar,slfilear,rcv=6)
ra =131215.1
dec=101112.
dist=cmpcordist(ra,dec,slar=slar)
; find all points within 30 arcminutes of ra,dec from projid a1199
; ..wrong..
ind=where(dist[2,*] lt 30.,count)
slarn=slar[ind]
ind=where(slarn.projid eq 'a1199',count)
npts=arch_getonoff(slarn,slfilear,ind,/scljy)
; .. correct:
ind=where((dist[2,*] lt 30.) and (slar.projid eq 'a1199'),count)
npts=arch_getonoff(slar,slfilear,ind,/scljy)
(See /pkg/rsi/local/libao/phil/Cor2/arch_getonoff.pro)
NAME:
arch_gettbl - input table of scans from cor archive
SYNTAX: nscans=arch_gettbl(yymmdd1,yymmdd2,slAr,slfileAr,rcvnum=rcvnum,
freq=freq,proj=proj,cor=cor)
ARGS:
yymmdd1 : long year,month,day of first day to get (ast)
for interim correlator data yy is two digits.
for was data, yy is 4 digits 2004mmdd
yymmdd2 : long year,month,day of last day to get (ast)
for interim correlator data yy is two digits.
for was data, yy is 4 digits 2004mmdd
KEYWORDS:
rcvnum : long .. receiver number to extract:
1=327,2=430,3=610,5=lbw,6=lbn,7=sbw,8=sbw,9=cb,$
10=xb,12=sbn,100=430ch
freq[2] : float Require cfr of band to be between freq[0] and freq[1]
in Mhz. At least 1 sbc of scan must match this.
proj : string project number to match (eg. 'a1473')
cor : If set then return the cor scanlist (see below)
RETURNS:
slAr[count] : {sl} or {corsl} scanlist array for each scan
slFileAr[n] : {slInd} one entry per file found
nscans : long number of scans found
DESCRIPTION:
This routine will scan the correlator archive and return a table
(an array of scanlist structures) for all scans in the archive that meet
the requested criteria. Two arrays are returned:
1. slAr[] holds 1 entry per scan
2. slfileAr[m] holds 1 entry per file found.
-------------------------------------------------------------------
For Interim correlator data:
Each slAr entry is an {sl} or {corsl} structure that contains:
scan : 0L, $; scannumber this entry
bytepos : 0L, $; byte pos start of this scan
fileindex : 0L, $; lets you point to a filename array
stat : 0B ,$; not used yet..
rcvnum : 0B ,$; receiver number 1-16
numfrq : 0B ,$; number of freq,cor boards used this scan
rectype : 0B ,$;type of pattern used to take the data. The
values are listed below (see getsl() under
generic idl routines for the most recent list).
1-calon,2-caloff,3-posOn,4-posOff,5-coron
6-cormap1,7-cormapdec,8-cordrift,
9-corcrossch (used by calibration routines)
10-x111auto (rfi monitoring),11-one(bml),
12-onoffbml
... see getsl() in GENEral idl
numrecs : 0L ,$; number of groups(records in scan)
freq : fltarr(4),$;topocentric freqMhz center each subband
julday : 0.D,$; julian date start of scan
srcname : ' ',$;source name (max 12 long)
procname : ' '};procedure name used.
If the /cor keyword is included then the slAr will be an {corsl}
structure and will contain the following additional fields:
projId : '' ,$; from the filename
patId : 0L ,$; groups scans beloging to a known pattern
secsPerrec : 0. ,$; seconds integration per record
channels : intarr(4),$; number output channels each sbc
bwNum : bytarr(4),$; bandwidth number 1=50Mhz,2=25Mhz...
lagConfig : bytarr(4),$; lag config each sbc
lag0 : fltarr(2,4),$; lag 0 power ratio (scan average)
blanking : 0B ,$; 0 or 1
azavg : 0. ,$; actual encoder azimuth average of scan
zaavg : 0. ,$; actual encoder za average of scan
raHrReq : 0.D,$; ra hours requested , start of scan
decDReq : 0.D,$; dec degrees requested, start of scan J2000
Delta ScanEnd-scanStart.. real angle for motion
raDelta : 0. ,$; delta ra Arcminutes real angle
decDelta : 0. ,$; delta dec Arcminutes
azErrAsec : 0. ,$; avg azErr Asecs great circle
zaErrAsec : 0. $; avg zaErr Asecs great circle
-------------------------------------------------------------------
For was2 data:
The return slar data structure is:
a={slwas ,$
scan : 0L, $; scannumber this entry
rowStart : 0L, $; row in fits file start of scan 0 based.
fileindex : 0L, $; lets you point to a filename array
stat : 0B ,$; not used yet..
rcvnum : 0B ,$; receiver number 1-16, 17=alfa
numfrq : 0B ,$; number of freq,cor boards used this scan
rectype : 0B ,$;1-calon,2-caloff,3-posOn,4-posOff
numrecs : 0L ,$; number of groups(records in scan)
freq : fltarr(8),$;topocentric freqMhz center each subband
julday : 0.D,$; julian day start of scan
srcname : ' ',$;source name (max 12 long)
procname : ' ',$;procedure name used.
stepName : ' ',$;name of step in procedure this scan
projId : '' ,$; from the filename
patId : 0L ,$; groups scans beloging to a known pattern
secsPerrec : 0. ,$; seconds integration per record
channels : intarr(8),$; number output channels each sbc
bw : fltarr(8),$; bandwidth used Mhz
backendmode: strarr(8),$; lag config each sbc
lag0 : fltarr(2,8),$; lag 0 power ratio (scan average)
blanking : 0B ,$; 0 or 1
azavg : 0. ,$; actual encoder azimuth average of scan
zaavg : 0. ,$; actual encoder za average of scan
encTmSt : 0. , $; secs Midnite ast when encoders read
; start of scan
raHrReq : 0.D,$; requested ra , start of scan J2000
decDReq : 0.D,$; requested dec, start of scan J2000.
; Delta end-start real angle for requested position
raDelta : 0. ,$; delta ra last-first recs. Amins real angle
decDelta : 0. ,$; delta dec (last-frist)Arcminutes real angle
pntErrAsec : 0. ,$; avg great circle pnt error
; alfa related
alfaAngle: 0. , $; alfa rotation angle used in deg
alfaCen :
You can use the slAr and slFileAr to then access the actual data files
without having to search through the file.
EXAMPLES
The following examples will read the scanlist array for all the data
beteen jan02 and apr02. The examples will then select different pieces
of information that are in this scanlist array. The final step will
extract some of the actual datascans.
;get tbl for all data for jan02->apr02.
nscans=arch_gettbl(020101,020430,slAr,slFileAr)
; work with this scanlist dataset..
; make a list of the unique source names
srcnames=slar[uniq(slar.srcname,sort(slar.srcname))].srcname
print,srcnames
;
; find all of the on/off position switch data for cband (rcv=9)
; (note corfindpat does not yet work for was2 data..)
n=corfindpat(slar,indar,pattype=1,rcv=9)
; find source NGC2264
indar=where(slar.srcname eq 'NGC2264',count)
now extract the scan averaged data. It will extract all data scans
that match the dataformat of the first entry of indar[0] (see arch_getdata().
n=arch_getdata(slar,slfilear,indar,b,type=3,/han,incompat=incompat)
SEE ALSO:getsl, arch_getdata,corfindpat.
NOTES:
22aug04. this is now implemented for was2 data. The cor= keyword is
ignored. The data stucture returned is a little different than the
interim corr structure.
(See /pkg/rsi/local/libao/phil/Cor2/arch_gettbl.pro)
NAME:
arch_openfile - open a file using archive arrays
SYNTAX: stat=arch_openfile(slAr,slfileAr,ind,lun)
ARGS:
slAr[l] : {sl} or {corsl} {slwas} returned by arch_gettbl
slFileAr[m] : {slInd} returned by arch_gettbl
ind : long index into slar for file to open
RETURNS:
stat: int 1 opened ok, 0 could not open
lun: int or {desc} for the open file
DESCRIPTION:
arch_openfile() will open a file that contains a given scan.
The user inputs slar[], slfilear[] that are returned by arch_gettbl()
and an index into slar[] for the scan that you want. The routine will
then find the file that this scan is in and open it. It returns the
lun (or descriptor if this is a was file) that can then be used to
do i/o on the file. The routine does not read any data from the file.
When done be sure and free the lun/desc using free_lun,lun or
wasclose,lun
EXAMPLES
For was data..
;get table for data august 20 thru august 31 2004.
nscans=arch_gettbl(20040820,20040831,slAr,slFileAr)
; find scan number=423314477L
scan=423314477L
ind=where(slar.scan eq scan,count)
if count eq 0 then begin
print,'scan:',scan,' not found'
goto,done
endif else begin
istat=slar_openfile(slar,slfilear,ind[0],desc)
if istat ne 1 then goto,done
endelse
;
; now input the scan
;
istat=corinpscan(desc,b,scan=scan)
(See /pkg/rsi/local/libao/phil/Cor2/arch_openfile.pro)
NAME:
chkcalonoff - check that hdrs are a valid cal on,off.
SYNTAX: istat=chkcalonoff(hdrOn,hdrOff)
ARGS:
lun: assigned to file we are reading
hdrOn: hdr from cal on.
hdrOff: hdr from cal off.
RETURNS:
istat: 1 it is a valid cal onoff pair
: -1 The hdrOn passed in is not a cal on
: -2 The hdrOff passed in is not a cal off
DESCRIPTION:
Check that the two headers passed in belong to a calonoff pair.
The routine checks the procedure name and that car[*,0] contains 'on' and
'off' respectively.
EXAMPLE:
input 2 records and check that they are from a cal on,off
print,posscan(lun,32500007L,1)
print,corgetm(lun,2,bb)
istat=chkcalonoff(bb[0].b1.h,bb[1].b1.h)
(See /pkg/rsi/local/libao/phil/Cor2/chkcalonoff.pro)
NAME:
coraccum - accumulate a record in a summary rec
SYNTAX: coraccum,binp,badd,new=new,scl=scl,array=array,brd=brd
ARGS: binp[]: {corget} input data
badd: {corget} accumulate data here
KEYWORDS:
new : if keyword set then this is the first call, alloc badd.
scl[] : float if provided, then scale the binp data by scl before
adding. This can be used to weight data by g/t.
if scl is an array it should equal the number of output
boards requested (see brd keyword). If it is a scalar
then that value will be used for all the boards.
array : if set then treat badd as an array and add element wise
the array binp[n] to badd[n]. Use with the new keyword.
By default if badd is an array then binp will be added
element wise.
brd[] : int if specified then only process the specified boards.
numbering is 1,2,3,4.
DESCRIPTION:
Accumulate 1 or more records worth of data into badd. If keyword
/new is set then allocate badd before adding. The header values in
badd will come from the first record added into badd.
The badd.b1.accum variable will be incremented by numrecs(binp)*scl for
each call. When corplot is called, it will scale the data by 1./badd.b1.accum.
When calling coraccum with the new keyword, you can include the
array keyword. This will allocate badd to be the same dimension as
binp. All future calls using badd will add binp element wise to badd.
This can be used when accumulating multiple maps.
Accumulated data must be of the same type (numlags, numbsbc, bw,etc..).
If you have observations to accumulate with only partial overlap of the
data types, you can use the brd keyword to specify which boards to accum. The
accumlated data must still be of the same type.
Example:
print,corget(lun,b)
coraccum,b,badd,/new
print,corget(lun,b)
coraccum,b,badd
corplot,badd
; Add n scans together element wise:
for i=0,n-1 do begin
print,corinpscan(lun,b)
coraccum,b,bsum,new=(i eq 0),/array
endfor
;
; input an entire scan and then plot the average of the records
; (this can also be done directly by corinpscan).
print,corinpscan(lun,b,scan=scan)
coraccum,b,bsum,/new
corplot,bsum
;
; let scan 1 be:4 brds,2pol, 1024 lags
; let scan 2 be:2 brds,2 pol,512 lags followed by 2 brds,2 pol, 1024 lags.
; To accumulate brds 1,2 of scan 1 with brds 3,4 of scan 2:
print,corinpscan(lun,b,scan=scan1)
coraccum,b,badd,/new,brd=[1,2]
print,corinpscan(lun,b,scan=scan2)
coraccum,b,badd,brd=[3,4]
help,badd.b1,/st
** Structure <3957c8>, 4 tags, length=9200, refs=2:
H STRUCT -> HDR Array[1]
P INT Array[2]
ACCUM DOUBLE 2.00000
D FLOAT Array[1024, 2]
(See /pkg/rsi/local/libao/phil/Cor2/coraccum.pro)
NAME:
coracf - compute the acf from the spectra
SYNTAX: bacf=coracf(b,/spectra)
ARGS:
b[n] : {corget} spectra to convert
KEYWORDS:
spectra: if set, then the input is acf's and they want spectra.
RETURNS:
bacf[n]: {corget} spectra converted to acf
DESCRIPTION:
Convert the spectra to acf's.
(See /pkg/rsi/local/libao/phil/Cor2/coracf.pro)
NAME:
corallocstr - allocate an array of structures.
SYNTAX: barr=corstostr(b,len)
ARGS: b: {corget} structure to duplicate
len: long length of the array to allocate
RETURNS:
barr[n]: {corget} array of structures to return.
DESCRIPTION:
Use corallocstr to allocate an array of {corget} structures. This routine
is necessary since each corget structure returned is an anonymous structure;
(See /pkg/rsi/local/libao/phil/Cor2/corallocstr.pro)
NAME:
coravg - average correlator data.
SYNTAX: bavg=coravg(binp,pol=pol,max=max,min=min,norm=norm,onlypol=onlypol)
ARGS: binp[]: {corget} input data
KEYWORDS:
pol : if set then average polarizations together
onlypol : only average polarizations. if binp[] is an array
do not average together records.
max : return max of each channel rather than average
min : return min of each channel rather than average
norm : if set then normalize each returned sbc to a mean of 1
(used when making bandpass correction spectra).
RETURNS:
bavg : {corget} averaged data.
DESCRIPTION:
coravg() will average correlator data. It has 3 basic functions:
1. average multiple records in a array together.
2. compute the average of an accumulated record (output of coraccum)
After doing this, it will set the bavg.b1.accum value to be the
negative of what was there (so corplot does not get confused).
3. average polarizations. This can be done on a single record, or
data from steps 1 or 2 above. If the /onlypol keyword is set then
it will average polarizations but it will not average together records.
The data is returned in bavg.
Polarization averaging will average the two polarizations on all boards
that have two polarization data. If polarizations are on separate correlator
boards then the routine will average all boards that have the same setup:
nlags, freq, and bandwidth. It will not combine boards that have two
polarizations per board.
If polarization averaging is used then the following header fields will
be modified.
b.(i).h.std.grptotrecs to the number of boards returned
b.(i).h.std.grpcurrec renumber 1..grptotrecs
b.(i).h.cor.numsbcout 2->1 for dual pol sbc
b.(i).h.cor.numbrdsused to the number of returned brds after averaging.
b.(i).h.cor.numsbcout to the number of returned sbc in each brd after
averaging
b.(i).h.dop.freqoffsets in case we reduce the number of boars from
pol averaging.
If binp contains multiple records then the max or min keyword can be
used to return the max or minimum of each channel. This will not work
on data from coraccum since the sum has already been performed.
Example:
; average polarizations:
print,corget(lun,b)
bavg=coravg(b,/pol)
;
; average a scan and then average polarizations
;
print,corinpscan(lun,b)
bavg=coravg(b,/pol)
;
; average the data accumulated by coraccum
; then average polarizations
;
print,corinpscan(lun,b)
coraccum,b,baccum,/new
print,corinpscan(lun,b)
coraccum,b,baccum
bavg=coravg(baccum,/pol)
;
; return max value from each channel
;
print,corinpscan(lun,b)
bmax=coravg(b,/max)
; average polarizations but do not average records. An
; example would be where pfposonoff() processes all of the
; on/off pairs in a file (possibly different sources), and you
; want to average the pols of each source.
n=pfposonoff(file,bar,tsys,cals) ; bar[n]
baravg=coravg(bar,/onlypol) ; baravg[n] with pols averaged
(See /pkg/rsi/local/libao/phil/Cor2/coravg.pro)
NAME:
coravgint - average multiple integrations.
SYNTAX: bavg=coravgint(bin)
ARGS: bin[] : {corget} input data
RETURNS:bavg : {corget} averaged data
DESCRIPTION:
If you have input an entire scans worth of data (eg corinpscan), you can
use coravgint to compute the average of all of the records of the scan. It
returns a {corget} structure that is the average of the N records input.
The header returned will be from the first record of the data (bin[0]).
Example:
assume scan 32500029 has 300 1 second integrations. Then:
print,corinpscan(lun,bin,scan=32500029L)
bavg=coravgint(bin)
bavg will then contain the average of the 300 records.
NOTES:
Many routines will do the averaging for you automatically if you request
it (eg corinpscan,coronoffpos,..)
SEE ALSO:
cormedian - to compute the median rather than the average of a scan.
(See /pkg/rsi/local/libao/phil/Cor2/coravgint.pro)
NAME:
corbl - baseline a correlator data set
SYNTAX: istat=corbl(bdat,blfit,maskused,deg=deg,mask=mask,$
edgefract=edgefract,auto=auto,sub=sub,svd=svd)
ARGS:
bdat[n]:{corget} input data to baseline. n can be >=1
KEYWORDS :
deg: int user can pass in polynomial degree for fit.
mask: {cormask} user can pass in cormask structure to use for mask
(see cormask()).
edgefract[n]: float The fraction of each edge of the bandpass to not use
in the fit. if n=1 then use the same fraction for both
sides. If n=2 then use [0] for the left edge and [1] for
the right edge. If both mask= and edgefr= are specified
then mask= will be used.
auto: if set then automatically do the baselining from the
user's input data, don't query the user for the parameters.
sub : if set then return the data-blfit rather than the fit.
svd : if set then use svdfit rather than poly_fit for the
fitting. This works better for higher order polynomials
but it is slower.
RETURNS :
blfit[n]:{corget} the baselined fits or bdat-bfit if /sub set.
maskused:{cormask} the mask used for the fits.
istat : int 1 if all the fits done, 0 if one or more sbc not fit.
DESCRIPTION:
corbl will baseline a correlator data set. bdat can be a single
corget structure or an array of corget structures. By default the routine
will interactively baseline the spectra for each correlator board. The user
supplies the board, mask, and fit order from a menu.
On exit, the routine will return the fit values in a corget structure.
If the /sub keyword is set then the difference bdat-bfit will be returned.
The mask that was used will be returned in maskused. If bdat is an array
then the averages will be used for the plotting (data and fits) but the
fits will be done separately for each record.
The auto keyword will do the fit without prompting the user. In this
case you must enter the keywords:
deg= as the polynomial degree for fit. It will use the same for all boards.
You must also specify the mask to use by passing in mask or setting
edgefract to the fraction of the bandpass to ignore on each edge (the same
value will be used for all boards).
For interactive use the program menu is:
KEY ARGS FUNCTION
m .. define mask
c b .. copy mask from board b to current brd
f n .. fit polynomial degree n
h h1 h2 .. change horizontal scale for plot to h1,h2
v v1 v2 .. change vertical scale for plot to v1,v2
b brd .. switch to board..1->nboards
q .. quit
current board: 1
brdsLeftToDo: Brd1 Brd2 Brd3 Brd4
The plots are color coded as:
White : first polarization of this board
Red : 2nd polarization of this board (if there is one).
Green : fits to the data.
yellow: the mask that is used for the fit
An example of using the program is:
print,corbl(bdat,bfit,/usesvd)
1. Adjust the horizontal and vertical range of the plot with:
h h1 h2
v v1 v2
2. define the mask for the current board:
m
then use the left button to start,stop a range and the left button to
quit.
3. try fitting various orders to the data and mask:
f 3
when the fit has been done, an extra line is output giving stats of the fit:
FitInfo deg: 3 usesvd:1 mask%: 75.7 Rms: 0.0602 0.0581
We used a 3rd order fit, the mask is 75.7% of the entire range, and the
residual rms's are .06 and .058 in whatever units your spectrum is.
Trying a 7th order fit would give:
f 7
FitInfo deg: 7 usesvd:1 mask%: 75.7 Rms: 0.0461 0.0455
Here the residuals have decreased to .046 .
4. Move to another board:
b 3
5. define the mask for this board then fit.
6. repeat for all boards. To exit enter q
If you do not used the /sub keyword, then the fits are returned. The
baseline can then be removed using:
bldat=cormath(bdat,blfit,/sub)
Idl has 2 polynomial fitting routines: poly_fit that uses matrix inversion
and svdfit that uses singular value decomposition. poly_fit is faster while
svdfit is more robust for larger order fits. The fits are done in double
precision using the xrange 0. to 1.
SEE ALSO: cormask, blmask (in general routines),cormath
(See /pkg/rsi/local/libao/phil/Cor2/corbl.pro)
NAME:
corblauto - automatic baselining of a correlator data set
SYNTAX: istat=corblauto(bdat,blfit,maskused,coef,$
edgefract=edgefract,masktouse=masktouse,$
nsig=nsig,ndel=ndel,maxloop=maxloop,$
deg=deg,sub=sub,verbose=verbose,fsin=fsin,
double=double,raw=raw,plver=plver,$
ARGS:
bdat:{corget} input data to baseline.
KEYWORDS :
edgefract[n]: float The fraction of each edge of the bandpass to not use
in the computations. if n=1 then use the same fraction for
both sides. If n=2 then use [0] for the left edge and [1] for
the right edge. The default value is .08.
The keyword masktouse will override this option.
masktouse:{cormask} if provided,then use this rather than edgefract for the
starting mask. This allows to you to force certain
areas to not be included in the fit. see cormask() for
creating the mask. If the raw keyword is used, then this
is just a float array instead of a cormask structure.
nsig: int The rms loop will throw out all outliers greater than
nsig*sigma. By default nsig is 3.
ndel: int The number of extra channels to remove whenever we
throw out a channel because of large rms. def=1
maxloop: int When removing outliers in the rms loop, it will continue
looping until no points are removed or maxloop iterations
is hit. By default maxloop is 15.
deg: int maximum degree of polynomial fit.
The routine starts at order 1 and iterates the
rmsloop,fit iterating the degree after each pass
until the last fitorder done is deg. By default deg is 1.
sub : if set then return the data-blfit rather than the fit.
verbose : 1 plot status after each fit, stop after last fit of plot
2 plot status each fit, stop after each plot
3 plot status each rms loop. stop each plot
-1 same as 1 but don't wait..
fsin : int the order of cos(nx),sin(nx) to also fit.
double : if set then fit using double precision
raw: if raw keyword set, then the input/output is a float
array instead of a correlator data struct. Raw should
be set to the number of elements the the bdata array.
The maskused argument will be returned as if 1 brd, 1 pol
was used. This assumes that the data in bdat are
evenly spaced.
plver[2]:float vertical range [min,max ]to use for plotting bandpasses
(before subtraction). the default is autoscaling.
RETURNS :
blfit :{corget} the baselined fits or bdat-bfit if /sub set.
maskused:{cormask} the mask used for the fits.
coefInfo:{} structure holding the coef's for the fit of each sbc
istat : int 1 if all the fits done, 0 if one or more sbc not fit.
DESCRIPTION:
corblauto will automatically create a mask and then use it to baseline
a correlator data set. The input is a {corget} structure of data. It will
return the masks created, the coefs of the fits, and the fits. It will
optionally remove the fits from the input data and return the difference
data-fit. The verb= option will plot the data as the fitting is done.
The routine can be used for mask creation (just ignore the fits)
or mask creation and baselining.
The /raw keyword lets you input an array of floats rather than a
{corget} structure. It will process the data and return the mask,fits
as if there was one correlator board with 1 polarization.
The parameters that control the algorithm are:
edgefract: This is the fraction of points on each edge of the data that
are not used. The default is .08 . For a sbc 1024 points this
would ignore 82 points on each edge.
masktouse: This is a mask with 0's and 1s. The 0 parts will not be
used in the fits. See cormask() on how to create the mask.
deg: The maximum degree for the fit. The default is 1. Values larger
than 12 or 13 may have pivot errors in the fitting.
ndel: The number of extra channels to throw out on the left and right
of each new channel we delete. This insures that skirts of
lines don't get included in the mask.
nsig : The rms computation removes all points greater than nsig*sigma.
By default nsig is set to 3.
maxLoop: The rms removal loops until there are no points within Nsigma*sigma
or maxloop iterations have occured. By default maxLoop is 15.
---------------------------------------------------------------------------
ALOGORITHM:
0. compute the starting points to use from the edge fraction or the masktouse.
The excluded points will never be included in the mask selection or fits.
let x0 by the x points, y0 be the y points after this selection.
1. set fitOrder=0, set yfit=mean(y0)
2. set yresiduals= y0 - yfit
the rms loop starts here
3. compute the rms of yresiduals
4. remove all points in y with yresiduals gt nsig*sigma .
For each point removed, also remove ndel adjacent points (determined
by the ndel keyword).
5. if the number of points removed in 4 is gt 0 and we've looped lt maxloop
times then goto 3.
fitting .. we now have a set of points within nsig of the mean.
6. fitdeg=fitdeg+1
7. fit a polynomial of fitDeg and a harmonic function of order fsin
to the points left in y.
8. evaluate the fitted function for the entire dataset y0
yfit=poly(x0,coef).
9. If fitdeg lt keyword deg goto 2.
---------------------------------------------------------------------------
The coef structure returns the fit results:
{corcoef}:
pol[2,4]: int 1-> this sbc/pol is used, 0--> it isn't
deg : int of polynomial fit
fsin : int of harmonic fit
nlags[4]: int number lags each sbc
coefAr[(deg+1)+fsin*2,2,4]: float coef for each fit. the polynomial
coef are followed by the cos,sin coef.
rms[2,4]: float of the fitted region within the mask
maskFract[2,4]: float npointsmask/nlags
To use the coefficients to recompute the fit, the x values should go
from -pi to +pi If Nchn is the number of channels in the spectra then
x=((findgen(Nchn)/(Nchn-1.) - .5)*2.*pi
EXAMPLES:
1. Use a polynomial fit of 3, remove the fit from the data, plot the
data and stop after each sbc done:
istat=corblauto(bdat,bfit,maskused,coefinfo,deg=3,/sub,verb=1)
; plot the datafit with the mask overplotted
corplot,bfit,cmask=maskused
2. Use a polynomial fit of 2 and a harmonic fit of 2. Create a mask
to use before calling this. Plot the results without stopping for keyboard
input
cormask,bdat,mask
istat=corblauto(bdat,bfit,maskused,coefinfo,deg=2,fsin=2,/sub,$
masktouse=masktouse,verb=-1)
On return :
coefinfo.coefar[0:2, are the c0,c1,c2 of the polynomial fit and
coefinfo.coefar[3:4] are the cos(x) sin(x) amplitudes while
coefinfo.coefar[5:6] are the cos(2x) sin(2x) amplitudes
SEE ALSO: cormask, corbl
(See /pkg/rsi/local/libao/phil/Cor2/corblauto.pro)
NAME:
corblautoeval - evaluate corblauto fit
SYNTAX: yfit=corblautoeval(x,coef,sbc=sbc,pol=pol)
ARGS:
x[n]: float xvalue to fit . Full scale should be 0. to 1.
coefI: {} coef struct returned by corblauto
KEYWORDS:
sbc: int sbc to eval. default is 1. count 1..n
pol int pol to eval. default is 1. count 1,2
RETURNS:
yfit[n]: float the fit evaluated at x
DESCRIPTION:
corblauto fits a polynomial and harmonic function returning the
fit coeficients in coef. This routine will evaluate the fit at the
specified points x[n]. x should go 0 to 1. for full scale.
SEE ALSO: corblauto
(See /pkg/rsi/local/libao/phil/Cor2/corblautoeval.pro)
NAME:
corcalcmp - given calon,off buffers, compute scale to kelvins
SYNTAX: istat=corcalcmp(calOn,calOff,calScl,datToScl,han=han,$
useinpm=useinpm,bla=bla,blda=blda,edgefract=edgefract,$
mask=mask,caldif=caldif,calValAr=calValAr,Tsys=Tsys
ARGS:
calOn[n]: {corget} cal on data. If n is greater than 1, then
the n calOn recs will be averaged together.
calOff[n]: {corget} cal off data. If n is greater than 1, then
the n caloff recs will be averaged together.
KEYWORDS:
han: if set then hanning smoooth the cal data
useinpm: if set then user passed mask in via mask keyword.
bla: If set then corblauto will use the caldeflection
to compute the mask of channels to use
blda: If set then corblauto will be passed
(calon-caloff)/caloff . It will fit and compute
the mask using this value.
edgefract[2]:float fraction of the bandpass to ignore on each
edge. Default is .08 . If 1 number is provided,
it is will be used for both edges.
extraParms: any extra parameters will be passed to
corblauto (if /bla is set). In particular you
can used deg=,fsin=,verb= .
RETURNS:
istat: int 1 ok, 0 trouble getting cal value
calScl[2,numbrds]: float factor that converts correlator counts to
kelvins for each pol of each board [pol,brd].
[0,n] is the first pol of each board (whether
it is polA or polB). If a board has only 1 pol
then calScl[1,n] will be 0.
See WARNIG below..
datToScl[m]: {corget} If provided, this data will be scaled to
kelvins for you.
mask: {cormask} The mask used when computing the cal. This is
computed in the following order:
1. if bla set, then mask comes from corblauto
2. if mask is provided then it is used as
passed in:
3. if edge= provided, then mask is computed
ignoring the edgefraction specified.
4. edgefraction is set to .08 and then use 3.
keyword.
caldif: {corget} return calon-caloff in correlator counts.
calValAr: fltar(2,nbrds) return calvalues for each board. index
[0,*] is the first spectra of board, [1,*] is the
2nd spectra of the board.
Tsys: fltar(2,nbrds,2) return Tsys as measured from the spectra
[pol,nbrds,calOn/calOff]
DESCRIPTION:
Given the cal on and off records, compute the scaling factors to go from
correlator counts to kelvins. It returns the results in calScl. If the
argument, datToScl is provided, then scale this data to kelvins using
the computed values (it should be {corget} structures).
The routine computes (where .i are the individual channels):
calDif.i= calon.i - calOff.i
It will then average the cal over a section (or mask) of the bandpass.
The mask is derived in the following order:
1. /useinpm is set. The user passes in the mask
2. /bla is set. corblauto will determine the mask from the
caldeflection.
3. edgefract is provided. It will make a mask ignoring this fraction
of channels from each edge.
It then computes (for each pol on each brd):
calScl[sbc,brd]= calK/mean(calDefl.i )
where the mean is summed over the non-zero mask elements.
If datToScl is provided then it will be converted to kelvins and returned.
The mask used in the cal computation is returned (if you do a bandpass
correction you will probably want to normalize to this set of
channels).
EXAMPLE:
This routine can be used if the cal onoff was taken with a
non stanard routine. Suppose you have a cal Off followed by a
cal on record:
print,corgetm(lun,b,2,/han) ; assume b[0] is cal off ,b[1] is cal on
istat=corcalcmp(b[1],b[0],calscl,b,/bla,verb=-1,deg=4,fsin=4,mask=mask)
Notes:
Calon,caloff are assumed to be from the same calon/off measurement.
It uses the first entry of calOn to find the rcvr, freq to get the
cal values.
datToScl should be the same setup as calOn,calOff (eg same number
of boards, and pols per board.
WARNING:
calscl is returned as a fltarr(2,nbrds) where the first index
is the sbc per board. If the board has only 1 sbc, then
calscl[1,x] is 0. Be careful passing calscl to the cormath() routine
to do the multiplication. Cormath wants 1 scaler for each sbc. It does
not want the extra zeros.
the WRONG way..
bk=cormath(b,smul=calscl) ... This will fail if a brd has 1 sbc.
the CORRECT way..
calsclOk=calscl[where(calscl ne 0.)] .. gets rid of the zeros
bk=cormath(b,smul=calsclok) .. this will work ok..
(See /pkg/rsi/local/libao/phil/Cor2/corcalcmp.pro)
NAME:
corcalib - intensity calibrate a single spectra
SYNTAX: istat=corcalib(lun,bdat,bcal,bfit,scan=scan,calscan=calscan,$
calinp=calinp,datinp=datinp,avg=avg,maxrecs=maxrecs,$
han=han,sl=sl,edgefract=edgefract,mask=mask,$
bpc=bpc,fitbpc=fitbpc,smobpc=smobpc,blrem=blrem,svd=svd,nochk=nochk
ARGS:
lun : int file descriptor to read from
KEYWORDS
scan : long scan number for data. default is current position
calscan : long scan number for cal on scan. def: scan following
the data scan.
calinp[2]:{corget} pass in the caldata rather than read it from lun
datinp[n]:{corget} pass in the data rather than read it from lun
avg: if set then return the averaged source record
maxrecs: long maximum number of recs to read in. default is 300.
han: if set then hanning smooth the data.
sl[]: {sl} array used to do direct access to scan.
edgefract[1/2]: float fraction of bandpass on each side to not use during
calibration. default .1
mask:{cormask} mask structure created for each brd via cormask routine
use this rather than edgefract.
note: currently mask.b1[1024,2] will use the mask for the
first pol of each brd for both entries of the board.
bpc: int 1 band pass correct with cal off
2 band pass correct with calon-caloff
3 band pass correct (smooth or fit) with data spectra
The default is no bandpass correction
fitbpc: int fit a polynomial of order fitbpc to the masked
version of the band pass correction and use the
polynomial rather than the data to do the bandpass
correction. This is only done if bpc is specified.
smobpc: int smooth the bandpass correction spectra by smobpc channels
before using it to correct the data. The number should be
an odd number of channels. This is only done if bpc is
specified.
blrem: int Remove a polynomial baseline of order blrem. Fit to the
masked portion of the spectra. This is done after
any bandpass correction or averaging.
svd : If baselining is done (blrem) then the default fitting
routine is poly_fit (matrix inversion). Setting svd
will use svdfit (single value decomposition) which is
more robust but slower.
nochk : if set then don't bother to check the cal records
to see if they are valid (in case they were written
with a non-standard program.
RETURNS:
bdat: {corget} intensity calibrated data spectra
bcal: {corget} intensity calibrated cal spectra
bfit: {corget} if supplied then return the smoothed or fit data that was
used for the bandpass correction.
istat: 1 ok
: 0 hiteof
:-1 no cal onoff recs
:-2 could not get cal value
:-3 cal,data scans different configs
:-5 sbc length does not match mask length
DESCRIPTION:
For those people who do not do position switching, corcalib allows you
to scale spectra from a src only scan to kelvins and optionally bandpass
correct it. The routine uses a src scan and a cal on,off pair. The data can
be read from disc or input directly to this routine (via calinp, datinp).
On output bdat and bcal will be in Kelvins. By default the individual records
of the scan are returned. If the /avg keyword is used, the average of all
of the src records will be returned. If the bfit argument is supplied then
the fit or smoothed version used for the bandpass correction will also be
returned. It will be scaled to the median bdat value so you can overplot
them.
If the data is input from disc then lun should be the file descriptor
from the open command. By default it will start reading the src scan from
the current file position. The scan=scan keyword lets you position to the
src scan before reading. By default the calscans will be the two scans
following the src scan. The calscan=calscan keyword lets you position
to the cal on scan before reading them. If the scans on disc have more than
300 records you need to use the maxrecs= keywords so the entire scan will
be input.
By default 10% of the bandpass on each edge is not used for the calibration
( when computing the conversion factor:
CalInKelvins/ <Calon[maskFrqChn]-calOff[maskFrqChn]>).
You can increase or decrease this with the edgefract keyword. The mask
keyword overrides the edgefract keyword and allows you to use a mask for
each sbc (use cormask to create the mask before calling corcalib).
The calibration will then only use the channels within the mask when computing
the gain calibration factors. This mask can be used to exclude rfi or
spectral lines.
Bandpass correction can be done with the cal off scan, the
calon-caloff difference spectrum, the data scan, or not at all. These
can be divided into the data scan as they are (although dividing the
data scan into itself is not very interesting!) or you can smooth or
fit a polynomial to the bandpass correction spectrum and then use the
fit/smooothed spectrum for the bandpass correction.
The keyword fitbpc=n will fit an n'th order polynomial to the bandpass
selected by the bpc keyword. Only the area within the mask (or edgefraction)
will be used for the fit.
The keyword smobpc=n will smooth the bandpass selected by the keyword
bpc and use it to do the bandpass correction (n should be an odd number >= 3).
You can pass in the data and/or calscans directly by using the
datinp, calinp keywords.
THE PROCESSING:
let Src be the src spectral data
let CalOn be the calOn spectra
let CalOff be the calOff spectra
let < > average over selected channels.The names will then have
Avg appended to them (eg calOnAvg=<calOn>)
let IndFrq be the set of frequency channels selected to use for the
calibration
The calibration consists of:
1. choose indFrq (the channels to use) in the following order:
a. The mask for each board from the mask keyword
b. Use edgefract to throw out this fraction of channels at each edge.
c. Use an edgefraction of .1 (10%)
2. compute CalOnAvg=<calOn[IndFrq]>,calOffAvg=<calOffAvg[IndFrq]>
3. Scale to Kelvins using: CntToK=CalValK/(calOnAvg-calOffAvg)
CalOnK =calOn*CntToK
CalOffK=calOff*CntToK
SrcK =Src*CntToK
4. If band pass correction is selected (bpc=1 or 2) then:
bpc=1: bpcN= calOff/<calOff[IndFrq]>
bpc=2: dif= calOn - calOff
bpcN= (dif)/<dif[IndFrq]>
If fitBpc > 0 then
bpcN=polyfit order fitbpc to bpcN[IndFrq]
else if smobpc > 2 then
bpcN=boxcar smooth (bcpN, smooth=smobpc)
then
SrcK=SrcK/bpcN
When deciding on the mask or edge fraction to use, you should have
a region where the calon-calOff is relatively flat (no filter rolloff and
no rfi).
EXAMPLE:
Suppose we have the following scans:
corlist,lun
SOURCE SCAN GRPS PROCEDURE STEP LST RCV
B1641+173 210200235 5 on on 17:42:08 5
B1641+173 210200236 1 calonoff on 17:47:10 5
B1641+173 210200237 1 calonoff off 17:47:21 5
40.42+0.70 210200238 5 on on 18:02:53 5
40.42+0.70 210200239 1 calonoff on 18:07:56 5
40.42+0.70 210200240 1 calonoff off 18:08:07 5
To process the first two sets:
--rew,lun
--print,corcalib(lun,bdat,bcal); will process the first set
--print,corcalib(lun,bdat,bcal,/han); will process the 2nd set with hanning
To process the 2nd set directly with an edgefraction=.12:
--print,corcalib(lun,bdat,bcal,scan=210200238L,edgefract=.12)
To input the data first, interactively create a mask, and then process
the data with a mask
--print,corinpscan(lun,bdatinp,scan=210200238L,/han)
--print,corgetm(lun,2 ,bcalinp,/han) ; reads the next two records
--cormask,bcalinp,mask ; interactively creates the mask
--print,corcalib(lun,bdat,bcal,calinp=bcalinp,datinp=bdatinp,mask=mask)
Use the same cal for multiple data scans:
--print,corgetm(lun,2 ,bcalinp,scan=210200236L/han);
--print,corcalib(lun,bdat1,bcal1,calinp=bcalinp,scan=210200235L)
--print,corcalib(lun,bdat2,bcal2,calinp=bcalinp,scan=210200238L)
Use the cal off for the bandpass correction. Use a 3rd order polynomial fit
to the cal off for the bandpass correction.
--print,corcalib(lun,bdat,bcal,scan=210200238L,bpc=1,fitbpc=3)
The bandpass correction is a bit tricky and depends on what type of
observations you are doing. The integration time for the off is usually
a lot less than the on positions so you need to use either the bandpass
fit or smoothing. It would probably be a good idea to add an option for
the user to input a bandpass to use for the correction (from an off src
position).
SEE ALSO:
cormask,corbl,corlist
(See /pkg/rsi/local/libao/phil/Cor2/corcalib.pro)
NAME:
corcalonoff - process a cal on/off pair
SYNTAX: istat=corcalonoff(lun,retDat,calOnrec,scan=scan,calval=calval,
sl=sl,swappol=swappol)
ARGS:
lun: unit number to read from (already opened).
retData[nbrds]: {cals} return data as array of structures.
calOnRec: {corget}optional parameter. If present it is the first
calon rec (in case you already read it).
KEYWORDS:
scan: if set, position to start of this scan before read
calval[2,nbrds]: cal value polA,polB for each board
sl[]: {sl} scan list array returned from getsl. If sl is
included then direct access positioning is available.
swappol: if set then swap polA,polB cal values. This can be
used to correct for the 1320 hipass cable reversal or
the use of a xfer switch in the iflo.
RETURNS:
istat - 1 ok, 0 did not complete
DESCRIPTION:
Read the calOnOff records and the cal values for the receiver
and frequency. Return the cal info in the retData[nbrds] array of
cal structs (one entry in retData for each correlator board).
The stucture format is:
retDat[i].h: {hdr} from first cal on
retDat[i].calval[2]:float calval in Kelvins for each sbc of this board
retDat[i].calscl[2]:float scale factor to use when scaling from
correlator counts to Kelvins: calT/(calon-caloff)
h.cor.calOn will contain the total power calOn (correlator units)
h.cor.calOff will contain the total power calOff (correlator units)
calon-caloff is the total power of the cal.
If a particular board has only 1 sbc, then the data will be in the
first entry [0] (whether or not it is polA or polB).
The routine will average multiple dumps if necessary.
You can scale your spectra to Kelvins via
b.bN.d[*,sbcn]= retDat[N-1].calscl[sbcn]*b.bN.d[*,sbcn]
(See /pkg/rsi/local/libao/phil/Cor2/corcalonoff.pro)
NAME:
corcalonoffm - process a cal on/off pair with a mask
SYNTAX: istat=corcalonoffm(lun,m,retDat,spc,scan=scan,calval=caltouse,sl=sl,$
edgefract=edgefract,han=han,swappol=swappol)
ARGS:
lun: unit number to read from (already opened).
m: {cormask} holds masks for the data 0 or 1.
It is created by the cormask routine.
Note. if m.b1[1024,2] has two entries per board, it uses
the first entry for both pols on the board.
retData[nbrds]: {cals} return data as array of structures.
spc[2]: {corget} if provided then return the calon (spc[0]) and
caloff spectra (spc[1]) after converting them to K.
KEYWORDS:
scan: if set, position to start of this scan before read
caltouse[2,nbrds]: cal value polA,polB for each board
sl[]: {sl} returned from getsl routine. If provided then
position to scan will be direct access.
edgefract: float if provided then create the mask ignoring m. Ignore
edgefract*lagsSbc lags from each side of the bandpass.
Return the mask in m (whatever is there gets overwritten).
han: if set then hanning smooth the spectra
swappol: if set then swap polA,polB cal values. This can be
used to correct for the 1320 hipass cable reversal or
the use of a xfer switch in the iflo.
RETURNS:
istat - 1 ok, 0 did not complete
DESCRIPTION:
Read the calOnOff records and the cal values for the receiver
and frequency. Compute the cal using the non zero lags in the mask m.
If edgefract is provided, then create the mask ignoring whatever is in m.
Zero edgefract*lagsbc lags from each side of the bandpass. Return this mask
in m (overwriting whatever was there).
The cal values and the scale factors (correlator counts to Kelvins) are
returned in the data structure retDat[n] (1 for each correlator board).
The data stucture format is:
retDat[i].h: {hdr} from first cal on
retDat[i].calval[2]:float calval in kelvin for each sbc of this board
retDat[i].calscl[2]:float to scale to kelvins calT/(calon-caloff)
h.cor.calOn will contain the total power calOn (correlator units)
h.cor.calOff will contain the total power calOff (correlator units)
calon-caloff is total power.
If a particular board has only 1 sbc, then the data will be in the
first entry [0] (whether or not it is polA or polB).
If the spc parameter is provided, the spectral data for the cal
on (spc[0]) and the cal off (spc[1]) will be returned after converting
them to Kelvins.
You can scale any other spectra to Kelvins using:
b.bN.d[*,sbcn]*= retDat[N-1].calscl[sbcn]
SEE ALSO:
corcalonoff
(See /pkg/rsi/local/libao/phil/Cor2/corcalonoffm.pro)
NAME:
corchkstr - check buffers for same structure
SYNTAX: istat=chbufstr(b1,b2)
ARGS:
b1: {corget} correlator data buf
b2: {corget} correlator data buf
RETURNS:
istat: 1 buffers are compatible
0 buffers are not compatible
DESCRIPTION:
Some routines store correlator data structures i{corget} in arrays
use corstostr. This is legal if the two buffers meet the following
requirements.
1. They have the same number of brds.
2. Each corresponding board has the same number of subcorrelators
3. Each sbc of each corresponding board has the same number of lags.
EXAMPLE:
istat=corget(lun,b1)
istat=corget(lun,b2)
same=corchkstr(b1,b2)
(See /pkg/rsi/local/libao/phil/Cor2/corchkstr.pro)
NAME:
corcmbsav - combine save files with multiple sources.
SYNTAX - npat=corcmbsav(savelistinp,b_all,b_sec,b_rms,b_info,usrcname,$
savefileout=savefileout,/filelist)
ARGS:
savelistinp[]: string list of save files to use as input
These files should contain:
bf[n]: an array of {corget} structures.
KEYWORDS:
filelist: if set then savelistinp is the data filelist used to input the
data. In this case this routine will generate the standard
save filenames to read:
corfile.ddmonyy.projid.n -> ddmonyy.projid.n.sav
It will assume that the .sav files are in the current directory.
savefileout: string name of save file to store the combined arrays in.
RETURNS:
npat : long number of patterns processed
b_all[]:{corget} the combined on/off -1 from the save files
b_sec[]:{corget} the combined secperchn from the save files
b_rms[]:{corget} the combined rms/mean from the save files
b_info[]:{} struct telling which kind of data is valid for
each b_rms/sec[m]
usrcname[l] string array holding the uniq source names present.
DESCRIPTION:
Process a set of idl save files each containing an array bf[n] of
corget structures. Combine these into a single array and return
them to the caller. optionally save them to an idl save file.
The routines pfposonoff(), pfposonoffrfi() process all of the
on/off patterns in a file (as long as they are the same data size).
They then save this info to an idl save file. This would typically be
1 days observations.
This routine's job is to combine all of the info from multiple
save files (days) concatentating them into 1 large array.
The input data in the save files are:
bf[n] {corget} holds the n on/off's processed for 1 file
arSecPerChn[n] {corget} secs/chan integration for each bf[n]. This is
output by pfposonoffrfi(). After rfi excision there may
not be the same integration time in each freq channel
arRmsByChn[n] {corget} rms/mean for each on/off after rfi excision
(actually after the rms fit, before the total power
test).
These will be combined into:
b_all[M] holds the combined bf[n]
b_sec[M] holds the combined arsecperchn[n]
b_rms[M] holds the combined arrmsbychn[n]
b_info[M] a struct holding source name for each entry and
whether or not there is b_sec and b_rms info for each
entry:
b_info[i].srcname ; source name
b_info[i].usesec ; 1 if it has secPerchn info from rfi processing
b_info[i].userms ; 1 if it had rms/mean info from rfi processing
The new arrays are returned to the caller. They can also be saved
to an idl save file with the savefileout keyword.
EXAMPLE:
Suppose project a1721 generated the following corfiles (either in
/share/olcor or /proj/a1721). You could process all the data with the
following code:
; The pfposonoff processing:
flist=['corfile.08may03.a1721.2',$
'corfile.08may03.a1721.3',$
'corfile.09may03.a1721.1']
dir=['/share/olcor/','/proj/a1721/'] ; directories to look for the files
voff=.03 ; plotting increment
han=1 ; hanning smooth
hor
ver
for i=0,n_elements(flist)-1 do begin &$
npairs=pfposonoff(flist[i],bf,tf,calsf,han=han,dir=dir,/scljy,/sav) &$
; plot the data so we can see what's up. offset each strip by voff
ver,-voff,(npairs+1)*voff &$
corplot,bf,off=voff &$
endfor
; The processed on/off-1 are now in the save files (one for each corfile)
; Use corcmbsav() to combine them:
; 1. input all the data from the daily save files.
; 2. save all of these arrays in the savbysrc filename.
; 3. the /filelist keyword flags that the string array has the
; original datafilenames and this routine should generate the
; input save file names.
;
savefileout='onoff_cmb.sav' ; place to store the new data
istat=corsavcmb(flist,b_all,b_sec,b_rms,b_info,$
savefileout=savefileout,/filelist)
;
WARNING:
This works as long as the data structure is the same type for
all of the on/offs. numboards, pol/brd, lags/pol
(See /pkg/rsi/local/libao/phil/Cor2/corcmbsav.pro)
NAME:
corcmpdist - compute distance (ra/dec) for scans
SYNTAX: dist=corcmpdist(ra,dec,b=b,slar=slar,hms=hms)
ARGS :
ra: float/double ra Hours J2000
dec: float/double dec Degrees J2000
KEYWORDS:
b[n]: {corget} array of correlator recs/scans to use for distance
sla[n]: {corsl} array of archive recs to use for distance
hms : if set then ra is in hms.s dec is in dms.s rather than hour,degrees
RETURNS:
dist[3,n]: float distance from ra,dec to each measurement.
[0,*] = ra - measured ra Arcminutes (greatcircle angle)
[1,*] = dec - measured dec arcminutes (angle)
[3,*] = total distance arcminutes (greatcircle angle)
DESCRIPTION:
Compute the distance from ra,dec to each measured entry. The measured
entries are passed in via b=b or slar=slar (you can use only one of these).
The requested ra,dec positions are used rather than the actual.
If b= is used then it is the requested position for each record.
If slar= is used then it is the average requested position of each scan.
The ra,dec arguments should also be in j2000 coordinates.
(See /pkg/rsi/local/libao/phil/Cor2/corcmpdist.pro)
NAME:
corcumfilter - cumfilter a correlator dataset.
SYNTAX: cmask=corcumfilter(b,limit,lengthfract,edgefract=edgefract)
ARGS:
b[n]: {corget} data to cumfilter
limit: float max deflection to allow
lengthfract:float fration of channels to use for filtering.
KEYWORDS:
edgefract: float fraction of channels to ignore on each side
(in case band pass still present)
RETURNS:
cmask[n]:{cormask} hold results of cumfiltering.
DESCRIPTION:
cumfilter is a routine written by carl heiles to remove outliers
from a dataset. The routine:
1. sorts the data by amplitude
2. selects length points (length=lengthfract*numpoints)
about the center of the sorted list
3. keeps all the points within limit*maxdeviation of the length
points from the center value.
The routine returns a cormask data structure that has a 0 for
points that were rejected and a 1 for points that were kept.
If b[n] is an array of structures then n cormask structures will
be returned (1 for each element of b).
EXAMPLE:
limit=2.
length=.5
print,corget(lun,b)
cmask=corcumfilter(b,rnage,length,edgefract=.08)
cmask is a single structure
You can then use cmask in basline fitting or operations:
istat=corbl(b,blfit,maskused,/auto,deg=1,mask=cmask,/sub)
print,corgetinpscan(lun,b)
cmask=corcumfilter(b,limit,length,edgefract=.08)
cmask is an array of cmask structures.
WARNING:
if b[] is an array , then corcumfilter will create and array of
masks. Some of the other cor routines (corplot,corstat, etc..) will
only accept single element cormasks so you need to be careful passing
any arrays of cmasks to other routines.
(See /pkg/rsi/local/libao/phil/Cor2/corcumfilter.pro)
NAME:
cordfbp - return the digital filter bandpasses.
SYNTAX: bp=cordfbp(b,force=force)
ARGS:
b: {corget} correlator data to make banpasses for.
KEYWORDS:
force: if set then force the recomputation of the bandpasses
rather than taking them from the common block if they are
already there.
RETURNS :
bp : {corget} holds the normalized digital filter bandpasses.
DESCRIPTION:
The normalized digital filter bandpasses are returned in a standard
{corget} data structure. If the input (b) is an array, then only
bandpasses for the first element of b are returned. The computed
bandpasses are stored in a local common block so that repeated calls with the
same type of data will go faster.
50 Mhz bandpasses are returned as all 1's (since there is no digital
filtering).
Stokes data returns the digital filter bandpasses for the
original PolA, polB data rather than I and Q. The last two polarized band
passes are returned as all 1's.
EXAMPLE:
; Input a data scan and divide each record by the digital filter bandpass.
; Note that the digital filter bandpasses are already normalized to unity.
;
istat=corinpscan(lun,b,scan=scan)
dfbp=cordfbp(b)
bpc=cormath(b,dfbp,/div)
The current filtering scheme seems to be:
filter on upconvert - all
filter 1 for interpolation after upconvert : double nyquist
filter 2 for interpolation after upnconvert: 12.5 Mhz and below
SEE ALSO: dfhbfilter, cormath
(See /pkg/rsi/local/libao/phil/Cor2/cordfbp.pro)
NAME:
corfindpat - get the indices for the start of a pattern
SYNTAX: nfound=corfindpat(sl,indar,pattype=pattype,rcv=rcv)
ARGS:
sl[]: {getsl} scan list array from getsl
KEYWORDS:
pattype: int Type of pattern we are looking for.
1 - on/off position switch with cal on/off
2 - on/off position switch whether or not cal there
3 - on followed by cal on ,off
4 - heiles calibrate scan two crosses
5 - heiles calibrate scan 1 or more crosses
6 - cal on,off
7 - x111auto with calonoff
8 - x111auto with or without cal
9 - heiles calibrate scan 4 crosses
If not provided than pattype 1 is the default.
dosl: lun if keyword dosl is set to a valid open file, then this
routine will do the sl=getsl(lun) call and return the
scan list in sl.
rcv : int if supplied then only find patterns for this receiver.
1..12 (see helpdt feeds) for a list of the feeds
RETURNS:
indar[npat]: long indices into sl[] for start of the pattern
npat : long number of patterns found
DESCRIPTION:
corfindpat() is used for the automatic processing of entire datafiles.
It processes a scanlist array (returned by getsl()) and returns an array of
pointers for the start of all the completed patterns of a particular type
that are located in the datafile. Patterns can be: on/off position switching,
heiles calibration runs, on with cal on/off,..etc.
The requirements for a completed scan depend on the pattern type:
type order needed requirements
1 posOn,posOff,calon,caloff. Number of on recs must equal number of
off recs
2 posOn,posOff. Number of on recs must equal number of
off recs.
3 posOn,calon,caloff.
4 calon,caloff,heiles calibrate scan with 2 crosses.
nrecs in each cross is the same.
5 calon,caloff,heiles calibrate scan with at least 1 cross.
nrecs in each cross is the same.
6 calon,caloff
7 x111auto (60recs) calon,caloff
8 x111auto (60recs)
9 calon,caloff,heiles calibrate scan with 4 crosses 120 samples/strip
nrecs in each cross is the same.
You can call sl=getsl(lun) once prior to this routine, or you can
set the keyword dosl=lun and it will create the sl array and return it.
It is also possible to limit the pattern to a particular receiver using
the rcv=rcvnum keyword.
EXAMPLE:
openr,lun,'/share/olcor/corfile.23aug02.x101.1',/get_lun
sl=getsl(lun)
; get poson,off,followed by cal on,off
npat=corfindpa(sl,indfound)
for i=0,npat-1 do begin
scan=sl[indfound[i]].scan
.. process this scan
; get ons followed by cal on,off. have the routine do the sl search.
; only get data for cband (rcv=9)
openr,lun,'/share/olcor/corfile.23aug02.x101.1',/get_lun
npat=corfindpa(sl,indfound,dosl=lun,pattype=3,rcv=9)
SEE ALSO: arch_gettbl,mmfindpattern,getsl (in general routines)
NOTE:
There is also a record type in the {sl} structure that lists the
name of the active pattern when the data was taken. It's coding and the
coding for this routine are not the same (sorry..). It may or may
not be accurate (for some test data, the pattern type is not set so the
last one supplied is used). Corfindpat differs in that it will
try and verify that the pattern is what it says it is and that it is
complete (eg calon followed by caloff).
(See /pkg/rsi/local/libao/phil/Cor2/corfindpat.pro)
NAME:
corfrq - compute the freq/vel array for a spectra
SYNTAX: retArr=corfrq(hdr,retvel=retvel,retrest=retrest)
ARGS:
hdr: header for this board
KEYWORDS:
retvel : if not equal to zero, then return velocity rather
than frequency array
retrest: if not equal to zero, then return the rest frequency
for the object being observered rather than the
topocentric frequency.
RETURNS:
retArr: Array of floating point frequencies or velocities.
DESCRIPTION:
Compute the topocentric frequency array (in Mhz) for the correlator board
corresponding to the header (hdr) passed in. If the keyword retvel
is set (not equal to zero) then return the velocity (optical definition).
If the keyword retrest is set, then return the rest frequency of the object
rather than the topocentric frequency.
The array returned (retArr) will have the same number of elements
as a single output sub correlator.
The order of the data assumes that the spectral channels are in
increasing frequency order (corget always returns the data this way).
If the spectra are spDat[2048] and then retAr[0] will be the lowest
frequecy or the highest velocity.
EXAMPLE:
.. assume 2 boards used, pola,b per board
corget,lun,b
frqArr=corfrq(b.b1.h)
frqArrRest=corfrq(b.b1.h,/rest)
velArr=corfrq(b.b1.h,/retvel)
plot,frqArr,b.b1.d[*,0]
(See /pkg/rsi/local/libao/phil/Cor2/corfrq.pro)
NAME:
corget - input next correlator record from disc
SYNTAX: istat=corget(lun,b,scan=scan,noscale=noscale,sl=sl,han=han)
ARGS:
lun: logical unit number to read from.
RETURNS:
b: structure holding the data input
istat: 1 ok
: 0 hiteof
:-1 i/o error..bad hdr, etc..
KEYWORDS:
noscale : if set, then do not scale each sub correlator to the
9level corrected 0 lag.
scan : if set, position to start of scan first..
han : if set, then hanning smooth the data
sl[] : {sl} array used to do direct access to scan.
This array is returned by the getsl procedure.
DESCRIPTION:
Read the next group of correlator data from the file pointed to
by lun. If keywords scan is present, position to scan before reading.
A group is the data from a single integration. Each correlator board
is written with a separate hdr and data array. The structure returned will
contain 1 to 4 elements depending on the number of boards being used.
b.b1
b.b2
b.b3
b.b4
each bN will have:
b.b3.h - complete header for this board
b.b3.p[2] int. 1-polA,2->polB, 0, no data this sbc
b.b3.accum double . accumulate scale factor (if used)
b.b3.d[nfreqchan,nsbc] the data. nsbc will be 1 or 2 depending on
how many sbc are using in this board.
use pol to determine what pol each sbc is. It will also tell you if
there is only 1 sbc pol[1] = 0. It will not compensate for
zeeman switching..
The header will contain:
.h.std - standard header
.h.cor - correlator portion of header
.h.pnt - pointing portion of header
.h.iflo- if,lo portion of header
.h.dop - doppler frequency/velocity portion of header
.h.proc- tcl procedure portion of header
The data is returned in increasing frequency order as floats.
If an i/o error occurs (hit eof) or the hdrid is incorrect (you are not
correctly positioned at the start of a header), then an error message
is output and the file is left positioned at the position on entry to the
routine.
EXAMPLE:
.. assume 2 boards used, pola,b per board (lagconfig 9)
istat=corget(lun,b)
b.b1.h - header first board
b.b2.d[*,0] - data from 2nd board, polA
SEE ALSO:
posscan,corgethdr
(See /pkg/rsi/local/libao/phil/Cor2/corget.pro)
NAME:
corgethdr - return the correlator header for the next group
SYNTAX: istat=corgethdr(lun,rethdr)
ARGS:
lun: int opened file to read
rethdr[nbrds]: {hdr} return array of headers brd 1.. thru brd n
(See /pkg/rsi/local/libao/phil/Cor2/corgethdr.pro)
NAME:
corgetm - input multiple correlator records to array.
SYNTAX: istat=corgetm(lun,ngrps,bb,scan=scan,sl=sl,han=han,noscale=noscale,$
check=check)
ARGS:
lun: int open file to read
KEYWORDS:
scan: long if provided, position to scan before reading.
sl[]: {sl} returned from getsl(). If provided then direct access
to scan is provided.
han: if set then hanning smooth
noscale: if set, then do not scale each sub correlator to the
9level corrected 0 lag.
check : if set, then verify that the records are compatible to store
in a single record.
RETURNS:
bb[]:{corget} array of corget structures that were read in
istat: 1 got all of the groups requested.
0 did not get all of the groups requested.
DESCRIPTION:
Read ngrp consecutive correlator records from disc and return them
in bb (an array of {corget} structures). The routine will read from the
current position on disc (after optional positioning to scan). It will
read across scan number boundaries if needed.
(See /pkg/rsi/local/libao/phil/Cor2/corgetm.pro)
NAME:
corhan - hanning smooth correlator data.
SYNTAX: corhan,b,bsmo
ARGS:
b[m]: {corget} data to hanning smooth
bsmo[m]: {corget} return smoothed data here. Note.. if only 1 argument
is passed to the routine then the data is smoothed in
place (it is returned in b).
DESCRIPTION:
corhan will hanning smooth the data in the structure b. If a single
argument is passed to the routine, then the smoothed data is returned
in place. If two arguments are passed (b,bsmo) then the data is returned
in the second argument.
EXAMPLE:
print,corget(lun,b)
corhan,b ; this smooths the data an returns it in b.
print,corget(lun,b)
corhan,b,bsmo ; this smooths the data an returns it in bsmo
(See /pkg/rsi/local/libao/phil/Cor2/corhan.pro)
NAME: corhcalrec- check if an input rec is a cal rec. SYNTAX: istat=corhcalrec(hdr) ARGS : hdr - header from 1 of the boards RETURNS: istat- 0 not a cal rec, 1-on, 2-off DESCRIPTION: Check if a record input is part of a cal record. EXAMPLE: corget(lun,b) istat=corhcalrec(b.b1.h)
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhcalval - return the pol A/B cal values for a sbc.
SYNTAX: stat=corhcalval(hdr,calval,date=date,swappol=swappol)
ARGS:
hdr: {hdr} header for board to check
KEYWORDS:
date[2]: intarray [year,dayNum] if provided, then compute the calvalues
that were valid at this epoch.
swappol: if set then swap the pola,polb cal values. This can
be used to correct for the 1320 hipass cable switch
problem or the use of a xfer switch in the iflo.
RETURNS:
calval[2]: float .. calValues in deg K for polA,polB
stat: int .. -1 error, 0 got the values ok.
DESCRIPTION:
Return the cal values in degrees K for the requested sbc. This
routine always returns 2 values (polA then polB) even if the header
is for a board that uses only one polarization.
The calvalues for the receiver in use are looked up and then the
values are interpolated to the observing frequency.
EXAMPLE:
input a correlator record and then get the calvalues for the
3rd correlator board:
print,corget(lun,b)
istat=corhcalval(b.b3.h,calval)
.. calval[2] now has the cal values in degrees K for polA and polB.
NOTE:
Some cals have measurements at a limited range of frequencies (in some
cases only 1 frequency). If the frequency is outside the range of measured
frequencies, then the closest measured calvalue is used (there is no
extrapolation).
The year daynum from the header is used to determine which set of
calvalue measurements to use (if the receiver has multiple timestamped
sets).
This routine computes the frequency of the sbc from hdr and then calls
calget().
SEE ALSO:
gen/calget gen/calval.pro, gen/calinpdata.pro
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhcfrrest - return the rest freq of band center.
SYNTAX: cfr=corhcfrrest(hdr)
ARGS:
hdr: {hdr} .. header for board to check
RETURNS:
cfr: double .. rest center freq of band in Mhz
DESCRIPTION:
Return the rf rest frequency for the band center of the requested board.
EXAMPLE:
get the rest freq of the 3nd board: cfrMhz=corhcfrrest(b.b3.h)
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhcfrtop - return the topocentric freq of band center.
SYNTAX: cfr=corhcfrtop(hdr)
ARGS:
hdr: {hdr} .. header for board to check
RETURNS:
cfr: double .. center freq of band in Mhz
DESCRIPTION:
Return the topocentric rf center of the band for the requested board.
hdr can be an array of hdrs but they must all come from the same
board. eg b[].b1.h
EXAMPLE:
get the freq of the 2nd board: cfrMhz=corhcfrtop(b.b2.h)
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhdnyquist - check if rec taken in double nyquist mode
SYNTAX: istat=corhdnyquist(hdr)
ARGS : hdr{} - header from 1 of the boards of the record
RETURNS: istat- 0 not taken in double nyquist mode.
1 taken in double nyquist mode.
DESCRIPTION:
Check if a record was taken in double nyquist mode.
EXAMPLE:
corget(lun,b)
istat=corhdnyquist(b.b1.h)
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhflipped - Obsolete..check if current data is flipped in freq.
SYNTAX: stat=corhflipped(corhdr)
ARGS:
corhdr[]: {corhdr} to check.
RETURNS: istat- 0 increasing frq, 1- decreasing freq order.
if corhdr[] is an array then istat will be an array of ints.
DESCRIPTION:
This routine has been replaced by corhflippedh(). corhflipped
was using the bit in the header that was not always being set correctly
for 1 ghz bandwidths.
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhflippedh - check if current data is flipped in freq.
SYNTAX: stat=corhflippedh(hdr,sbcnum)
ARGS:
hdr[]: {hdr} to check.
sbcnum:long sbc to check 1 thru 4.
RETURNS: istat- 0 increasing frq, 1- decreasing freq order.
if corhdr[] is an array then istat will be an array of ints.
DESCRIPTION:
Check if the correlator data for this sbc is stored in increasing
or decreasing frequency order (even or odd number of high side lo's).
This routine replaces Corhflipped. Corhflipped gave incorrect results for
1 ghz if's if the LO was below 1500 Mhz and it was being used as hi side, or
the lo was above 1500 Mhz and it was being used as a low side lo.
corhflippedh should work in all cases.
EXAMPLE:
sbcnum=1
check first board: istat=corhflippedh(b.b1.h,sbcnum)
(See /pkg/rsi/local/libao/phil/Cor2/corhflippedh.pro)
NAME:
corhgainget - return the gain given a header
SYNTAX: stat=corhgainget(hdr,gainval,date=date,az=az,za=za,onlyza=onlyza)
ARGS:
hdr[n]: {hdr} header for board to check
KEYWORDS:
date[2]: intarray [year,dayNum] if provided, then compute the gain value
at this epoch.
az[n]: fltarray If provided, use this for the azimuth value rather than
the header values
za[n]: fltarray If provided, use this for the zenith angle values rather
than the header values
onlyza: If set then return the za dependence (average of az)
RETURNS:
gainval: float .. gainvalue in K/Jy
stat: int -1 --> error, no data returned
0 --> requested freq is outside freq range of fits.
Return gain of the closed frequency.
1 --> frequency interpolated gain value returned.
DESCRIPTION:
Return the telescope gain value in K/Jy for the requested sbc.
The gain fits for the receiver in use are input and then the
values are interpolated to the az, za and observing frequency.
If hdr[] is an array then the following restrictions are:
1. each element must be from the same receiver and at the same
frequency (eg. all the records from a single scan).
2. If the az,za keywords are provided, they must be dimensioned the same
as hdr
EXAMPLE:
input a correlator record and then get the gain value for the
3rd correlator board:
print,corget(lun,b)
istat=corhgainget(b.b3.h,gain)
.. gain now has the gain value in K/Jy
input an entire scan and compute the gain for all
records of 1 sbc. assume 300 records..
print,corinpscan(lun,bar)
istat=corhgainget(bar.b3.h,gain)
gain is now a array of 300 records
NOTE:
Some receivers have measurements at a limited range of frequencies (in some
cases only 1 frequency). If the frequency is outside the range of measured
frequencies, then the closest measured gain is used (there is no
extrapolation in frequency).
The year daynum from the header is used to determine which set of
gain fits to use (if the receiver has multiple timestamped sets).
This routine takes the az,za, date, and frequency from the
header and then calls gainget().
If you input an array of corget recs , then they must all be from the
same sbc.
SEE ALSO:
gen/gainget gen/gaininpdata.pro
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhintrec - return integration time for a record
SYNTAX: secs=corhintrec(hdr)
ARGS : hdr{} - header from 1 of the boards of the record
RETURNS: secs:float seconds of integration for this record
DESCRIPTION:
Return the integration time for a record.
EXAMPLE:
corget(lun,b)
istat=corhintrec(b.b1.h)
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhsefdget - return the sefd given a header
SYNTAX: stat=corhsefdget(hdr,sefdval,date=date,az=az,za=za,onlyza=onlyza)
ARGS:
hdr[n]: {hdr} header for board to check
KEYWORDS:
date[2]: intarray [year,dayNum] if provided, then compute the gain value
at this epoch.
az[n]: fltarray If provided, use this for the azimuth value rather than
the header values
za[n]: fltarray If provided, use this for the zenith angle values rather
than the header values
onlyza: If set then return the za dependence (average of az)
RETURNS:
sefdval: float .. sefd value in Jy
stat: int -1 --> error, no data returned
0 --> requested freq is outside freq range of fits.
Return sefd of the closed frequency.
1 --> frequency interpolated sefd value returned.
DESCRIPTION:
Return the telescope sefd value in Jy for the requested sbc.
The sefd fits for the receiver in use are input and then the
values are interpolated to the az, za and observing frequency.
If hdr[] is an array then the following restrictions are:
1. each element must be from the same receiver and at the same
frequency (eg. all the records from a single scan).
2. If the az,za keywords are provided, they must be dimensioned the same
as hdr
EXAMPLE:
input a correlator record and then get the sefd value for the
3rd correlator board:
print,corget(lun,b)
istat=corhsefdget(b.b3.h,sefdval)
.. sefdval now has the sefd value in Jy
input an entire scan and compute the sefd for all
records of 1 sbc. assume 300 records..
print,corinpscan(lun,bar)
istat=corhsefdget(bar.b3.h,sefdAr)
sefdAr is now a array of 300 records
NOTE:
Some receivers have measurements at a limited range of frequencies (in some
cases only 1 frequency). If the frequency is outside the range of measured
frequencies, then the closest measured sefd is used (there is no
extrapolation in frequency).
The year daynum from the header is used to determine which set of
sefd fits to use (if the receiver has multiple timestamped sets).
This routine takes the az,za, date, and frequency from the
header and then calls sefdget().
If you input an array of corget recs , then they must all be from the
same sbc.
SEE ALSO:
gen/sefdget gen/sefdinpdata.pro
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhstate - decode status words for correlator header
SYNTAX: statInfo=corhstate(corhdr)
ARGS:
corhdr[]:{hdrcor}
RETURNS:
statInfo[]{corstate} decoded status structure
DESCRIPTION:
The correlator header contains various info that is encoded in bitmaps
(h.cor.state) This routine decodes this
bitmap and returns it in a structure. The input can be 1 or more
cor headers.
Example:
print,corget(lun,b)
corstate=corhstate(b.b1.h.cor)
help,corstate,/st ; to print it out
The returned structure format is:
a={corstate,$
dnyquist : 0 ,$; 1--> double nyquist
chiptest : 0 ,$; 1--> chip testmoded
blankOn : 0 ,$; 1--> radar blanking on
pwrcntInc : 0 ,$; 1--> power counter data included
relbitshift: 0 ,$; relative bit shift used
relbitshiftsgn: 0 ,$; 1-->relbitshift neg, 0--> relbitshift pos
rawacf : 0 ,$; 1--> raw acf (not bias removal
cmbacf : 0 ,$; 1--> corrected acfs
spc : 0 ,$; 1--> spectra
pack : 0 ,$; 1--> packed data (for decoder)
startImdSw: 0 ,$; 1--> start immediate software
startImdHw: 0 ,$; 1--> start immediate hardware
startTick1: 0 ,$; 1--> start 1 sec tick
startTick10: 0 ,$; 1--> start 10 sec tick
dmpDelTrgInt: 0 ,$; 1--> dump delay trig from interrupt
dmpDelTrgTick: 0 ,$; 1--> dump delay trig from tick
adjPwrStartScn: 0 ,$; 1--> adjpower start of scan (automatic)
adjPwrStartRec: 0 ,$; 1--> adjpower start of rece (automatic)
spcFlipped: 0 ,$; 1--> spectra is flipped on disc
isACor: 0 ,$; 1--> running as a correlator (0-->decoder)
totCntIncluded: 0 ,$; 1--> total counts included
calOff : 0 ,$; 1--> cal is off this rec
calOn : 0 ,$; 1--> cal is on this rec
complexDat: 0 ,$; 1--> complex data taken
level9: 0 ,$; 1--> 9 level data. 0--> 3 level
stokes: 0 ,$; 1--> stokes data
pwrCntI: 0 ,$; 1--> power count I included
pwrCntQ: 0 }; 1--> power count Q included
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corhstokes - check if record taken in stokes mode
SYNTAX: istat=corhstokes(hdr)
ARGS : hdr{} - header from 1 of the boards of the record
RETURNS: istat- 0 not taken in stokes mode. 1-->taken in stokes mode
DESCRIPTION:
Check if a record was taken stokes (polarization) mode.
EXAMPLE:
corget(lun,b)
istat=corhstokes(b.b1.h)
(See /pkg/rsi/local/libao/phil/Cor2/corhquery.pro)
NAME:
corimg - create an image of freq vs time for 1 sbc.
SYNTAX: d=corimg(b,brd,sbc,bscl=bscl,cmpscl=cmpscl,edgchn=edgchn)
ARGS:
b[n] : {corget} input data
brd : 1-4 which board to use
sbc : 1,2 in this board to use
KEYWORDS:
bscl : {corget} divide into each entry in b[n]. Usually comes from
coravgint().
cmpscl: int 1 compute mean ,use to flatten image (default)
2 compute median,use to flatten image
3 donot flatten image.
edg