m_parameters.f90 8.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322
  1. # 0 "<stdin>"
  2. # 0 "<built-in>"
  3. # 0 "<command-line>"
  4. # 1 "/usr/include/stdc-predef.h" 1 3 4
  5. # 17 "/usr/include/stdc-predef.h" 3 4
  6. # 2 "<command-line>" 2
  7. # 1 "<stdin>"
  8. # 10 "<stdin>"
  9. ! File: m_parameters.F90
  10. !
  11. ! Created: 6 August 2010
  12. !
  13. ! Last modified: 6/8/2010
  14. !
  15. ! Author: Pavel Sakov
  16. ! NERSC
  17. !
  18. ! Purpose: Provide a simpl nml list-based parameter input into EnKF.
  19. !
  20. ! Description: Provides code for reading parameters from a specified
  21. ! parameter file.
  22. !
  23. ! Modifications: none
  24. module m_parameters
  25. use qmpi
  26. implicit none
  27. integer, parameter, private :: STRLEN = 512
  28. integer, parameter, private :: FID = 101
  29. character(STRLEN), public :: PRMFNAME = "NONE"
  30. integer, public :: ENSSIZE = 0
  31. namelist /ensemble/ ENSSIZE
  32. character(STRLEN), public :: METHODTAG = "NONE"
  33. namelist /method/ METHODTAG
  34. real, public :: LOCRAD = 0.0d0
  35. character(STRLEN), public :: LOCFUNTAG = "GASPARI-COHN"
  36. namelist /localisation/ LOCRAD, LOCFUNTAG
  37. real, public :: INFL = 1.0d0
  38. real, public :: RFACTOR1 = 1.0d0
  39. real, public :: RFACTOR2 = 1.0d0
  40. real, public :: KFACTOR = 2.0d0
  41. namelist /moderation/ INFL, RFACTOR1, RFACTOR2, KFACTOR
  42. character(STRLEN), public :: JMAPFNAME = "NONE"
  43. character(STRLEN), public :: POINTFNAME = "NONE"
  44. character(STRLEN), public :: MEANSSHFNAME = "NONE"
  45. namelist /files/ JMAPFNAME, POINTFNAME, MEANSSHFNAME
  46. integer, parameter, private :: NPRMESTMAX = 10
  47. integer :: nprmest = 0
  48. character(STRLEN), dimension(NPRMESTMAX), public :: PRMESTNAME
  49. real, dimension(NPRMESTMAX), public :: PRMINFL
  50. namelist /prmest/ PRMESTNAME, PRMINFL
  51. public prm_read, prm_describe, prm_print, prm_getinfl, prm_prmestexists, ucase
  52. contains
  53. subroutine prm_read
  54. integer :: ios, i
  55. call getarg(1, PRMFNAME)
  56. if (master) then
  57. print *, 'EnKF: reading parameters from "', trim(PRMFNAME), '":'
  58. end if
  59. open(unit = FID, file = trim(PRMFNAME), form = "formatted",&
  60. status = "old", iostat = ios)
  61. if (ios /= 0) then
  62. if (master) then
  63. print *, 'ERROR: could not open "', trim(PRMFNAME), '", iostatus =', ios
  64. stop
  65. end if
  66. end if
  67. read(unit = FID, nml = method, iostat = ios)
  68. if (ios /= 0) then
  69. if (master) then
  70. print *, 'ERROR: "', trim(PRMFNAME), '": could not read namelist "method"'
  71. end if
  72. stop
  73. end if
  74. rewind(FID)
  75. read(unit = FID, nml = ensemble, iostat = ios)
  76. if (ios /= 0) then
  77. if (master) then
  78. print *, 'ERROR: "', trim(PRMFNAME), '": could not read namelist "ensemble"'
  79. end if
  80. stop
  81. end if
  82. rewind(FID)
  83. read(unit = FID, nml = localisation, iostat = ios)
  84. if (ios /= 0) then
  85. if (master) then
  86. print *, 'ERROR: "', trim(PRMFNAME), '": could not read namelist "localisation"'
  87. end if
  88. stop
  89. end if
  90. rewind(FID)
  91. read(unit = FID, nml = moderation, iostat = ios)
  92. if (ios /= 0) then
  93. if (master) then
  94. print *, 'ERROR: "', trim(PRMFNAME), '": could not read namelist "moderation"'
  95. end if
  96. stop
  97. end if
  98. rewind(FID)
  99. read(unit = FID, nml = files, iostat = ios)
  100. if (ios /= 0) then
  101. if (master) then
  102. print *, 'ERROR: "', trim(PRMFNAME), '": could not read namelist "files"'
  103. end if
  104. stop
  105. end if
  106. rewind(FID)
  107. do i = 1, NPRMESTMAX
  108. PRMESTNAME(i) = ""
  109. end do
  110. read(unit = FID, nml = prmest, iostat = ios)
  111. if (ios /= 0) then
  112. if (master) then
  113. print *, 'ERROR: "', trim(PRMFNAME), '": could not read namelist "prmest"'
  114. end if
  115. stop
  116. end if
  117. do i = 1, NPRMESTMAX
  118. if (PRMESTNAME(i) == "") then
  119. nprmest = i - 1
  120. exit
  121. end if
  122. end do
  123. rewind(FID)
  124. close(FID)
  125. call ucase(METHODTAG)
  126. call ucase(LOCFUNTAG)
  127. end subroutine prm_read
  128. subroutine prm_describe
  129. if (.not. master) then
  130. return
  131. end if
  132. print '(a)', ' Example of EnKF parameter file:'
  133. print *
  134. print '(a)', '&method'
  135. print '(a)', ' methodtag = "DEnKF"'
  136. print '(a)', '/'
  137. print '(a)', '&ensemble'
  138. print '(a)', ' enssize = 0'
  139. print '(a)', '/'
  140. print '(a)', '&localisation'
  141. print '(a)', ' locfuntag = "Gaspari-Cohn"'
  142. print '(a)', ' locrad = 300.0'
  143. print '(a)', '/'
  144. print '(a)', '&moderation'
  145. print '(a)', ' infl = 1.01 (<number>)'
  146. print '(a)', ' rfactor1 = 1.0 (<number>)'
  147. print '(a)', ' rfactor2 = 2.0 (<number>)'
  148. print '(a)', ' kfactor = 2.0 (<number>)'
  149. print '(a)', '/'
  150. print '(a)', '&files'
  151. print '(a)', ' jmapfname = "jmap.txt" (<file name>)'
  152. print '(a)', ' pointfname = "point2nc.txt" (<file name>)'
  153. print '(a)', ' meansshfname = "meanssh.uf" (<file name>)'
  154. print *
  155. print '(a)', 'Parameter options:'
  156. print '(a)', ' method = "EnKF" | "DEnKF"*'
  157. print '(a)', ' enssize = <number> (0* to use all available states)'
  158. print '(a)', ' locfuntag = "Gaspari-Cohn"* | "Step" | "None"'
  159. print '(a)', ' locrad = <support radius in km>'
  160. print '(a)', ' infl = <multiple, for ensemble anomalies> (* 1.0)'
  161. print '(a)', ' rfactor1 = <obs. error variance multiple> (* 1.0)'
  162. print '(a)', ' rfactor2 = <additional multiple for updating ens. anomalies> (* 1.0)'
  163. print '(a)', ' kfactor = <max. allowed increment in terms of ensemble spread> (* 2.0)'
  164. print '(a)', ' jmapfname* = <file with j remapping> (* none)'
  165. print '(a)', ' pointfname* = <file with point coordinates> (* none)'
  166. print '(a)', ' meansshfname* = <file with mean SSH> (* none)'
  167. end subroutine prm_describe
  168. subroutine prm_print
  169. integer :: i
  170. if (.not. master) then
  171. return
  172. end if
  173. print '(a)', ' EnKF parameters:'
  174. print '(a)', ' method:'
  175. print '(a, a, a)', ' methodtag = "', trim(METHODTAG), '"'
  176. print '(a)', ' ensemble:'
  177. print '(a, i0)', ' enssize = ', ENSSIZE
  178. print '(a)', ' localisation:'
  179. print '(a, f5.0)', ' locrad = ', LOCRAD
  180. print '(a, a, a)', ' locfuntag = "', trim(LOCFUNTAG), '"'
  181. print '(a)', ' moderation:'
  182. print '(a, f5.3)', ' infl = ', INFL
  183. print '(a, f3.1)', ' rfactor1 = ', RFACTOR1
  184. print '(a, f3.1)', ' rfactor2 = ', RFACTOR2
  185. print '(a, f3.1)', ' kfactor = ', KFACTOR
  186. print '(a)', ' files:'
  187. print '(a, a, a)', ' jmapfname = "', trim(JMAPFNAME), '"'
  188. print '(a, a, a)', ' pointfname = "', trim(POINTFNAME), '"'
  189. print '(a, a, a)', ' meansshfname = "', trim(MEANSSHFNAME), '"'
  190. print '(a, i0, a)', ' prmest: ', nprmest, ' fields'
  191. do i = 1, nprmest
  192. print '(a, a, a, f5.3)', ' prmestname = "', trim(PRMESTNAME(i)), '", infl = ', PRMINFL(i)
  193. end do
  194. print *
  195. end subroutine prm_print
  196. function prm_getinfl(fldname)
  197. real :: prm_getinfl
  198. character(*), intent(in) :: fldname
  199. integer :: i
  200. prm_getinfl = INFL
  201. do i = 1, nprmest
  202. if (trim(fldname) == PRMESTNAME(i)) then
  203. prm_getinfl = PRMINFL(i)
  204. print '(a, a, a, f5.3)', ' "', trim(fldname), '": using inflation = ', prm_getinfl
  205. return
  206. end if
  207. end do
  208. end function prm_getinfl
  209. function prm_prmestexists(varname)
  210. logical :: prm_prmestexists
  211. character(*), intent(in) :: varname
  212. integer :: i
  213. prm_prmestexists = .false.
  214. do i = 1, nprmest
  215. if (trim(varname) == PRMESTNAME(i)) then
  216. prm_prmestexists = .true.
  217. return
  218. end if
  219. end do
  220. end function prm_prmestexists
  221. ! Shift a character string to upper case.
  222. !
  223. subroutine ucase(string)
  224. character(*) :: string
  225. integer :: i
  226. do i = 1, len(string)
  227. if (string(i:i) >= 'a' .and. string(i:i) <= 'z') then
  228. string(i:i) = achar (iachar ( string(i:i) ) - 32)
  229. end if
  230. end do
  231. end subroutine ucase
  232. end module m_parameters