;+
; NAME: TMAGZEN
;	
; PURPOSE:
;	To produce an array that masks the pixels that would be found
;       on an image within N degrees of the magnetic zenith at Eureka.
;
; CATEGORY:
;	
; CALLING SEQUENCE:
;       TMAGZEN,YEAR,MONTH,DAY,CAM,FILT,DEG,MASK,CCP,CRP [,SMALL=SMALL] $
;	      [,HEIGHT=HEIGHT]	
; INPUTS:
;       YEAR,MONTH,DAY: The year,month and day the image that the mask 
;		will be used on was created.
;       CAM,FILT: The camera and filter number of the image that the 
;		mask will be used on.
;       DEG: The number of degrees from the magnetic zenith at Eureka.
; KEYWORD PARAMETERS:
;       SMALL:	If set, then the mask will be a 256 by 256 BYTE array. 
;		Otherwise, a 512 by 512 BYTE array will be returned.
;       HEIGHT: If set, then the calculations will be performed using 
;		the parameters from the magnetic zenith vector at 300 
;		km. Otherwise the magnetic vector at a height of 1 km 
;		will be used.
; OUTPUTS:
;       MASK:	A BYTE array, that contains 1's were the pixels of an 
;		image would be withing DEG degrees of the magnetic zenith 
;		and 0's elsewhere.
;       CCP:	Returns the column number of the pixel that corresponds 
;		to the magnetic zenith.
;       CRP:	Returns the row number of the pixel that corresponds 
;		to the magnetic zenith.
; COMMON BLOCKS:
;       None.
; SIDE EFFECTS:
;       None.
; RESTRICTIONS:
;       None.
; PROCEDURE:
;       ....
; MODIFICATION HISTORY:
;       Written August 1994, by T A Oliynyk.	
;-

PRO tmagzen,year,month,day,cam,filt,deg,mask,ccp,crp,small=small,height=height

loc=pocaloc(year,month,day)
IF KEYWORD_SET(height) THEN loc(2)=height
igrf,loc(0),loc(1),loc(2),year,bx,by,bz,bhor,dip,dec
hmag=bhor*1E5	; convert gauss to nT
vmag=bz*1E5
rot=dec

IF KEYWORD_SET(small) THEN edge=256 ELSE edge=512

col=INDGEN(edge)#REPLICATE(1,1,edge)
row=TRANSPOSE(col)

IF NOT KEYWORD_SET(small) THEN small=0
; transform row and column numbers to zenith and azimuth angles
cr2az,year,month,day,cam,filt,col,row,az,za,small=small,site=3 $
     ,/degrees

lambda=(90*!DTOR-ATAN(hmag/vmag))  ; latitude of zenith magnetic vector
phi=(rot+180)*!DTOR                ; east longtitude of zenith magnetic vector

beta=(90-za)*!DTOR  ; zenith angle converted to latitude 
alpha=az*!DTOR      ; azmuth angle is already equal to east longtitude

; transformation produced angles, in degrees, that the points labeled 
; in the zenith array, za, and the azimuth array, az, lie , along great 
; circles, from the magnetic zenith vector.

del=ACOS(SIN(lambda)*SIN(beta)+COS(phi-alpha)*COS(lambda)*COS(beta))*!RADEG
mask=BYTARR(edge,edge)
; produce array to locate pixels in an image that lie DEG degrees from 
; the magnetic zenith.
mask(WHERE(del LE deg))=1B

cpix=MIN(WHERE(del EQ MIN(del)))
crp=cpix/edge
ccp=cpix MOD edge

RETURN
END