PRO ieasegrid,iopt,alon,alat,thelon,thelat,ascale
;
;	COMPUTES THE INVERSE "EASE" GRID TRANSFORM
;
;	GIVEN THE IMAGE TRANSFORMATION COORDINATES (THELON,THELAT) AND
;	THE SCALE (ASCALE) THE CORRESPONDING LON,LAT (ALON,ALAT) IS COMPUTED
;	USING THE "EASE GRID" (VERSION 1.0) TRANSFORMATION GIVEN IN FORTRAN
;	SOURCE CODE SUPPLIED BY NSIDC
;
;	IOPT IS EASE TYPE: IOPT=11=NORTH, IOPT=12=SOUTH, IOPT=13=CYLINDRICAL
;
;	THE RADIUS OF THE EARTH USED IN THIS PROJECTION IS IMBEDDED INTO
;	ASCALE WHILE THE PIXEL DIMENSION IN KM IS IMBEDDED IN BSCALE
;	THE BASE VALUES ARE: RADIUS EARTH= 6371.228 KM
;			     PIXEL DIMEN =25.067525 KM
;	THEN, BSCALE = BASE_PIXEL_DIMEN
;	      ASCALE = RADIUS_EARTH/BASE_PIXEL_DIMEN
;

; revised and corrected by dgl 23 july 2005
	PI2=1.57079633D0		; Pi/2
	DTR=PI2/90.0D0
	X1=THELON
	Y1=THELAT
	IF IOPT EQ 11 THEN BEGIN	; EASE GRID NORTH
		ALON=ATAN(X1,-Y1)/DTR
		IF ABS(SIN(DTR*ALON)) GT ABS(COS(ALON*DTR)) THEN $
		  TEMP=(X1/SIN(ALON*DTR))/ASCALE $
		ELSE
		  TEMP=(-Y1/COS(ALON*DTR))/ASCALE
		IF ABS(TEMP) LE 1.0 THEN $
			ALAT=90.0D0-2.*ASIN(TEMP)/DTR $
		ELSE BEGIN
			ALAT=-90.0D0+2.*ASIN(1./TEMP)/DTR
			PRINT,'*** ERROR IN EASE INVERSE SINE ***',TEMP
		ENDELSE
	ENDIF ELSE IF IOPT EQ 12 THEN BEGIN	; EASE GRID SOUTH
		ALON=ATAN(X1,Y1)/DTR
		IF ABS(COS(ALON*DTR)) GT ABS(SIN(ALON*DTR)) THEN $
		    TEMP=(Y1/COS(ALON*DTR))/ASCALE $
		ELSE $
		    TEMP=(X1/SIN(ALON*DTR))/ASCALE
		IF ABS(TEMP) LE 1.0 THEN $
			ALAT=90.0D0-2.*ACOS(TEMP)/DTR $
		ELSE BEGIN
			PRINT,'*** ERROR IN EASE INVERSE COSINE ***',TEMP
			ALAT=0.0
		ENDELSE
	ENDIF ELSE IF IOPT EQ 13 THEN BEGIN	; EASE CYLINDRICAL
		ALON=((X1/ASCALE)/COS(30.0*DTR))*90.0D0/PI2
		TEMP=(Y1*COS(30.0*DTR))/ASCALE
		IF ABS(TEMP) LE 1.0 THEN $
			ALAT=ASIN(TEMP)/DTR $
		ELSE BEGIN
			PRINT,'*** ERROR IN EASE INVERSE SINE ***',TEMP
			ALAT=SIGN(90.0,TEMP)
		ENDELSE
	ENDIF
	END