Category : Miscellaneous Language Source Code
Archive   : LAT-LON.ZIP
Filename : MAPPRO.FOR

C M A P P R O
C
C MAP PROjection returns the projected cartesian co-ordinate
C positions given lat and long co-ordinates.
C
C CALL MAPPRO (PROTYP,RPX,RPY,RPLON,RPLAT,RPSCAL,RLCNT,RLBND,RMGSZ)
C
C PROTYP - (Int) Projection type [I]
C 1 - Equirectangular
C 2 - Stereographic
C 3 - Polar stereographic
C RPX - (R*4) Projected X co-ord [O]
C RPY - (R*4) Projected Y co-ord [O]
C RPLON - (R*4) Longitude in real degrees [I]
C RPLAT - (R*4) Latitude in real degrees [I]
C RPSCAL - (R*4) Projection scaling parameters [I]
C (1) - X multiplicative factor
C (2) - X additive factor
C (3) - Y multiplicative factor
C (4) - Y additive factor
C RLCNT - (R*4) Long/Lat at centre of input window [I]
C (1) - Lon centre in radian measure
C (2) - Lat centre in radian measure
C RLBND - (R*4) Input Lon/Lat window bounds in [I]
C radian measure -
C (1) - Lower longitude of area
C (2) - Upper longitude of area
C (3) - Lower latitude of area
C (4) - Upper latitude of area
C RMGSZ - (R*4) Output map-grid dimensions, [I]
C (pixels,lines) - a 2-element array.
C NOTE:- Parameter RPSCAL scales the
C output coordinates RPX & RPY to fit
C into the map-grid dimensions.
C For proper use of this module: if RPSCAL
C has been pre-determined using module MAPSCL,
C then RMGSZ(1) = XYBND(3),
C RMGSZ(2) = XYBND(4).
C Refer to module MAPSCL for a definition
C of parameter XYBND.
C-

SUBROUTINE MAPPRO (PROTYP,RPX,RPY,RPLON,RPLAT,
& RPSCAL,RLCNT,RLBND,RMGSZ)
C ---------------------------------------------------------
C - Parameters -
C ---------------------------------------------------------
INTEGER PROTYP
REAL RPX,RPY,RPLON,RPLAT,RPSCAL(4),RLCNT(2)
REAL RLBND(4),RMGSZ(2)
C ---------------------------------------------------------
C - Local Definitions -
C ---------------------------------------------------------
REAL CDLON,COSLAT,COSLT0,SINLAT,SINLT0,RDEN
REAL RPNT(2),RLOCLT,RLOCLN,RTAN
REAL RU,RV,R90,R270

C ---------------------------------------------------------
C - Functions -
C ---------------------------------------------------------
REAL RADEG
C ---------------------------------------------------------
C - Data -
C ---------------------------------------------------------
DATA R90,R270 /1.5707963,4.712389/
C ---------------------------------------------------------
RLOCLT = RADEG (RPLAT)
RLOCLN = RADEG (RPLON)
IF (PROTYP.EQ.1) GO TO 100
IF (PROTYP.EQ.2) GO TO 200
IF (PROTYP.EQ.3) GO TO 300
RETURN
100 IF (RLBND(1).LT.0) RPNT(1) = RLOCLN - RLBND(1)
IF (RLBND(1).GE.0) RPNT(1) = RLBND(1) - RLOCLN
IF (RLBND(4).LT.0) RPNT(2) = RLOCLT - RLBND(4)
IF (RLBND(4).GE.0) RPNT(2) = RLBND(4) - RLOCLT
RPX = RPNT(1)*RPSCAL(1) + RPSCAL(2)
RPY = RPNT(2)*RPSCAL(3) + RPSCAL(4)
RETURN
200 CDLON = COS(RLOCLN-RLCNT(1))
COSLAT = COS(RLOCLT)
COSLT0 = COS(RLCNT(2))
SINLAT = SIN(RLOCLT)
SINLT0 = SIN(RLCNT(2))
RDEN = 1.0 + SINLAT*SINLT0 + COSLAT*COSLT0*CDLON
RPX = (COSLAT*SIN(RLOCLN-RLCNT(1)))/RDEN
RPY = (SINLAT*COSLT0 - SINLT0*COSLAT*CDLON)/RDEN
RPX = RPX*RPSCAL(1) + RPSCAL(2)
RPY = RPY*RPSCAL(3) + RPSCAL(4)
RETURN
300 RTAN = 2*TAN((R90 - RLOCLT)/2)
RU = RPSCAL(1) * (RTAN * COS(RLOCLN))
RV = RPSCAL(3) * (RTAN * SIN(RLOCLN))
RPX = RU*COS(R270-RLCNT(1)) - RV*SIN(R270-RLCNT(1)) + RPSCAL(2)
RPY = RU*SIN(R270-RLCNT(1)) + RV*COS(R270-RLCNT(1)) + RPSCAL(4)
RPY = RMGSZ(2)-RPY
RETURN
END