da_set_lc.inc

1 subroutine da_set_lc(proj)
2 
3    !-----------------------------------------------------------------------
4    ! Purpose: Initialize the remaining items in the proj structure for a
5    ! lambert conformal grid.
6    !-----------------------------------------------------------------------
7 
8    implicit none
9     
10    type(proj_info), intent(inout)     :: proj
11 
12    real                               :: arg
13    real                               :: deltalon1
14    real                               :: tl1r
15    real                               :: ctl1r
16 
17    ! Compute cone factor
18    call da_lc_cone(proj%truelat1, proj%truelat2, proj%cone)
19    if (print_detail_map) then
20       print '(A,F8.6)', 'Computed cone factor: ', proj%cone
21    end if
22    ! Compute longitude differences and ensure we stay out of the
23    ! forbidden "cut zone"
24    deltalon1 = proj%lon1 - proj%stdlon
25    if (deltalon1 > +180.) deltalon1 = deltalon1 - 360.
26    if (deltalon1 < -180.) deltalon1 = deltalon1 + 360.
27 
28    ! Convert truelat1 to radian and compute COS for later use
29    tl1r = proj%truelat1 * rad_per_deg
30    ctl1r = COS(tl1r)
31 
32    ! Compute the radius to our known lower-left (SW) corner
33    proj%rsw = proj%rebydx * ctl1r/proj%cone * &
34            (TAN((90.*proj%hemi-proj%lat1)*rad_per_deg/2.) / &
35             TAN((90.*proj%hemi-proj%truelat1)*rad_per_deg/2.))**proj%cone
36 
37    ! Find pole point
38    arg = proj%cone*(deltalon1*rad_per_deg)
39    proj%polei = 1. - proj%hemi * proj%rsw * Sin(arg)
40    proj%polej = 1. + proj%rsw * COS(arg)  
41    if (print_detail_map) then
42       print '(A,2F10.3)', 'Computed pole i/j = ', proj%polei, proj%polej
43    end if
44 
45 end subroutine da_set_lc                             
46 
47