1module verticalModes_mod
2 ! MODULE verticalModes_mod (prefix='vms' category='4. Data Object transformations')
3 !
4 !:Purpose: To 1) compute empirical orthogonal functions (EOFs) from either ensemble-derived
5 ! vertical background-error covariances matrices or a prescribed vertical correlation
6 ! function (i.e., the so-called vertical modes) and to 2) project back or forth
7 ! ensemble pertubations onto these modes.
8 ! Therefore, capablity #2 behaves like a spectral transform but in the vertical
9 ! dimension.
10 !
11 use earthConstants_mod
12 use utilities_mod
13 use varNameList_mod
14 use horizontalCoord_mod
15 use verticalCoord_mod
16 use localizationFunction_mod
17 use calcHeightAndPressure_mod
18 use gridStatevector_mod
19 use ensembleStatevector_mod
20 use timeCoord_mod
21 use midasMpi_mod
22 implicit none
23 save
24 private
25
26 ! Public derived type
27 public :: struct_vms
28 ! Public subroutines
29 public :: vms_computeModesFromEns, vms_computeModesFromFunction
30 public :: vms_writeModes, vms_transform
31
32 type :: struct_oneVar3d
33 integer :: nLev
34 character(len=4) :: varName
35 real(8), allocatable :: autoCovariance(:,:)
36 real(8), allocatable :: eigenVectors(:,:)
37 real(8), allocatable :: eigenValues(:)
38 real(8), allocatable :: eigenVectorsInv(:,:)
39 end type struct_oneVar3d
40
41 type :: struct_vms
42 private
43 integer :: nVar3d
44 type(struct_oneVar3d), allocatable :: allVar3d(:)
45 logical :: initialized = .false.
46 end type struct_vms
47
48contains
49
50 !--------------------------------------------------------------------------
51 ! vms_setup
52 !--------------------------------------------------------------------------
53 subroutine vms_setup(varNamesList, vco, vModes)
54 !
55 !:Purpose: To setup the vModes structure
56 !
57 implicit none
58
59 ! Arguments:
60 character(len=4), intent(in) :: varNamesList(:)
61 type(struct_vco), pointer, intent(in) :: vco
62 type(struct_vms), intent(inout) :: vModes
63
64 ! Locals:
65 integer :: nLev, varNameIndex, var3dIndex
66
67 if (vModes%initialized) return
68
69 !
70 !- Set parameters
71 !
72 vModes%nVar3d = 0
73 do varNameIndex = 1, size(varNamesList)
74 if (vnl_varLevelFromVarname(varNamesList(varNameIndex)) /= 'SF' ) then
75 vModes%nVar3d = vModes%nVar3d + 1
76 end if
77 end do
78
79 allocate(vModes%allVar3d(vModes%nVar3d))
80
81 var3dIndex = 0
82 do varNameIndex = 1, size(varNamesList)
83 if (vnl_varLevelFromVarname(varNamesList(varNameIndex)) /= 'SF' ) then
84 var3dIndex = var3dIndex + 1
85 vModes%allVar3d(var3dIndex)%varName = varNamesList(varNameIndex)
86 if (vnl_varLevelFromVarName(vModes%allVar3d(var3dIndex)%varName) == 'MM') then
87 nLev = vco_getNumLev(vco,'MM')
88 else
89 nLev = vco_getNumLev(vco,'TH')
90 end if
91 allocate(vModes%allVar3d(var3dIndex)%autoCovariance(nLev,nLev))
92 allocate(vModes%allVar3d(var3dIndex)%eigenVectors(nLev,nLev))
93 allocate(vModes%allVar3d(var3dIndex)%eigenValues(nLev))
94 allocate(vModes%allVar3d(var3dIndex)%eigenVectorsInv(nLev,nLev))
95 vModes%allVar3d(var3dIndex)%nLev = nLev
96 end if
97 end do
98
99 !
100 !- Ending
101 !
102 vModes%initialized = .true.
103
104 end subroutine vms_setup
105
106 !--------------------------------------------------------------------------
107 ! vms_computeModesFromEns
108 !--------------------------------------------------------------------------
109 subroutine vms_computeModesFromEns(ensPerts,vModes)
110 !
111 !:Purpose: To compute vertical modes from ensemble-derived
112 ! vertical background-error covariances matrices
113 !
114 implicit none
115
116 ! Arguments:
117 type(struct_ens), intent(inout) :: ensPerts
118 type(struct_vms), intent(inout) :: vModes
119
120 ! Locals:
121 type(struct_hco), pointer :: hco_ens
122 type(struct_vco), pointer :: vco_ens
123 real(4), pointer :: ptr4d_1_r4(:,:,:,:)
124 real(4), pointer :: ptr4d_2_r4(:,:,:,:)
125 real(8), allocatable :: latWeight(:) ! Weight given to grid point in the statistic computation
126 real(8), allocatable :: sumWeight(:,:)
127 integer :: varLevIndex1, varLevIndex2, levIndex1, levIndex2
128 integer :: lonIndex, latIndex, nEns, memberIndex
129 integer :: var3dIndex, myLonBeg, myLonEnd, myLatBeg, myLatEnd
130 character(len=4), pointer :: varNamesList(:)
131
132 !
133 !- Structure initialization
134 !
135 nullify(varNamesList)
136 call ens_varNamesList(varNamesList, ensPerts)
137
138 vco_ens => ens_getVco(ensPerts)
139
140 call vms_setup(varNamesList, vco_ens, & ! IN
141 vModes) ! OUT
142
143 !
144 !- Compute the vertical auto-covariance matrices
145 !
146 nEns = ens_getNumMembers(ensPerts)
147
148 hco_ens => ens_getHco(ensPerts)
149 allocate(latWeight(hco_ens%nj))
150 do latIndex = 1, hco_ens%nj
151 latWeight(latIndex) = cos(hco_ens%lat(latIndex))
152 if (mmpi_myid == 0) then
153 write(*,*) latIndex, hco_ens%lat(latIndex), cos(hco_ens%lat(latIndex))
154 end if
155 end do
156
157 call ens_getLatLonBounds(ensPerts, myLonBeg, myLonEnd, myLatBeg, myLatEnd)
158
159 do var3dIndex = 1, vModes%nVar3d
160
161 allocate(sumWeight(vModes%allVar3d(var3dIndex)%nLev,vModes%allVar3d(var3dIndex)%nLev))
162 sumWeight(:,:) = 0.d0
163
164 vModes%allVar3d(var3dIndex)%autoCovariance(:,:) = 0.d0
165
166 do memberIndex = 1, nEns
167
168 !$OMP PARALLEL DO PRIVATE (varLevIndex1,varLevIndex2,levIndex2,levIndex1,ptr4d_1_r4,ptr4d_2_r4,latIndex,lonIndex)
169 do levIndex2 = 1, vModes%allVar3d(var3dIndex)%nLev
170 varLevIndex2 = levIndex2 + ens_getOffsetFromVarName(ensPerts,vModes%allVar3d(var3dIndex)%varName)
171 do levIndex1 = 1, vModes%allVar3d(var3dIndex)%nLev
172 varLevIndex1 = levIndex1 + ens_getOffsetFromVarName(ensPerts,vModes%allVar3d(var3dIndex)%varName)
173
174 ptr4d_1_r4 => ens_getOneLev_r4(ensPerts,varLevIndex1)
175 ptr4d_2_r4 => ens_getOneLev_r4(ensPerts,varLevIndex2)
176
177 do latIndex = myLatBeg, myLatEnd
178 do lonIndex = myLonBeg, myLonEnd
179
180 vModes%allVar3d(var3dIndex)%autoCovariance(levIndex1,levIndex2) = &
181 vModes%allVar3d(var3dIndex)%autoCovariance(levIndex1,levIndex2) &
182 + real(ptr4d_1_r4(memberIndex,1,lonIndex,latIndex),8) &
183 * real(ptr4d_2_r4(memberIndex,1,lonIndex,latIndex),8) &
184 * latWeight(latIndex)
185
186 sumWeight(levIndex1,levIndex2) = sumWeight(levIndex1,levIndex2) + latWeight(latIndex)
187
188 end do
189 end do
190
191 end do
192 end do
193 !$OMP END PARALLEL DO
194
195 end do ! Loop on Ensemble
196
197 !- Communication
198 call mmpi_allreduce_sumR8_2d(vModes%allVar3d(var3dIndex)%autoCovariance, "GRID")
199 call mmpi_allreduce_sumR8_2d(sumWeight, "GRID")
200
201 !- Apply the appropriate scaling
202 do levIndex2 = 1, vModes%allVar3d(var3dIndex)%nLev
203 do levIndex1 = 1, vModes%allVar3d(var3dIndex)%nLev
204 vModes%allVar3d(var3dIndex)%autoCovariance(levIndex1,levIndex2) = &
205 vModes%allVar3d(var3dIndex)%autoCovariance(levIndex1,levIndex2) &
206 / (sumWeight(levIndex1,levIndex2) / dble(nEns) * dble(nEns-1))
207 end do
208 end do
209
210 deallocate(sumWeight)
211
212 end do ! Loop on 3d variables
213
214 !
215 !- Compute the eigenvectors and eigenvalue
216 !
217 call vms_computeModes(vModes)
218
219 end subroutine vms_computeModesFromEns
220
221 !--------------------------------------------------------------------------
222 ! vms_computeModesFromFunction
223 !--------------------------------------------------------------------------
224 subroutine vms_computeModesFromFunction(vco, lengthScaleTop, lengthScaleBot, &
225 vModes)
226 !
227 !:Purpose: To compute vertical modes from a prescribed correlation function
228 !
229 implicit none
230
231 ! Arguments:
232 type(struct_vco), pointer, intent(in) :: vco
233 real(8) , intent(in) :: lengthScaleTop
234 real(8) , intent(in) :: lengthScaleBot
235 type(struct_vms), intent(inout) :: vModes
236
237 ! Locals:
238 real(8), pointer :: vertLocation_MM(:)
239 real(8), pointer :: vertLocation_TH(:)
240 real(8), pointer :: vertLocation(:)
241 real(8) :: pSurfRef, distance, correl
242 real(8) :: lengthScaleGradient
243 real(8) :: lengthScaleLev1, lengthScaleLev2, lengthScaleAvg
244 integer :: levIndex, levIndex1, levIndex2
245 integer :: var3dIndex
246 character(len=4) :: varNamesList(2)
247
248 !
249 !- Structure initialization
250 !
251 varNamesList = (/'P_M','P_T'/)
252 call vms_setup(varNamesList, vco, & ! IN
253 vModes) ! OUT
254
255 !
256 !- Compute the vertical auto-covariance matrices
257 !
258 pSurfRef = 101000.D0
259 call czp_fetch1DLevels(vco, pSurfRef, & ! IN
260 profM_opt=vertLocation_MM) ! OUT
261 call czp_fetch1DLevels(vco, pSurfRef, & ! IN
262 profT_opt=vertLocation_TH) ! OUT
263
264 do levIndex = 1, vco%nLev_M
265 vertLocation_MM(levIndex) = log(vertLocation_MM(levIndex))
266 end do
267 do levIndex = 1, vco%nLev_T
268 vertLocation_TH(levIndex) = log(vertLocation_TH(levIndex))
269 end do
270
271 call lfn_setup('FifthOrder') ! IN
272
273 do var3dIndex = 1, vModes%nVar3d
274
275 if (vnl_varLevelFromVarName(vModes%allVar3d(var3dIndex)%varName) == 'MM') then
276 vertLocation => vertLocation_MM
277 else
278 vertLocation => vertLocation_TH
279 end if
280
281 lengthScaleGradient = (lengthScaleTop-lengthScaleBot) / &
282 (vertLocation(1)-vertLocation(vModes%allVar3d(var3dIndex)%nLev))
283
284 do levIndex1 = 1, vModes%allVar3d(var3dIndex)%nLev
285 do levIndex2 = 1, vModes%allVar3d(var3dIndex)%nLev
286 distance = abs(vertLocation(levIndex2) - vertLocation(levIndex1))
287
288 lengthScaleLev1 = lengthScaleGradient * &
289 (vertLocation(levIndex1)-vertLocation(vModes%allVar3d(var3dIndex)%nLev)) + &
290 lengthScaleBot
291 lengthScaleLev2 = lengthScaleGradient * &
292 (vertLocation(levIndex2)-vertLocation(vModes%allVar3d(var3dIndex)%nLev)) + &
293 lengthScaleBot
294 lengthScaleAvg = (lengthScaleLev1+lengthScaleLev2)/2.d0
295
296 if (levIndex1 == 1 .and. mmpi_myid == 0) then
297 write(*,*) 'vms_computeModesFromFunction: function length scale (', &
298 levIndex2,') = ', lengthScaleAvg
299 end if
300
301 correl = lfn_response(distance,lengthScaleAvg)
302 vModes%allVar3d(var3dIndex)%autoCovariance(levIndex1,levIndex2) = correl
303 end do
304 end do
305
306 end do
307
308 !
309 !- Compute the eigenvectors and eigenvalue
310 !
311 call vms_computeModes(vModes)
312
313 end subroutine vms_computeModesFromFunction
314
315 !--------------------------------------------------------------------------
316 ! vms_computeModes
317 !--------------------------------------------------------------------------
318 subroutine vms_computeModes(vModes)
319 !
320 !:Purpose: Compute vertical modes from the vertical covariances matrices
321 ! contained in the vModes structure.
322 !
323 implicit none
324
325 ! Arguments:
326 type(struct_vms), intent(inout) :: vModes
327
328 ! Locals:
329 real(8) :: tolerance
330 integer :: levIndex1, levIndex2, var3dIndex
331 integer :: matrixRank
332
333 !
334 !- Compute the eigenvectors and eigenvalue
335 !
336 tolerance = 1.0D-50
337 do var3dIndex = 1, vModes%nVar3d
338
339 !- Compute the eigenvectors and eigenvalue
340 call utl_eigenDecomp(vModes%allVar3d(var3dIndex)%autoCovariance, & ! IN
341 vModes%allVar3d(var3dIndex)%eigenValues, & ! OUT
342 vModes%allVar3d(var3dIndex)%eigenVectors, & ! OUT
343 tolerance, & ! IN
344 matrixRank) ! OUT
345 if ( matrixRank < vModes%allVar3d(var3dIndex)%nLev ) then
346 write(*,*) 'varName =', vModes%allVar3d(var3dIndex)%varName
347 write(*,*) 'matrixRank =', matrixRank
348 write(*,*) 'nLev =', vModes%allVar3d(var3dIndex)%nLev
349 call utl_abort('vms_computeModes: vertical matrix is rank deficient')
350 end if
351
352 !- Compute the inverse of the eigenVector matrix to go from grid point space
353 ! to vertical mode space
354 ! Note: The inverse of an orthogonal matrix is simply its transpose
355 do levIndex2 = 1, vModes%allVar3d(var3dIndex)%nLev
356 do levIndex1 = 1, vModes%allVar3d(var3dIndex)%nLev
357 vModes%allVar3d(var3dIndex)%eigenVectorsInv(levIndex2,levIndex1) = &
358 vModes%allVar3d(var3dIndex)%eigenVectors(levIndex1,levIndex2)
359 end do
360 end do
361
362 end do
363
364 end subroutine vms_computeModes
365
366 !--------------------------------------------------------------------------
367 ! vms_transform
368 !--------------------------------------------------------------------------
369 subroutine vms_transform(vModes, VertModesState, GridState, &
370 TransformDirection, lonBeg, lonEnd, &
371 latBeg, latEnd, nLev, varName)
372 !
373 !:Purpose: To project back or forth ensemble pertubations onto the vertical
374 ! modes contained in the vModes structure.
375 !
376 implicit none
377
378 ! Arguments:
379 type(struct_vms), intent(in) :: vModes
380 integer, intent(in) :: lonBeg
381 integer, intent(in) :: lonEnd
382 integer, intent(in) :: latBeg
383 integer, intent(in) :: latEnd
384 integer, intent(in) :: nLev
385 real(8), intent(inout) :: VertModesState(lonBeg:lonEnd,latBeg:latEnd,nLev) ! 3D vertical modes coefficients
386 real(8), intent(inout) :: GridState(lonBeg:lonEnd,latBeg:latEnd,nLev) ! 3D field in grid point space
387 character(len=*), intent(in) :: TransformDirection ! VertModesToGridPoint or GridPointToVertModes
388 character(len=*), intent(in) :: varName
389
390 ! Locals:
391 integer :: latIndex, lonIndex, nMode
392 integer :: varIndexAssociated, varIndex
393
394 !
395 !- 1. Find the appropriate modes for the provided varName
396 !
397 varIndexAssociated = -1
398
399 ! Search for a match in varName
400 do varIndex = 1, vModes%nVar3d
401 if (vModes%allVar3d(varIndex)%varName == varName) then
402 varIndexAssociated = varIndex
403 exit
404 end if
405 end do
406
407 ! Search for a match in level type (MM or TH)
408 if (varIndexAssociated == -1) then
409 if (vnl_varLevelFromVarName(varName) == 'MM') then
410 do varIndex = 1, vModes%nVar3d
411 if (vModes%allVar3d(varIndex)%varName == 'P_M') then
412 varIndexAssociated = varIndex
413 exit
414 end if
415 end do
416 else if (vnl_varLevelFromVarName(varName) == 'TH') then
417 do varIndex = 1, vModes%nVar3d
418 if (vModes%allVar3d(varIndex)%varName == 'P_T') then
419 varIndexAssociated = varIndex
420 exit
421 end if
422 end do
423 else
424 write(*,*) trim(vModes%allVar3d(varIndex)%varName), vnl_varLevelFromVarName(varName)
425 call utl_abort('vms_transform: varName is not on momentum or thermodynamic levels')
426 end if
427 end if
428
429 if (varIndexAssociated == -1) then
430 write(*,*)
431 write(*,*) 'varName in input', trim(varName)
432 write(*,*) 'varNames found in vModes'
433 do varIndex = 1, vModes%nVar3d
434 write(*,*) ' ... ', trim(vModes%allVar3d(varIndex)%varName)
435 end do
436 call utl_abort('vms_transform: could not match varName with any modes')
437 end if
438
439 if (nLev /= vModes%allVar3d(varIndexAssociated)%nLev) then
440 call utl_abort('vms_transform: the number of levels are not consistent')
441 end if
442 nMode = nLev
443
444 !
445 !- 2. Do the transform
446 !
447 select case (trim(TransformDirection))
448 case ('GridPointToVertModes')
449
450 !$OMP PARALLEL DO PRIVATE (latIndex,lonIndex)
451 do latIndex = latBeg, latEnd
452 do lonIndex = lonBeg, lonEnd
453 vertModesState(lonIndex,latIndex,1:nMode) = &
454 matmul(vModes%allVar3d(varIndexAssociated)%eigenVectorsInv(1:nMode,1:nLev), &
455 gridState(lonIndex,latIndex,1:nLev))
456 end do
457 end do
458 !$OMP END PARALLEL DO
459
460 case ('VertModesToGridPoint')
461
462 !$OMP PARALLEL DO PRIVATE (latIndex,lonIndex)
463 do latIndex = latBeg, latEnd
464 do lonIndex = lonBeg, lonEnd
465 gridState(lonIndex,latIndex,1:nLev) = &
466 matmul(vModes%allVar3d(varIndexAssociated)%eigenVectors(1:nLev,1:nMode), &
467 vertModesState(lonIndex,latIndex,1:nMode))
468 end do
469 end do
470 !$OMP END PARALLEL DO
471
472 case default
473 write(*,*)
474 write(*,*) 'Error: TranformDirection Unknown ', trim(TransformDirection)
475 call utl_abort('vms_transform')
476 end select
477
478 end subroutine vms_transform
479
480 !--------------------------------------------------------------------------
481 ! vms_writeModes
482 !--------------------------------------------------------------------------
483 subroutine vms_writeModes(vModes)
484 !
485 !:Purpose: To write the content of the provided vModes structure
486 !
487 implicit none
488
489 ! Arguments:
490 type(struct_vms), intent(in) :: vModes
491
492 ! Locals:
493 integer :: var3dIndex, levIndex1, levIndex2
494 integer :: fstouv, fnom, fstfrm, fclos, iunstats, errorID
495 character(len=4), allocatable :: varName3d(:)
496 character(len=128) :: outfilename
497
498 if (.not. vModes%initialized) then
499 call utl_abort('vms_writeModes: The vModes structure is not initialized')
500 end if
501
502 if (mmpi_myid == 0) then
503 !
504 !- Write in FST file format
505 !
506 allocate(varName3d(vModes%nVar3d))
507 iunstats = 0
508 errorID = fnom(iunstats,'./vertAutoCov.fst','RND',0)
509 errorID = fstouv(iunstats,'RND')
510 do var3dIndex = 1, vModes%nVar3d
511 call writeMatrix2d_r8(vModes%allVar3d(var3dIndex)%autoCovariance, vModes%allVar3d(var3dIndex)%nlev, &
512 iunstats, vModes%allVar3d(var3dIndex)%varName,'AUTOCOVAR')
513 call writeMatrix2d_r8(vModes%allVar3d(var3dIndex)%eigenVectors, vModes%allVar3d(var3dIndex)%nlev, &
514 iunstats, vModes%allVar3d(var3dIndex)%varName,'EIGENVECTORS')
515 call writeArray1d_r8(vModes%allVar3d(var3dIndex)%eigenValues, vModes%allVar3d(var3dIndex)%nlev, &
516 iunstats, vModes%allVar3d(var3dIndex)%varName,'EIGENVALUES')
517 varName3d(var3dIndex) = vModes%allVar3d(var3dIndex)%varName
518 end do
519 call writeArray1d_c4(varName3d, vModes%nVar3d, iunstats, 'VAR', 'VARNAME')
520 errorID = fstfrm(iunstats)
521 errorID = fclos (iunstats)
522
523 !
524 !- Write in txt format (for plotting purposes)
525 !
526 do var3dIndex = 1, vModes%nVar3d
527 outfilename = "./vModes_autoCovar_"//trim(vModes%allVar3d(var3dIndex)%varName)//".txt"
528 open (unit=99,file=outfilename,action="write",status="new")
529 do levIndex2 = 1, vModes%allVar3d(var3dIndex)%nLev
530 do levIndex1 = 1, vModes%allVar3d(var3dIndex)%nLev
531 if (levIndex1 == vModes%allVar3d(var3dIndex)%nLev) then
532 write(99,'(2X,E10.4)') vModes%allVar3d(var3dIndex)%autoCovariance(levIndex1,levIndex2)
533 else
534 write(99,'(2X,E10.4,$)') vModes%allVar3d(var3dIndex)%autoCovariance(levIndex1,levIndex2)
535 end if
536 end do
537 end do
538 close(unit=99)
539 outfilename = "./vModes_eigenVectors_"//trim(vModes%allVar3d(var3dIndex)%varName)//".txt"
540 open (unit=99,file=outfilename,action="write",status="new")
541 do levIndex2 = 1, vModes%allVar3d(var3dIndex)%nLev
542 do levIndex1 = 1, vModes%allVar3d(var3dIndex)%nLev
543 if (levIndex1 == vModes%allVar3d(var3dIndex)%nLev) then
544 write(99,'(2X,E10.4)') vModes%allVar3d(var3dIndex)%eigenVectors(levIndex1,levIndex2)
545 else
546 write(99,'(2X,E10.4,$)') vModes%allVar3d(var3dIndex)%eigenVectors(levIndex1,levIndex2)
547 end if
548 end do
549 end do
550 close(unit=99)
551 outfilename = "./vModes_eigenValues_"//trim(vModes%allVar3d(var3dIndex)%varName)//".txt"
552 open (unit=99,file=outfilename,action="write",status="new")
553 do levIndex1 = 1, vModes%allVar3d(var3dIndex)%nLev
554 write(99,'(2X,E10.4)') vModes%allVar3d(var3dIndex)%eigenValues(levIndex1)
555 end do
556 close(unit=99)
557 end do
558 end if
559
560 end subroutine vms_writeModes
561
562 !--------------------------------------------------------------------------
563 ! writeMatrix2d_r8
564 !--------------------------------------------------------------------------
565 subroutine writeMatrix2d_r8(matrix2d,rank,iun,nomvar_in,etiket_in)
566 !
567 !:Purpose: To write a 2D matrix rank x rank
568 !
569 implicit none
570
571 ! Arguments:
572 integer, intent(in) :: rank
573 real(8), intent(in) :: matrix2d(rank,rank)
574 integer, intent(in) :: iun
575 character(len=*), intent(in) :: nomvar_in
576 character(len=*), intent(in) :: etiket_in
577
578 ! Locals:
579 real(4), allocatable :: workecr(:,:)
580 real(4) :: work
581 integer :: errorID, fstecr
582 integer :: dateo, npak, ni, nj, nk
583 integer :: ip1, ip2, ip3, deet, npas, datyp
584 integer :: ig1 ,ig2 ,ig3 ,ig4
585 character(len=1) :: grtyp
586 character(len=4) :: nomvar
587 character(len=2) :: typvar
588 character(len=12) :: etiket
589
590 allocate(workecr(rank, rank))
591
592 npak = -32
593 dateo = 0
594 deet = 0
595 npas = 0
596 ni = rank
597 nj = rank
598 nk = 1
599 ip1 = 0
600 ip2 = 0
601 ip3 = 0
602 typvar = 'XX'
603 nomvar = nomvar_in
604 etiket = etiket_in
605 grtyp = 'X'
606 ig1 = 0
607 ig2 = 0
608 ig3 = 0
609 ig4 = 0
610 datyp = 5
611
612 !- Convert to real 4
613 workecr(:,:) = real(matrix2d(:,:),4)
614
615 !- Writing
616 errorID = fstecr(workecr, work, npak, iun, dateo, deet, npas, ni, nj, &
617 nk, ip1, ip2, ip3, typvar, nomvar, etiket, grtyp, &
618 ig1, ig2, ig3, ig4, datyp, .true.)
619
620 deallocate(workecr)
621
622 end subroutine writeMatrix2d_r8
623
624 !--------------------------------------------------------------------------
625 ! writeArray1d_r8
626 !--------------------------------------------------------------------------
627 subroutine writeArray1d_r8(array,size,iun,nomvar_in,etiket_in)
628 !
629 !:Purpose: To write a 1D array in real 8
630 !
631 implicit none
632
633 ! Arguments:
634 integer, intent(in) :: size
635 real(8), intent(in) :: array(size)
636 integer, intent(in) :: iun
637 character(len=*), intent(in) :: nomvar_in
638 character(len=*), intent(in) :: etiket_in
639
640 ! Locals:
641 real(4), allocatable :: workecr(:)
642 real(4) :: work
643 integer :: errorID, fstecr
644 integer :: dateo, npak, ni, nj, nk
645 integer :: ip1, ip2, ip3, deet, npas, datyp
646 integer :: ig1 ,ig2 ,ig3 ,ig4
647 character(len=1) :: grtyp
648 character(len=4) :: nomvar
649 character(len=2) :: typvar
650 character(len=12) :: etiket
651
652 allocate(workecr(size))
653
654 npak = -32
655 dateo = 0
656 deet = 0
657 npas = 0
658 ni = size
659 nj = 1
660 nk = 1
661 ip1 = 0
662 ip2 = 0
663 ip3 = 0
664 typvar = 'XX'
665 nomvar = nomvar_in
666 etiket = etiket_in
667 grtyp = 'X'
668 ig1 = 0
669 ig2 = 0
670 ig3 = 0
671 ig4 = 0
672 datyp = 5
673
674 !- Convert to real 4
675 workecr(:) = real(array(:),4)
676
677 !- Writing
678 errorID = fstecr(workecr, work, npak, iun, dateo, deet, npas, ni, nj, &
679 nk, ip1, ip2, ip3, typvar, nomvar, etiket, grtyp, &
680 ig1, ig2, ig3, ig4, datyp, .true.)
681
682 deallocate(workecr)
683
684 end subroutine writeArray1d_r8
685
686 !--------------------------------------------------------------------------
687 ! writeArray1d_c4
688 !--------------------------------------------------------------------------
689 subroutine writeArray1d_c4(array,size,iun,nomvar_in,etiket_in)
690 !
691 !:Purpose: To write a 1D array in character with len = 4
692 !
693 implicit none
694
695 ! Arguments:
696 integer, intent(in) :: size
697 character(len=4), intent(in) :: array(size)
698 integer, intent(in) :: iun
699 character(len=*), intent(in) :: nomvar_in
700 character(len=*), intent(in) :: etiket_in
701
702 ! Locals:
703 real(4) :: work
704 integer :: errorID, fstecr_s
705 integer :: dateo, npak, ni, nj, nk
706 integer :: ip1, ip2, ip3, deet, npas, datyp
707 integer :: ig1 ,ig2 ,ig3 ,ig4
708 character(len=1) :: grtyp
709 character(len=4) :: nomvar
710 character(len=2) :: typvar
711 character(len=12) :: etiket
712
713 npak = -32
714 dateo = 0
715 deet = 0
716 npas = 0
717 ni = 4 ! 4 Characters
718 nj = size
719 nk = 1
720 ip1 = 0
721 ip2 = 0
722 ip3 = 0
723 typvar = 'XX'
724 nomvar = nomvar_in
725 etiket = etiket_in
726 grtyp = 'X'
727 ig1 = 0
728 ig2 = 0
729 ig3 = 0
730 ig4 = 0
731 datyp = 7 ! Character
732
733 !- Writing
734 errorID = fstecr_s(array, work, npak, iun, dateo, deet, npas, ni, nj, &
735 nk, ip1, ip2, ip3, typvar, nomvar, etiket, grtyp, &
736 ig1, ig2, ig3, ig4, datyp, .true.)
737
738 end subroutine writeArray1d_c4
739
740end module verticalModes_mod