1 | MODULE function_ana |
---|
2 | ! |
---|
3 | IMPLICIT NONE |
---|
4 | ! |
---|
5 | CONTAINS |
---|
6 | ! |
---|
7 | !**************************************************************************************** |
---|
8 | SUBROUTINE function_ana1(ni, nj, xcoor, ycoor, fnc_ana) |
---|
9 | !**************************************************************************************** |
---|
10 | |
---|
11 | !**** *function ana* - calculate analytical function |
---|
12 | |
---|
13 | INTEGER, PARAMETER :: wp = SELECTED_REAL_KIND(12,307) ! double |
---|
14 | |
---|
15 | INTEGER, INTENT(IN) :: ni, nj |
---|
16 | REAL(kind=wp), DIMENSION(ni,nj), INTENT(IN) :: xcoor, ycoor |
---|
17 | REAL(kind=wp), DIMENSION(ni,nj), INTENT(OUT) :: fnc_ana |
---|
18 | |
---|
19 | REAL (kind=wp), PARAMETER :: dp_pi=3.14159265359 |
---|
20 | REAL (kind=wp), PARAMETER :: dp_conv = dp_pi/180. |
---|
21 | REAL(kind=wp) :: dp_length, coef, coefmult |
---|
22 | INTEGER :: i,j |
---|
23 | CHARACTER(LEN=7) :: cl_anaftype="fcos" |
---|
24 | |
---|
25 | DO j=1,nj |
---|
26 | DO i=1,ni |
---|
27 | |
---|
28 | SELECT CASE (cl_anaftype) |
---|
29 | CASE ("fcos") |
---|
30 | dp_length = 1.2*dp_pi |
---|
31 | coef = 2. |
---|
32 | coefmult = 1. |
---|
33 | fnc_ana(i,j) = coefmult*(coef - COS( dp_pi*(ACOS( COS(xcoor(i,j)*dp_conv)*COS(ycoor(i,j)*dp_conv) )/dp_length)) ) |
---|
34 | |
---|
35 | CASE ("fcossin") |
---|
36 | dp_length = 1.d0*dp_pi |
---|
37 | coef = 21.d0 |
---|
38 | coefmult = 3.846d0 * 20.d0 |
---|
39 | fnc_ana(i,j) = coefmult*(coef - COS( dp_pi*(ACOS( COS(ycoor(i,j)*dp_conv)*COS(ycoor(i,j)*dp_conv) )/dp_length)) * & |
---|
40 | SIN( dp_pi*(ASIN( SIN(xcoor(i,j)*dp_conv)*SIN(ycoor(i,j)*dp_conv) )/dp_length)) ) |
---|
41 | END SELECT |
---|
42 | |
---|
43 | ENDDO |
---|
44 | ENDDO |
---|
45 | |
---|
46 | END SUBROUTINE function_ana1 |
---|
47 | |
---|
48 | !!**************************************************************************************** |
---|
49 | !SUBROUTINE function_ana2(lon, lat, mask, fnc_ana) |
---|
50 | SUBROUTINE function_ana2(ni, nj, lon, lat, fnc_ana) |
---|
51 | !!**************************************************************************************** |
---|
52 | |
---|
53 | !**** *function ana* - calculate analytical function |
---|
54 | |
---|
55 | INTEGER, PARAMETER :: wp = SELECTED_REAL_KIND(12,307) ! double |
---|
56 | |
---|
57 | INTEGER, INTENT(IN) :: ni, nj |
---|
58 | REAL (kind=wp), DIMENSION(ni,nj), INTENT(IN) :: lon, lat |
---|
59 | REAL(kind=wp), DIMENSION(ni,nj), INTENT(OUT) :: fnc_ana |
---|
60 | |
---|
61 | REAL (kind=wp), PARAMETER :: coef=2., dp_pi=3.14159265359 |
---|
62 | REAL (kind=wp) :: dp_length, dp_conv |
---|
63 | INTEGER :: i, j |
---|
64 | |
---|
65 | ! Analytical Gulf Stream |
---|
66 | REAL (kind=wp) :: gf_coef, gf_ori_lon, gf_ori_lat, & |
---|
67 | & gf_end_lon, gf_end_lat, gf_dmp_lon, gf_dmp_lat |
---|
68 | REAL (kind=wp) :: gf_per_lon |
---|
69 | REAL (kind=wp) :: dx, dy, dr, dth, dc, dr0, dr1 |
---|
70 | |
---|
71 | ! Parameters for analytical function |
---|
72 | dp_length = 1.2*dp_pi |
---|
73 | dp_conv = dp_pi/180. |
---|
74 | gf_coef = 1.0 ! Coefficient for Gulf Stream term (0.0 = no Gulf Stream) |
---|
75 | gf_ori_lon = -80.0 ! Origin of the Gulf Stream (longitude in deg) |
---|
76 | gf_ori_lat = 25.0 ! Origin of the Gulf Stream (latitude in deg) |
---|
77 | gf_end_lon = -1.8 ! End point of the Gulf Stream (longitude in deg) |
---|
78 | gf_end_lat = 50.0 ! End point of the Gulf Stream (latitude in deg) |
---|
79 | gf_dmp_lon = -25.5 ! Point of the Gulf Stream decrease (longitude in deg) |
---|
80 | gf_dmp_lat = 55.5 ! Point of the Gulf Stream decrease (latitude in deg) |
---|
81 | |
---|
82 | ! ALLOCATE(fnc_ana(0:xy-1)) |
---|
83 | |
---|
84 | dr0 = SQRT(((gf_end_lon - gf_ori_lon)*dp_conv)**2 + & |
---|
85 | & ((gf_end_lat - gf_ori_lat)*dp_conv)**2) |
---|
86 | dr1 = SQRT(((gf_dmp_lon - gf_ori_lon)*dp_conv)**2 + & |
---|
87 | & ((gf_dmp_lat - gf_ori_lat)*dp_conv)**2) |
---|
88 | |
---|
89 | DO j=1,nj |
---|
90 | DO i=1,ni |
---|
91 | |
---|
92 | ! Original OASIS fcos analytical test function |
---|
93 | fnc_ana(i,j)=(coef-COS(dp_pi*(ACOS(COS(lat(i,j)*dp_conv)*& |
---|
94 | & COS(lon(i,j)*dp_conv))/dp_length))) |
---|
95 | gf_per_lon = lon(i,j) |
---|
96 | IF (gf_per_lon > 180.0) gf_per_lon = gf_per_lon - 360.0 |
---|
97 | IF (gf_per_lon < -180.0) gf_per_lon = gf_per_lon + 360.0 |
---|
98 | dx = (gf_per_lon - gf_ori_lon)*dp_conv |
---|
99 | dy = (lat(i,j) - gf_ori_lat)*dp_conv |
---|
100 | dr = SQRT(dx*dx + dy*dy) |
---|
101 | dth = ATAN2(dy, dx) |
---|
102 | dc = 1.3*gf_coef |
---|
103 | IF (dr > dr0) dc = 0.0 |
---|
104 | IF (dr > dr1) dc = dc * COS(dp_pi*0.5*(dr-dr1)/(dr0-dr1)) |
---|
105 | fnc_ana(i,j) = fnc_ana(i,j) + & |
---|
106 | & (MAX(1000.0*SIN(0.4*(0.5*dr+dth) + & |
---|
107 | & 0.007*COS(50.0*dth) + 0.37*dp_pi),999.0) - 999.0) * dc |
---|
108 | |
---|
109 | ENDDO |
---|
110 | ENDDO |
---|
111 | |
---|
112 | END SUBROUTINE function_ana2 |
---|
113 | |
---|
114 | !**************************************************************************************** |
---|
115 | SUBROUTINE function_vortex(ni, nj, xcoor, ycoor, fnc_ana) |
---|
116 | !**************************************************************************************** |
---|
117 | !!! VORTEX FROM TEMPEST-REMAP (as in XIOS) |
---|
118 | ! |
---|
119 | INTEGER, PARAMETER :: wp = SELECTED_REAL_KIND(12,307) ! double |
---|
120 | ! |
---|
121 | INTEGER, INTENT(IN) :: ni, nj |
---|
122 | REAL(kind=wp), DIMENSION(ni,nj), INTENT(IN) :: xcoor, ycoor |
---|
123 | REAL(kind=wp), DIMENSION(ni,nj), INTENT(OUT) :: fnc_ana |
---|
124 | |
---|
125 | REAL(kind=wp), PARAMETER :: dp_pi = 3.14159265359 |
---|
126 | REAL(kind=wp), PARAMETER :: dLon0 = 5.5 |
---|
127 | REAL(kind=wp), PARAMETER :: dLat0 = 0.2 |
---|
128 | REAL(kind=wp), PARAMETER :: dR0 = 3.0 |
---|
129 | REAL(kind=wp), PARAMETER :: dD = 5.0 |
---|
130 | REAL(kind=wp), PARAMETER :: dT = 6.0 |
---|
131 | REAL(kind=wp) :: dp_length, dp_conv |
---|
132 | |
---|
133 | REAL(kind=wp) :: dSinC, dCosC, dCosT, dSinT |
---|
134 | REAL(kind=wp) :: dTrm, dX, dY, dZ |
---|
135 | REAL(kind=wp) :: dlon, dlat |
---|
136 | REAL(kind=wp) :: dRho, dVt, dOmega |
---|
137 | |
---|
138 | INTEGER :: i,j |
---|
139 | CHARACTER(LEN=7) :: cl_anaftype="vortex" |
---|
140 | |
---|
141 | dp_conv = dp_pi/180. |
---|
142 | dSinC = SIN( dLat0 ) |
---|
143 | dCosC = COS( dLat0 ) |
---|
144 | |
---|
145 | DO j=1,nj |
---|
146 | DO i=1,ni |
---|
147 | ! Find the rotated longitude and latitude of a point on a sphere |
---|
148 | ! with pole at (dLon0, dLat0). |
---|
149 | dCosT = COS( ycoor(i,j)*dp_conv ) |
---|
150 | dSinT = SIN( ycoor(i,j)*dp_conv ) |
---|
151 | |
---|
152 | dTrm = dCosT * COS( xcoor(i,j)*dp_conv - dLon0 ) |
---|
153 | dX = dSinC * dTrm - dCosC * dSinT |
---|
154 | dY = dCosT * SIN( xcoor(i,j)*dp_conv - dLon0 ) |
---|
155 | dZ = dSinC * dSinT + dCosC * dTrm |
---|
156 | |
---|
157 | dlon = ATAN2( dY, dX ) |
---|
158 | IF( dlon < 0.0 ) dlon = dlon + 2.0 * dp_pi |
---|
159 | dlat = ASIN( dZ ) |
---|
160 | |
---|
161 | dRho = dR0 * COS(dlat) |
---|
162 | dVt = 3.0 * SQRT(3.0)/2.0/COSH(dRho)/COSH(dRho)*TANH(dRHo) |
---|
163 | IF (dRho == 0.0) THEN |
---|
164 | dOmega = 0.0 |
---|
165 | ELSE |
---|
166 | dOmega = dVt / dRho |
---|
167 | END IF |
---|
168 | |
---|
169 | fnc_ana(i,j) = 2.0 * ( 1.0 + TANH( dRho / dD * SIN( dLon - dOmega * dT ) ) ) |
---|
170 | |
---|
171 | END DO |
---|
172 | END DO |
---|
173 | ! |
---|
174 | END SUBROUTINE function_vortex |
---|
175 | |
---|
176 | END MODULE function_ana |
---|