forked from NASA-LIS/LISF
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathpolfixs.F90
110 lines (108 loc) · 2.57 KB
/
polfixs.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
!-----------------------BEGIN NOTICE -- DO NOT EDIT-----------------------
! NASA Goddard Space Flight Center
! Land Information System Framework (LISF)
! Version 7.4
!
! Copyright (c) 2022 United States Government as represented by the
! Administrator of the National Aeronautics and Space Administration.
! All Rights Reserved.
!-------------------------END NOTICE -- DO NOT EDIT-----------------------
!BOP
!
! !ROUTINE: polfixs
! \label{polfixs}
!
! !REVISION HISTORY:
! 04-10-96 Mark Iredell; Initial Specification
!
! !INTERFACE:
subroutine polfixs(nm,nx,km,rlat,rlon,ib,lo,go)
implicit none
! !ARGUMENTS:
integer :: nm
integer :: nx
integer :: km
real :: rlat(nm)
real :: rlon(nm)
integer :: ib
logical*1 :: lo(nx)
real :: go(nx)
!
! !DESCRIPTION:
! This subroutine averages multiple pole scalar values
! on a latitude/longitude grid. bitmaps may be averaged too.
!
! The arguments are:
! \begin{description}
! \item[nm]
! number of grid points
! \item[nx]
! leading dimension of fields
! \item[rlat]
! latitudes in degrees
! \item[rlon]
! longitudes in degrees
! \item[ib]
! integer bitmap flags
! \item[lo]
! logical bitmaps
! \item[go]
! returned scalar value
! \end{description}
!
!EOP
integer :: n, k
real :: tsp, gnp, gsp, wsp, tnp, wnp
real, PARAMETER :: rlatnp=89.9995, rlatsp=-89.9995
! do k=1,km
wnp=0.0
gnp=0.0
tnp=0.0
wsp=0.0
gsp=0.0
tsp=0.0
! average multiple pole values
do n=1,nm
if(rlat(n).ge.rlatnp) then
wnp=wnp+1
if(ib.eq.0.or.lo(n)) then
gnp=gnp+go(n)
tnp=tnp+1
endif
elseif(rlat(n).le.rlatsp) then
wsp=wsp+1
if(ib.eq.0.or.lo(n)) then
gsp=gsp+go(n)
tsp=tsp+1
endif
endif
enddo
! distribute average values back to multiple poles
if(wnp.gt.1) then
if(tnp.ge.wnp/2) then
gnp=gnp/tnp
else
gnp=0.
endif
do n=1,nm
if(rlat(n).ge.rlatnp) then
if(ib.ne.0) lo(n)=tnp.ge.wnp/2
go(n)=gnp
endif
enddo
endif
if(wsp.gt.1) then
if(tsp.ge.wsp/2) then
gsp=gsp/tsp
else
gsp=0.
endif
do n=1,nm
if(rlat(n).le.rlatsp) then
if(ib.ne.0) lo(n)=tsp.ge.wsp/2
go(n)=gsp
endif
enddo
endif
! enddo
end subroutine polfixs