00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040 module ice_history
00041
00042
00043
00044 use ice_kinds_mod
00045 use ice_broadcast
00046 use ice_communicate, only: my_task, master_task
00047 use ice_blocks
00048 use ice_grid
00049 use ice_read_write
00050 use ice_fileunits
00051 use ice_history_fields
00052 use ice_history_write
00053
00054
00055
00056 implicit none
00057 save
00058
00059 character (len=16) :: vname_in
00060 character (len=55) :: vdesc_in
00061 character (len=55) :: vcomment_in
00062
00063
00064
00065
00066
00067 character (len=16), parameter ::
00068 tstr = 'TLON TLAT time',
00069 ustr = 'ULON ULAT time',
00070 tcstr = 'area: tarea' ,
00071 ucstr = 'area: uarea'
00072
00073 real (kind=dbl_kind) ::
00074 avgct(max_nstrm)
00075
00076
00077
00078
00079
00080 logical (kind=log_kind) ::
00081 f_tmask = .true.,
00082 f_tarea = .true., f_uarea = .true.,
00083 f_dxt = .true., f_dyt = .true.,
00084 f_dxu = .true., f_dyu = .true.,
00085 f_HTN = .true., f_HTE = .true.,
00086 f_ANGLE = .true., f_ANGLET = .true.
00087
00088 character (len=max_nstrm) ::
00089
00090 f_hi = 'mxxxx', f_hs = 'mxxxx', &
00091 f_fs = 'mxxxx', &
00092 f_Tsfc = 'mxxxx', f_aice = 'mxxxx', &
00093 f_uvel = 'mxxxx', f_vvel = 'mxxxx', &
00094 f_transix = 'mxxxx', f_transiy = 'mxxxx', &
00095 f_fswdn = 'mxxxx', f_fswup = 'mxxxx', &
00096 f_flwdn = 'mxxxx', &
00097 f_snow = 'mxxxx', f_snow_ai = 'mxxxx', &
00098 f_rain = 'mxxxx', f_rain_ai = 'mxxxx', &
00099 f_faero_atm = 'mxxxx', f_faero_ocn = 'mxxxx', &
00100 f_sst = 'mxxxx', f_sss = 'mxxxx', &
00101 f_uocn = 'mxxxx', f_vocn = 'mxxxx', &
00102 f_frzmlt = 'mxxxx', &
00103 f_fswfac = 'mxxxx', &
00104 f_fswabs = 'mxxxx', f_fswabs_ai = 'mxxxx', &
00105 f_alvdr = 'mxxxx', f_alidr = 'mxxxx', &
00106 f_albice = 'mxxxx', f_albsno = 'mxxxx', &
00107 f_albpnd = 'mxxxx', f_coszen = 'mxxxx', &
00108 f_flat = 'mxxxx', f_flat_ai = 'mxxxx', &
00109 f_fsens = 'mxxxx', f_fsens_ai = 'mxxxx', &
00110 f_flwup = 'mxxxx', f_flwup_ai = 'mxxxx', &
00111 f_evap = 'mxxxx', f_evap_ai = 'mxxxx', &
00112 f_qi = 'mxxxx', f_qs = 'mxxxx', &
00113 f_Tair = 'mxxxx', &
00114 f_Tref = 'mxxxx', f_Qref = 'mxxxx', &
00115 f_congel = 'mxxxx', f_frazil = 'mxxxx', &
00116 f_snoice = 'mxxxx', &
00117 f_meltt = 'mxxxx', f_melts = 'mxxxx', &
00118 f_meltb = 'mxxxx', f_meltl = 'mxxxx', &
00119 f_fresh = 'mxxxx', f_fresh_ai = 'mxxxx', &
00120 f_fsalt = 'mxxxx', f_fsalt_ai = 'mxxxx', &
00121 f_fhocn = 'mxxxx', f_fhocn_ai = 'mxxxx', &
00122 f_fswthru = 'mxxxx', f_fswthru_ai = 'mxxxx', &
00123 f_strairx = 'mxxxx', f_strairy = 'mxxxx', &
00124 f_strtltx = 'mxxxx', f_strtlty = 'mxxxx', &
00125 f_strcorx = 'mxxxx', f_strcory = 'mxxxx', &
00126 f_strocnx = 'mxxxx', f_strocny = 'mxxxx', &
00127 f_strintx = 'mxxxx', f_strinty = 'mxxxx', &
00128 f_strength = 'mxxxx', f_opening = 'mxxxx', &
00129 f_divu = 'mxxxx', f_shear = 'mxxxx', &
00130 f_sig1 = 'mxxxx', f_sig2 = 'mxxxx', &
00131 f_dvidtt = 'mxxxx', f_dvidtd = 'mxxxx', &
00132 f_daidtt = 'mxxxx', f_daidtd = 'mxxxx', &
00133 f_mlt_onset = 'mxxxx', f_frz_onset = 'mxxxx', &
00134 f_dardg1dt = 'mxxxx', f_dardg2dt = 'mxxxx', &
00135 f_dvirdgdt = 'mxxxx', f_iage = 'mxxxx', &
00136 f_FY = 'mxxxx', &
00137 f_aeron = 'xxxxx', f_aero = 'xxxxx', &
00138 f_apond = 'xxxxx', f_apondn = 'xxxxx', &
00139 f_hisnap = 'mxxxx', f_aisnap = 'mxxxx', &
00140 f_aicen = 'mxxxx', f_vicen = 'mxxxx', &
00141 f_trsig = 'mxxxx', f_icepresent = 'mxxxx', &
00142 f_fsurf_ai = 'mxxxx', f_fcondtop_ai= 'mxxxx', &
00143 f_fmeltt_ai = 'xxxxx', &
00144 f_fsurfn_ai = 'xxxxx',f_fcondtopn_ai= 'xxxxx', &
00145 f_fmelttn_ai= 'xxxxx', f_flatn_ai = 'xxxxx'
00146
00147
00148
00149
00150
00151 namelist / icefields_nml / &
00152 f_tmask , &
00153 f_tarea , f_uarea , &
00154 f_dxt , f_dyt , &
00155 f_dxu , f_dyu , &
00156 f_HTN , f_HTE , &
00157 f_ANGLE , f_ANGLET , &
00158 f_bounds , &
00159
00160
00161 f_hi, f_hs , &
00162 f_fs, &
00163 f_Tsfc, f_aice , &
00164 f_uvel, f_vvel , &
00165 f_transix, f_transiy , &
00166 f_fswdn, f_fswup , &
00167 f_flwdn, &
00168 f_snow, f_snow_ai , &
00169 f_rain, f_rain_ai , &
00170 f_faero_atm, f_faero_ocn, &
00171 f_sst, f_sss , &
00172 f_uocn, f_vocn , &
00173 f_frzmlt , &
00174 f_fswfac , &
00175 f_fswabs, f_fswabs_ai, &
00176 f_alvdr, f_alidr , &
00177 f_albice, f_albsno , &
00178 f_albpnd, f_coszen , &
00179 f_flat, f_flat_ai , &
00180 f_fsens, f_fsens_ai , &
00181 f_flwup, f_flwup_ai , &
00182 f_evap, f_evap_ai , &
00183 f_qi, f_qs , &
00184 f_Tair , &
00185 f_Tref, f_Qref , &
00186 f_congel, f_frazil , &
00187 f_snoice, &
00188 f_meltt, f_melts , &
00189 f_meltb, f_meltl , &
00190 f_fresh, f_fresh_ai , &
00191 f_fsalt, f_fsalt_ai , &
00192 f_fhocn, f_fhocn_ai , &
00193 f_fswthru, f_fswthru_ai,&
00194 f_strairx, f_strairy , &
00195 f_strtltx, f_strtlty , &
00196 f_strcorx, f_strcory , &
00197 f_strocnx, f_strocny , &
00198 f_strintx, f_strinty , &
00199 f_strength, f_opening , &
00200 f_divu, f_shear , &
00201 f_sig1, f_sig2 , &
00202 f_dvidtt, f_dvidtd , &
00203 f_daidtt, f_daidtd , &
00204 f_mlt_onset, f_frz_onset, &
00205 f_dardg1dt, f_dardg2dt , &
00206 f_dvirdgdt , &
00207 f_hisnap, f_aisnap , &
00208 f_aicen, f_vicen , &
00209 f_aeron, f_aero , &
00210 f_iage , &
00211 f_FY , &
00212 f_apond, f_apondn , &
00213 f_trsig, f_icepresent,&
00214 f_fsurf_ai, f_fcondtop_ai,&
00215 f_fmeltt_ai, &
00216 f_fsurfn_ai,f_fcondtopn_ai,&
00217 f_fmelttn_ai,f_flatn_ai
00218
00219
00220
00221 contains
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231 subroutine init_hist (dt)
00232
00233
00234
00235
00236
00237
00238
00239
00240
00241
00242
00243
00244
00245
00246
00247 use ice_constants
00248 use ice_calendar, only: yday, days_per_year, histfreq, &
00249 histfreq_n, nstreams
00250 use ice_flux, only: mlt_onset, frz_onset
00251 use ice_restart, only: restart
00252 use ice_state, only: tr_aero, tr_iage, tr_FY, tr_pond
00253 use ice_exit
00254
00255
00256
00257 real (kind=dbl_kind), intent(in) ::
00258 dt
00259
00260
00261
00262 integer (kind=int_kind) :: n, k, ns, ns1, ns2, lenf
00263 integer (kind=int_kind) :: hfreqn
00264 integer (kind=int_kind), dimension(max_nstrm) ::
00265 ntmp
00266 integer (kind=int_kind) :: nml_error
00267
00268 character (len=3) :: nchar
00269 character (len=40) :: stmp
00270
00271
00272
00273
00274
00275 call get_fileunit(nu_nml)
00276 if (my_task == master_task) then
00277 open (nu_nml, file=nml_filename, status='old',iostat=nml_error)
00278 if (nml_error /= 0) then
00279 nml_error = -1
00280 else
00281 nml_error = 1
00282 endif
00283 do while (nml_error > 0)
00284 read(nu_nml, nml=icefields_nml,iostat=nml_error)
00285 if (nml_error > 0) read(nu_nml,*)
00286 end do
00287 if (nml_error == 0) close(nu_nml)
00288 endif
00289 call release_fileunit(nu_nml)
00290
00291 call broadcast_scalar(nml_error, master_task)
00292 if (nml_error /= 0) then
00293 close (nu_nml)
00294 call abort_ice('ice: error reading icefields_nml')
00295 endif
00296
00297
00298 nstreams = 0
00299 do ns = 1, max_nstrm
00300 if (histfreq(ns) == '1' .or. histfreq(ns) == 'h' .or. &
00301 histfreq(ns) == 'd' .or. histfreq(ns) == 'm' .or. &
00302 histfreq(ns) == 'y') then
00303 nstreams = nstreams + 1
00304 else if (histfreq(ns) /= 'x') then
00305 call abort_ice('ice: histfreq contains illegal element')
00306 endif
00307 enddo
00308 if (nstreams == 0) write (nu_diag,*) 'WARNING: No history output'
00309 do ns1 = 1, nstreams
00310 do ns2 = 1, nstreams
00311 if (histfreq(ns1) == histfreq(ns2) .and. ns1/=ns2 &
00312 .and. my_task == master_task) then
00313 call abort_ice('ice: histfreq elements must be unique')
00314 endif
00315 enddo
00316 enddo
00317
00318 if (.not. tr_iage) f_iage = 'xxxxx'
00319 if (.not. tr_FY) f_FY = 'xxxxx'
00320 if (.not. tr_pond) then
00321 f_apond = 'xxxxx'
00322 f_apondn = 'xxxxx'
00323 endif
00324 if (.not. tr_aero) then
00325 f_faero_atm = 'xxxxx'
00326 f_faero_ocn = 'xxxxx'
00327 f_aero = 'xxxxx'
00328 f_aeron = 'xxxxx'
00329 endif
00330
00331
00332
00333 if (f_albsno(1:1) /= 'x') f_albsno = f_albice
00334 if (f_albpnd(1:1) /= 'x') f_albpnd = f_albice
00335 if (f_coszen(1:1) /= 'x') f_coszen = f_albice
00336
00337
00338 if (f_fmeltt_ai(1:1) /= 'x') f_fmelttn_ai = f_fmeltt_ai
00339
00340 #ifndef ncdf
00341 f_bounds = .false.
00342 #endif
00343
00344 call broadcast_scalar (f_tmask, master_task)
00345 call broadcast_scalar (f_tarea, master_task)
00346 call broadcast_scalar (f_uarea, master_task)
00347 call broadcast_scalar (f_dxt, master_task)
00348 call broadcast_scalar (f_dyt, master_task)
00349 call broadcast_scalar (f_dxu, master_task)
00350 call broadcast_scalar (f_dyu, master_task)
00351 call broadcast_scalar (f_HTN, master_task)
00352 call broadcast_scalar (f_HTE, master_task)
00353 call broadcast_scalar (f_ANGLE, master_task)
00354 call broadcast_scalar (f_ANGLET, master_task)
00355 call broadcast_scalar (f_bounds, master_task)
00356
00357
00358 call broadcast_scalar (f_hi, master_task)
00359 call broadcast_scalar (f_hs, master_task)
00360 call broadcast_scalar (f_fs, master_task)
00361 call broadcast_scalar (f_Tsfc, master_task)
00362 call broadcast_scalar (f_aice, master_task)
00363 call broadcast_scalar (f_uvel, master_task)
00364 call broadcast_scalar (f_vvel, master_task)
00365 call broadcast_scalar (f_transix, master_task)
00366 call broadcast_scalar (f_transiy, master_task)
00367 call broadcast_scalar (f_fswdn, master_task)
00368 call broadcast_scalar (f_fswup, master_task)
00369 call broadcast_scalar (f_flwdn, master_task)
00370 call broadcast_scalar (f_snow, master_task)
00371 call broadcast_scalar (f_snow_ai, master_task)
00372 call broadcast_scalar (f_rain, master_task)
00373 call broadcast_scalar (f_rain_ai, master_task)
00374 call broadcast_scalar (f_faero_atm, master_task)
00375 call broadcast_scalar (f_faero_ocn, master_task)
00376 call broadcast_scalar (f_sst, master_task)
00377 call broadcast_scalar (f_sss, master_task)
00378 call broadcast_scalar (f_uocn, master_task)
00379 call broadcast_scalar (f_vocn, master_task)
00380 call broadcast_scalar (f_frzmlt, master_task)
00381 call broadcast_scalar (f_fswfac, master_task)
00382 call broadcast_scalar (f_fswabs, master_task)
00383 call broadcast_scalar (f_fswabs_ai, master_task)
00384 call broadcast_scalar (f_alvdr, master_task)
00385 call broadcast_scalar (f_alidr, master_task)
00386 call broadcast_scalar (f_albice, master_task)
00387 call broadcast_scalar (f_albsno, master_task)
00388 call broadcast_scalar (f_albpnd, master_task)
00389 call broadcast_scalar (f_coszen, master_task)
00390 call broadcast_scalar (f_flat, master_task)
00391 call broadcast_scalar (f_flat_ai, master_task)
00392 call broadcast_scalar (f_fsens, master_task)
00393 call broadcast_scalar (f_fsens_ai, master_task)
00394 call broadcast_scalar (f_flwup, master_task)
00395 call broadcast_scalar (f_flwup_ai, master_task)
00396 call broadcast_scalar (f_evap, master_task)
00397 call broadcast_scalar (f_evap_ai, master_task)
00398 call broadcast_scalar (f_qi, master_task)
00399 call broadcast_scalar (f_qs, master_task)
00400 call broadcast_scalar (f_Tair, master_task)
00401 call broadcast_scalar (f_Tref, master_task)
00402 call broadcast_scalar (f_Qref, master_task)
00403 call broadcast_scalar (f_congel, master_task)
00404 call broadcast_scalar (f_frazil, master_task)
00405 call broadcast_scalar (f_snoice, master_task)
00406 call broadcast_scalar (f_meltt, master_task)
00407 call broadcast_scalar (f_meltb, master_task)
00408 call broadcast_scalar (f_meltl, master_task)
00409 call broadcast_scalar (f_melts, master_task)
00410 call broadcast_scalar (f_fresh, master_task)
00411 call broadcast_scalar (f_fresh_ai, master_task)
00412 call broadcast_scalar (f_fsalt, master_task)
00413 call broadcast_scalar (f_fsalt_ai, master_task)
00414 call broadcast_scalar (f_fhocn, master_task)
00415 call broadcast_scalar (f_fhocn_ai, master_task)
00416 call broadcast_scalar (f_fswthru, master_task)
00417 call broadcast_scalar (f_fswthru_ai, master_task)
00418 call broadcast_scalar (f_strairx, master_task)
00419 call broadcast_scalar (f_strairy, master_task)
00420 call broadcast_scalar (f_strtltx, master_task)
00421 call broadcast_scalar (f_strtlty, master_task)
00422 call broadcast_scalar (f_strcorx, master_task)
00423 call broadcast_scalar (f_strcory, master_task)
00424 call broadcast_scalar (f_strocnx, master_task)
00425 call broadcast_scalar (f_strocny, master_task)
00426 call broadcast_scalar (f_strintx, master_task)
00427 call broadcast_scalar (f_strinty, master_task)
00428 call broadcast_scalar (f_strength, master_task)
00429 call broadcast_scalar (f_opening, master_task)
00430 call broadcast_scalar (f_divu, master_task)
00431 call broadcast_scalar (f_shear, master_task)
00432 call broadcast_scalar (f_sig1, master_task)
00433 call broadcast_scalar (f_sig2, master_task)
00434 call broadcast_scalar (f_dvidtt, master_task)
00435 call broadcast_scalar (f_dvidtd, master_task)
00436 call broadcast_scalar (f_daidtt, master_task)
00437 call broadcast_scalar (f_daidtd, master_task)
00438 call broadcast_scalar (f_mlt_onset, master_task)
00439 call broadcast_scalar (f_frz_onset, master_task)
00440 call broadcast_scalar (f_dardg1dt, master_task)
00441 call broadcast_scalar (f_dardg2dt, master_task)
00442 call broadcast_scalar (f_dvirdgdt, master_task)
00443 call broadcast_scalar (f_aisnap, master_task)
00444 call broadcast_scalar (f_hisnap, master_task)
00445 call broadcast_scalar (f_aicen, master_task)
00446 call broadcast_scalar (f_vicen, master_task)
00447 call broadcast_scalar (f_trsig, master_task)
00448 call broadcast_scalar (f_icepresent, master_task)
00449 call broadcast_scalar (f_fsurf_ai, master_task)
00450 call broadcast_scalar (f_fcondtop_ai, master_task)
00451 call broadcast_scalar (f_fmeltt_ai, master_task)
00452 call broadcast_scalar (f_fsurfn_ai, master_task)
00453 call broadcast_scalar (f_fcondtopn_ai, master_task)
00454 call broadcast_scalar (f_fmelttn_ai, master_task)
00455 call broadcast_scalar (f_flatn_ai, master_task)
00456
00457 call broadcast_scalar (f_aero, master_task)
00458 call broadcast_scalar (f_aeron, master_task)
00459 call broadcast_scalar (f_iage, master_task)
00460 call broadcast_scalar (f_FY, master_task)
00461 call broadcast_scalar (f_apond, master_task)
00462 call broadcast_scalar (f_apondn, master_task)
00463
00464 do ns1 = 1, nstreams
00465
00466
00467
00468
00469
00470
00471
00472
00473
00474 if (f_hi(1:1) /= 'x') &
00475 call define_hist_field(n_hi,"hi","m",tstr, tcstr, &
00476 "grid cell mean ice thickness", &
00477 "ice volume per unit grid cell area", c1, c0, &
00478 ns1, f_hi)
00479
00480 if (f_hs(1:1) /= 'x') &
00481 call define_hist_field(n_hs,"hs","m",tstr, tcstr, &
00482 "grid cell mean snow thickness", &
00483 "snow volume per unit grid cell area", c1, c0, &
00484 ns1, f_hs)
00485
00486 if (f_fs(1:1) /= 'x') &
00487 call define_hist_field(n_fs,"fs"," ",tstr, tcstr, &
00488 "grid cell mean snow fraction", &
00489 "none", c1, c0, &
00490 ns1, f_fs)
00491
00492 if (f_Tsfc(1:1) /= 'x') &
00493 call define_hist_field(n_Tsfc,"Tsfc","degC",tstr, tcstr, &
00494 "snow/ice surface temperature", &
00495 "averaged with Tf if no ice is present", c1, c0, &
00496 ns1, f_Tsfc)
00497
00498 if (f_aice(1:1) /= 'x') &
00499 call define_hist_field(n_aice,"aice","%",tstr, tcstr, &
00500 "ice area (aggregate)", &
00501 "none", c100, c0, &
00502 ns1, f_aice)
00503
00504 if (f_qi(1:1) /= 'x') &
00505 call define_hist_field(n_qi,"qi","J",tstr, tcstr, &
00506 "internal ice heat content", &
00507 "none", c1, c0, &
00508 ns1, f_qi)
00509
00510 if (f_qs(1:1) /= 'x') &
00511 call define_hist_field(n_qs,"qs","J",tstr, tcstr, &
00512 "internal snow heat content", &
00513 "none", c1, c0, &
00514 ns1, f_qs)
00515
00516 if (f_uvel(1:1) /= 'x') &
00517 call define_hist_field(n_uvel,"uvel","cm/s",ustr, ucstr, &
00518 "ice velocity (x)", &
00519 "positive is x direction on U grid", m_to_cm, c0, &
00520 ns1, f_uvel)
00521
00522 if (f_vvel(1:1) /= 'x') &
00523 call define_hist_field(n_vvel,"vvel","cm/s",ustr, ucstr, &
00524 "ice velocity (y)", &
00525 "positive is y direction on U grid", m_to_cm, c0, &
00526 ns1, f_vvel)
00527
00528 if (f_transix(1:1) /= 'x') &
00529 call define_hist_field(n_transix,"transix","kg/s",tstr, tcstr, &
00530 "ice mass transport (x) on East side", &
00531 "positive is x direction on U grid", c1, c0, &
00532 ns1, f_transix)
00533
00534 if (f_transiy(1:1) /= 'x') &
00535 call define_hist_field(n_transiy,"transiy","kg/s",tstr, tcstr, &
00536 "ice mass transport (y) on North side", &
00537 "positive is y direction on U grid", c1, c0, &
00538 ns1, f_transiy)
00539
00540 if (f_fswdn(1:1) /= 'x') &
00541 call define_hist_field(n_fswdn,"fswdn","W/m^2",tstr, tcstr, &
00542 "down solar flux", &
00543 "positive downward", c1, c0, &
00544 ns1, f_fswdn)
00545
00546 if (f_fswup(1:1) /= 'x') &
00547 call define_hist_field(n_fswup,"fswup","W/m^2",tstr, tcstr, &
00548 "upward solar flux", &
00549 "positive downward", c1, c0, &
00550 ns1, f_fswup)
00551
00552 if (f_flwdn(1:1) /= 'x') &
00553 call define_hist_field(n_flwdn,"flwdn","W/m^2",tstr, tcstr, &
00554 "down longwave flux", &
00555 "positive downward", c1, c0, &
00556 ns1, f_flwdn)
00557
00558 if (f_snow(1:1) /= 'x') &
00559 call define_hist_field(n_snow,"snow","cm/day",tstr, tcstr, &
00560 "snowfall rate (cpl)", &
00561 "none", mps_to_cmpdy/rhofresh, c0, &
00562 ns1, f_snow)
00563
00564 if (f_snow_ai(1:1) /= 'x') &
00565 call define_hist_field(n_snow_ai,"snow_ai","cm/day",tstr, tcstr, &
00566 "snowfall rate", &
00567 "weighted by ice", mps_to_cmpdy/rhofresh, c0, &
00568 ns1, f_snow_ai)
00569
00570 if (f_rain(1:1) /= 'x') &
00571 call define_hist_field(n_rain,"rain","cm/day",tstr, tcstr, &
00572 "rainfall rate (cpl)", &
00573 "none", mps_to_cmpdy/rhofresh, c0, &
00574 ns1, f_rain)
00575
00576 if (f_rain_ai(1:1) /= 'x') &
00577 call define_hist_field(n_rain_ai,"rain_ai","cm/day",tstr, tcstr, &
00578 "rainfall rate", &
00579 "weighted by ice", mps_to_cmpdy/rhofresh, c0, &
00580 ns1, f_rain_ai)
00581
00582 if (f_sst(1:1) /= 'x') &
00583 call define_hist_field(n_sst,"sst","C",tstr, tcstr, &
00584 "sea surface temperature", &
00585 "none", c1, c0, &
00586 ns1, f_sst)
00587
00588 if (f_sss(1:1) /= 'x') &
00589 call define_hist_field(n_sss,"sss","ppt",tstr, tcstr, &
00590 "sea surface salinity", &
00591 "none", c1, c0, &
00592 ns1, f_sss)
00593
00594 if (f_uocn(1:1) /= 'x') &
00595 call define_hist_field(n_uocn,"uocn","cm/s",ustr, ucstr, &
00596 "ocean current (x)", &
00597 "positive is x direction on U grid", m_to_cm, c0, &
00598 ns1, f_uocn)
00599
00600 if (f_vocn(1:1) /= 'x') &
00601 call define_hist_field(n_vocn,"vocn","cm/s",ustr, ucstr, &
00602 "ocean current (y)", &
00603 "positive is y direction on U grid", m_to_cm, c0, &
00604 ns1, f_vocn)
00605
00606 if (f_fswfac(1:1) /= 'x') &
00607 call define_hist_field(n_fswfac,"fswfac","1",tstr, tcstr, &
00608 "shortwave scaling factor", &
00609 "ratio of netsw new:old", c1, c0, &
00610 ns1, f_fswfac)
00611
00612 if (f_frzmlt(1:1) /= 'x') &
00613 call define_hist_field(n_frzmlt,"frzmlt","W/m^2",tstr, tcstr, &
00614 "freeze/melt potential", &
00615 "if >0, new ice forms; if <0, ice melts", c1, c0, &
00616 ns1, f_frzmlt)
00617
00618 if (f_fswabs(1:1) /= 'x') &
00619 call define_hist_field(n_fswabs,"fswabs","W/m^2",tstr, tcstr, &
00620 "snow/ice/ocn absorbed solar flux (cpl)", &
00621 "positive downward", c1, c0, &
00622 ns1, f_fswabs)
00623
00624 if (f_fswabs_ai(1:1) /= 'x') &
00625 call define_hist_field(n_fswabs_ai,"fswabs_ai","W/m^2",tstr, tcstr, &
00626 "snow/ice/ocn absorbed solar flux", &
00627 "weighted by ice area", c1, c0, &
00628 ns1, f_fswabs_ai)
00629
00630 if (f_alvdr(1:1) /= 'x') &
00631 call define_hist_field(n_alvdr,"alvdr","%",tstr, tcstr, &
00632 "visible direct albedo", &
00633 "none", c100, c0, &
00634 ns1, f_alvdr)
00635
00636 if (f_alidr(1:1) /= 'x') &
00637 call define_hist_field(n_alidr,"alidr","%",tstr, tcstr, &
00638 "near IR direct albedo", &
00639 "none", c100, c0, &
00640 ns1, f_alidr)
00641
00642 if (f_albice(1:1) /= 'x') &
00643 call define_hist_field(n_albice,"albice","%",tstr, tcstr, &
00644 "bare ice albedo", &
00645 "averaged for coszen>0, weighted by aice", c100, c0, &
00646 ns1, f_albice)
00647
00648 if (f_albsno(1:1) /= 'x') &
00649 call define_hist_field(n_albsno,"albsno","%",tstr, tcstr, &
00650 "snow albedo", &
00651 "averaged for coszen>0, weighted by aice", c100, c0, &
00652 ns1, f_albsno)
00653
00654 if (f_albpnd(1:1) /= 'x') &
00655 call define_hist_field(n_albpnd,"albpnd","%",tstr, tcstr, &
00656 "melt pond albedo", &
00657 "averaged for coszen>0, weighted by aice", c100, c0, &
00658 ns1, f_albpnd)
00659
00660 if (f_coszen(1:1) /= 'x') &
00661 call define_hist_field(n_coszen,"coszen","radian",tstr, tcstr, &
00662 "cosine of the zenith angle", &
00663 "negative below horizon", c1, c0, &
00664 ns1, f_coszen)
00665
00666 if (f_flat(1:1) /= 'x') &
00667 call define_hist_field(n_flat,"flat","W/m^2",tstr, tcstr, &
00668 "latent heat flux (cpl)", &
00669 "positive downward", c1, c0, &
00670 ns1, f_flat)
00671
00672 if (f_flat_ai(1:1) /= 'x') &
00673 call define_hist_field(n_flat_ai,"flat_ai","W/m^2",tstr, tcstr, &
00674 "latent heat flux", &
00675 "weighted by ice area", c1, c0, &
00676 ns1, f_flat_ai)
00677
00678 if (f_fsens(1:1) /= 'x') &
00679 call define_hist_field(n_fsens,"fsens","W/m^2",tstr, tcstr, &
00680 "sensible heat flux (cpl)", &
00681 "positive downward", c1, c0, &
00682 ns1, f_fsens)
00683
00684 if (f_fsens_ai(1:1) /= 'x') &
00685 call define_hist_field(n_fsens_ai,"fsens_ai","W/m^2",tstr, tcstr, &
00686 "sensible heat flux", &
00687 "weighted by ice area", c1, c0, &
00688 ns1, f_fsens_ai)
00689
00690 if (f_flwup(1:1) /= 'x') &
00691 call define_hist_field(n_flwup,"flwup","W/m^2",tstr, tcstr, &
00692 "upward longwave flux (cpl)", &
00693 "positive downward", c1, c0, &
00694 ns1, f_flwup)
00695
00696 if (f_flwup_ai(1:1) /= 'x') &
00697 call define_hist_field(n_flwup_ai,"flwup_ai","W/m^2",tstr, tcstr, &
00698 "upward longwave flux", &
00699 "weighted by ice area", c1, c0, &
00700 ns1, f_flwup_ai)
00701
00702 if (f_evap(1:1) /= 'x') &
00703 call define_hist_field(n_evap,"evap","cm/day",tstr, tcstr, &
00704 "evaporative water flux (cpl)", &
00705 "none", mps_to_cmpdy/rhofresh, c0, &
00706 ns1, f_evap)
00707
00708 if (f_evap_ai(1:1) /= 'x') &
00709 call define_hist_field(n_evap_ai,"evap_ai","cm/day",tstr, tcstr, &
00710 "evaporative water flux", &
00711 "weighted by ice area", mps_to_cmpdy/rhofresh, c0, &
00712 ns1, f_evap_ai)
00713
00714 if (f_Tair(1:1) /= 'x') &
00715 call define_hist_field(n_Tair,"Tair","C",tstr, tcstr, &
00716 "air temperature", &
00717 "none", c1, -tffresh, &
00718 ns1, f_Tair)
00719
00720 if (f_Tref(1:1) /= 'x') &
00721 call define_hist_field(n_Tref,"Tref","C",tstr, tcstr, &
00722 "2m reference temperature", &
00723 "none", c1, c0, &
00724 ns1, f_Tref)
00725
00726 if (f_Qref(1:1) /= 'x') &
00727 call define_hist_field(n_Qref,"Qref","g/kg",tstr, tcstr, &
00728 "2m reference specific humidity", &
00729 "none", kg_to_g, c0, &
00730 ns1, f_Qref)
00731
00732 if (f_congel(1:1) /= 'x') &
00733 call define_hist_field(n_congel,"congel","cm/day",tstr, tcstr, &
00734 "congelation ice growth", &
00735 "none", mps_to_cmpdy/dt, c0, &
00736 ns1, f_congel)
00737
00738 if (f_frazil(1:1) /= 'x') &
00739 call define_hist_field(n_frazil,"frazil","cm/day",tstr, tcstr, &
00740 "frazil ice growth", &
00741 "none", mps_to_cmpdy/dt, c0, &
00742 ns1, f_frazil)
00743
00744 if (f_snoice(1:1) /= 'x') &
00745 call define_hist_field(n_snoice,"snoice","cm/day",tstr, tcstr, &
00746 "snow-ice formation", &
00747 "none", mps_to_cmpdy/dt, c0, &
00748 ns1, f_snoice)
00749
00750 if (f_meltt(1:1) /= 'x') &
00751 call define_hist_field(n_meltt,"meltt","cm/day",tstr, tcstr, &
00752 "top ice melt", &
00753 "none", mps_to_cmpdy/dt, c0, &
00754 ns1, f_meltt)
00755
00756 if (f_meltb(1:1) /= 'x') &
00757 call define_hist_field(n_meltb,"meltb","cm/day",tstr, tcstr, &
00758 "basal ice melt", &
00759 "none", mps_to_cmpdy/dt, c0, &
00760 ns1, f_meltb)
00761
00762 if (f_meltl(1:1) /= 'x') &
00763 call define_hist_field(n_meltl,"meltl","cm/day",tstr, tcstr, &
00764 "lateral ice melt", &
00765 "none", mps_to_cmpdy/dt, c0, &
00766 ns1, f_meltl)
00767
00768 if (f_melts(1:1) /= 'x') &
00769 call define_hist_field(n_melts,"melts","cm/day",tstr, tcstr, &
00770 "snow melt", &
00771 "none", mps_to_cmpdy/dt, c0, &
00772 ns1, f_melts)
00773
00774 if (f_fresh(1:1) /= 'x') &
00775 call define_hist_field(n_fresh,"fresh","cm/day",tstr, tcstr, &
00776 "freshwtr flx ice to ocn (cpl)", &
00777 "if positive, ocean gains fresh water", &
00778 mps_to_cmpdy/rhofresh, c0, &
00779 ns1, f_fresh)
00780
00781 if (f_fresh_ai(1:1) /= 'x') &
00782 call define_hist_field(n_fresh_ai,"fresh_ai","cm/day",tstr, tcstr, &
00783 "freshwtr flx ice to ocn", &
00784 "weighted by ice area", mps_to_cmpdy/rhofresh, c0, &
00785 ns1, f_fresh_ai)
00786
00787 if (f_fsalt(1:1) /= 'x') &
00788 call define_hist_field(n_fsalt,"fsalt","kg/m^2/day",tstr, tcstr, &
00789 "salt flux ice to ocn (cpl)", &
00790 "if positive, ocean gains salt", secday, c0, &
00791 ns1, f_fsalt)
00792
00793 if (f_fsalt_ai(1:1) /= 'x') &
00794 call define_hist_field(n_fsalt_ai,"fsalt_ai","kg/m^2/day",tstr, tcstr,&
00795 "salt flux ice to ocean", &
00796 "weighted by ice area", secday, c0, &
00797 ns1, f_fsalt_ai)
00798
00799 if (f_fhocn(1:1) /= 'x') &
00800 call define_hist_field(n_fhocn,"fhocn","W/m^2",tstr, tcstr, &
00801 "heat flux ice to ocn (cpl)", &
00802 "if positive, ocean gains heat", c1, c0, &
00803 ns1, f_fhocn)
00804
00805 if (f_fhocn_ai(1:1) /= 'x') &
00806 call define_hist_field(n_fhocn_ai,"fhocn_ai","W/m^2",tstr, tcstr, &
00807 "heat flux ice to ocean", &
00808 "weighted by ice area", c1, c0, &
00809 ns1, f_fhocn_ai)
00810
00811 if (f_fswthru(1:1) /= 'x') &
00812 call define_hist_field(n_fswthru,"fswthru","W/m^2",tstr, tcstr, &
00813 "SW thru ice to ocean (cpl)", &
00814 "if positive, ocean gains heat", c1, c0, &
00815 ns1, f_fswthru)
00816
00817 if (f_fswthru_ai(1:1) /= 'x') &
00818 call define_hist_field(n_fswthru_ai,"fswthru_ai","W/m^2",tstr, tcstr,&
00819 "SW flux thru ice to ocean", &
00820 "weighted by ice area", c1, c0, &
00821 ns1, f_fswthru_ai)
00822
00823 if (f_strairx(1:1) /= 'x') &
00824 call define_hist_field(n_strairx,"strairx","N/m^2",ustr, ucstr, &
00825 "atm/ice stress (x)", &
00826 "positive is x direction on U grid", c1, c0, &
00827 ns1, f_strairx)
00828
00829 if (f_strairy(1:1) /= 'x') &
00830 call define_hist_field(n_strairy,"strairy","N/m^2",ustr, ucstr, &
00831 "atm/ice stress (y)", &
00832 "positive is y direction on U grid", c1, c0, &
00833 ns1, f_strairy)
00834
00835 if (f_strtltx(1:1) /= 'x') &
00836 call define_hist_field(n_strtltx,"strtltx","N/m^2",ustr, ucstr, &
00837 "sea sfc tilt stress (x)", &
00838 "positive is x direction on U grid", c1, c0, &
00839 ns1, f_strtltx)
00840
00841 if (f_strtlty(1:1) /= 'x') &
00842 call define_hist_field(n_strtlty,"strtlty","N/m^2",ustr, ucstr, &
00843 "sea sfc tilt stress (y)", &
00844 "positive is y direction on U grid", c1, c0, &
00845 ns1, f_strtlty)
00846
00847 if (f_strcorx(1:1) /= 'x') &
00848 call define_hist_field(n_strcorx,"strcorx","N/m^2",ustr, ucstr, &
00849 "coriolis stress (x)", &
00850 "positive is x direction on U grid", c1, c0, &
00851 ns1, f_strcorx)
00852
00853 if (f_strcory(1:1) /= 'x') &
00854 call define_hist_field(n_strcory,"strcory","N/m^2",ustr, ucstr, &
00855 "coriolis stress (y)", &
00856 "positive is y direction on U grid", c1, c0, &
00857 ns1, f_strcory)
00858
00859 if (f_strocnx(1:1) /= 'x') &
00860 call define_hist_field(n_strocnx,"strocnx","N/m^2",ustr, ucstr, &
00861 "ocean/ice stress (x)", &
00862 "positive is x direction on U grid", c1, c0, &
00863 ns1, f_strocnx)
00864
00865 if (f_strocny(1:1) /= 'x') &
00866 call define_hist_field(n_strocny,"strocny","N/m^2",ustr, ucstr, &
00867 "ocean/ice stress (y)", &
00868 "positive is y direction on U grid", c1, c0, &
00869 ns1, f_strocny)
00870
00871 if (f_strintx(1:1) /= 'x') &
00872 call define_hist_field(n_strintx,"strintx","N/m^2",ustr, ucstr, &
00873 "internal ice stress (x)", &
00874 "positive is x direction on U grid", c1, c0, &
00875 ns1, f_strintx)
00876
00877 if (f_strinty(1:1) /= 'x') &
00878 call define_hist_field(n_strinty,"strinty","N/m^2",ustr, ucstr, &
00879 "internal ice stress (y)", &
00880 "positive is y direction on U grid", c1, c0, &
00881 ns1, f_strinty)
00882
00883 if (f_strength(1:1) /= 'x') &
00884 call define_hist_field(n_strength,"strength","N/m",tstr, tcstr, &
00885 "compressive ice strength", &
00886 "none", c1, c0, &
00887 ns1, f_strength)
00888
00889 if (f_opening(1:1) /= 'x') &
00890 call define_hist_field(n_opening,"opening","%/day",tstr, tcstr, &
00891 "lead area opening rate", &
00892 "none", secday*c100, c0, &
00893 ns1, f_opening)
00894
00895 if (f_divu(1:1) /= 'x') &
00896 call define_hist_field(n_divu,"divu","%/day",tstr, tcstr, &
00897 "strain rate (divergence)", &
00898 "none", secday*c100, c0, &
00899 ns1, f_divu)
00900
00901 if (f_shear(1:1) /= 'x') &
00902 call define_hist_field(n_shear,"shear","%/day",tstr, tcstr, &
00903 "strain rate (shear)", &
00904 "none", secday*c100, c0, &
00905 ns1, f_shear)
00906
00907 if (f_sig1(1:1) /= 'x') &
00908 call define_hist_field(n_sig1,"sig1"," ",ustr, ucstr, &
00909 "norm. principal stress 1", &
00910 "sig1 is instantaneous", c1, c0, &
00911 ns1, f_sig1)
00912
00913 if (f_sig2(1:1) /= 'x') &
00914 call define_hist_field(n_sig2,"sig2"," ",ustr, ucstr, &
00915 "norm. principal stress 2", &
00916 "sig2 is instantaneous", c1, c0, &
00917 ns1, f_sig2)
00918
00919 if (f_dvidtt(1:1) /= 'x') &
00920 call define_hist_field(n_dvidtt,"dvidtt","cm/day",tstr, tcstr, &
00921 "volume tendency thermo", &
00922 "none", mps_to_cmpdy, c0, &
00923 ns1, f_dvidtt)
00924
00925 if (f_dvidtd(1:1) /= 'x') &
00926 call define_hist_field(n_dvidtd,"dvidtd","cm/day",tstr, tcstr, &
00927 "volume tendency dynamics", &
00928 "none", mps_to_cmpdy, c0, &
00929 ns1, f_dvidtd)
00930
00931 if (f_daidtt(1:1) /= 'x') &
00932 call define_hist_field(n_daidtt,"daidtt","%/day",tstr, tcstr, &
00933 "area tendency thermo", &
00934 "none", secday*c100, c0, &
00935 ns1, f_daidtt)
00936
00937 if (f_daidtd(1:1) /= 'x') &
00938 call define_hist_field(n_daidtd,"daidtd","%/day",tstr, tcstr, &
00939 "area tendency dynamics", &
00940 "none", secday*c100, c0, &
00941 ns1, f_daidtd)
00942
00943 if (f_mlt_onset(1:1) /= 'x') &
00944 call define_hist_field(n_mlt_onset,"mlt_onset","day of year", &
00945 tstr, tcstr,"melt onset date", &
00946 "midyear restart gives erroneous dates", c1, c0, &
00947 ns1, f_mlt_onset)
00948
00949 if (f_frz_onset(1:1) /= 'x') &
00950 call define_hist_field(n_frz_onset,"frz_onset","day of year", &
00951 tstr, tcstr,"freeze onset date", &
00952 "midyear restart gives erroneous dates", c1, c0, &
00953 ns1, f_frz_onset)
00954
00955 if (f_dardg1dt(1:1) /= 'x') &
00956 call define_hist_field(n_dardg1dt,"dardg1dt","%/day",tstr, tcstr, &
00957 "ice area ridging rate", &
00958 "none", secday*c100, c0, &
00959 ns1, f_dardg1dt)
00960
00961 if (f_dardg2dt(1:1) /= 'x') &
00962 call define_hist_field(n_dardg2dt,"dardg2dt","%/day",tstr, tcstr, &
00963 "ridge area formation rate", &
00964 "none", secday*c100, c0, &
00965 ns1, f_dardg2dt)
00966
00967 if (f_dvirdgdt(1:1) /= 'x') &
00968 call define_hist_field(n_dvirdgdt,"dvirdgdt","cm/day",tstr, tcstr, &
00969 "ice volume ridging rate", &
00970 "none", mps_to_cmpdy, c0, &
00971 ns1, f_dvirdgdt)
00972
00973 if (f_hisnap(1:1) /= 'x') &
00974 call define_hist_field(n_hisnap,"hisnap","m",tstr, tcstr, &
00975 "ice volume snapshot", &
00976 "none", c1, c0, &
00977 ns1, f_hisnap)
00978
00979 if (f_aisnap(1:1) /= 'x') &
00980 call define_hist_field(n_aisnap,"aisnap"," ",tstr, tcstr, &
00981 "ice area snapshot", &
00982 "none", c1, c0, &
00983 ns1, f_aisnap)
00984
00985 if (f_trsig(1:1) /= 'x') &
00986 call define_hist_field(n_trsig,"trsig","N/m^2",tstr, tcstr, &
00987 "internal stress tensor trace", &
00988 "ice strength approximation", c1, c0, &
00989 ns1, f_trsig)
00990
00991 if (f_icepresent(1:1) /= 'x') &
00992 call define_hist_field(n_icepresent,"ice_present","1",tstr, tcstr, &
00993 "fraction of time-avg interval that any ice is present", &
00994 "ice extent flag", c1, c0, &
00995 ns1, f_icepresent)
00996
00997 if (f_fsurf_ai(1:1) /= 'x') &
00998 call define_hist_field(n_fsurf_ai,"fsurf_ai","W/m^2",tstr, tcstr, &
00999 "net surface heat flux", &
01000 "positive downward, excludes conductive flux, weighted by ice area", &
01001 c1, c0, ns1, f_fsurf_ai)
01002
01003 if (f_fcondtop_ai(1:1) /= 'x') &
01004 call define_hist_field(n_fcondtop_ai,"fcondtop_ai","W/m^2", &
01005 tstr, tcstr,"top surface conductive heat flux", &
01006 "positive downwards, weighted by ice area", c1, c0, &
01007 ns1, f_fcondtop_ai)
01008
01009 if (f_fmeltt_ai(1:1) /= 'x') &
01010 call define_hist_field(n_fmeltt_ai,"fmeltt_ai","W/m^2",tstr, tcstr, &
01011 "net surface heat flux causing melt", &
01012 "always >= 0, weighted by ice area", c1, c0, &
01013 ns1, f_fmeltt_ai)
01014
01015
01016 if (f_aicen(1:1) /= 'x') then
01017 do n=1,ncat_hist
01018 write(nchar,'(i3.3)') n
01019 write(vname_in,'(a,a)') 'aicen', trim(nchar)
01020 stmp = 'ice area, category '
01021 write(vdesc_in,'(a,2x,a)') trim(stmp), trim(nchar)
01022 call define_hist_field(n_aicen(n,:),vname_in,"%",tstr, tcstr, &
01023 vdesc_in, "Ice range:", c100, c0, &
01024 ns1, f_aicen)
01025 enddo
01026 endif
01027
01028 if (f_vicen(1:1) /= 'x') then
01029 do n=1,ncat_hist
01030 write(nchar,'(i3.3)') n
01031 write(vname_in,'(a,a)') 'vicen', trim(nchar)
01032 stmp = 'ice volume, category '
01033 write(vdesc_in,'(a,2x,a)') trim(stmp), trim(nchar)
01034 call define_hist_field(n_vicen(n,:),vname_in,"m",tstr, tcstr, &
01035 vdesc_in, "none", c1, c0, &
01036 ns1, f_vicen)
01037 enddo
01038 endif
01039
01040 if (f_fsurfn_ai(1:1) /= 'x') then
01041 do n=1,ncat_hist
01042 write(nchar,'(i3.3)') n
01043 write(vname_in,'(a,a)') 'fsurfn_ai', trim(nchar)
01044 stmp = 'net surface heat flux, category '
01045 write(vdesc_in,'(a,2x,a)') trim(stmp), trim(nchar)
01046 call define_hist_field(n_fsurfn_ai(n,:),vname_in,"W/m^2", &
01047 tstr, tcstr, vdesc_in, "weighted by ice area", c1, c0, &
01048 ns1, f_fsurfn_ai)
01049 enddo
01050 endif
01051
01052 if (f_fcondtopn_ai(1:1) /= 'x') then
01053 do n=1,ncat_hist
01054 write(nchar,'(i3.3)') n
01055 write(vname_in,'(a,a)') 'fcondtopn_ai', trim(nchar)
01056 stmp = 'top sfc conductive heat flux, cat '
01057 write(vdesc_in,'(a,2x,a)') trim(stmp), trim(nchar)
01058 call define_hist_field(n_fcondtopn_ai(n,:),vname_in,"W/m^2", &
01059 tstr, tcstr, vdesc_in, "weighted by ice area", c1, c0, &
01060 ns1, f_fcondtopn_ai)
01061 enddo
01062 endif
01063
01064 if (f_fmelttn_ai(1:1) /= 'x') then
01065 do n=1,ncat_hist
01066 write(nchar,'(i3.3)') n
01067 write(vname_in,'(a,a)') 'fmelttn_ai', trim(nchar)
01068 stmp = 'net sfc heat flux causing melt, cat '
01069 write(vdesc_in,'(a,2x,a)') trim(stmp), trim(nchar)
01070 call define_hist_field(n_fmelttn_ai(n,:),vname_in,"W/m^2", &
01071 tstr, tcstr, vdesc_in, "weighted by ice area", c1, c0, &
01072 ns1, f_fmelttn_ai)
01073 enddo
01074 endif
01075
01076 if (f_flatn_ai(1:1) /= 'x') then
01077 do n=1,ncat_hist
01078 write(nchar,'(i3.3)') n
01079 write(vname_in,'(a,a)') 'flatn_ai', trim(nchar)
01080 stmp = 'latent heat flux, category '
01081 write(vdesc_in,'(a,2x,a)') trim(stmp), trim(nchar)
01082 call define_hist_field(n_flatn_ai(n,:),vname_in,"W/m^2", &
01083 tstr, tcstr, vdesc_in, "weighted by ice area", c1, c0, &
01084 ns1, f_flatn_ai)
01085 enddo
01086 endif
01087
01088
01089
01090
01091 if (f_iage(1:1) /= 'x') then
01092 call define_hist_field(n_iage,"iage","years",tstr, tcstr, &
01093 "sea ice age", &
01094 "none", c1/(secday*days_per_year), c0, &
01095 ns1, f_iage)
01096 endif
01097
01098
01099 if (f_FY(1:1) /= 'x') then
01100 call define_hist_field(n_FY,"FYarea"," ",tstr, tcstr, &
01101 "first-year ice area", &
01102 "weighted by ice area", c1, c0, &
01103 ns1, f_FY)
01104 endif
01105
01106
01107 if (f_aero(1:1) /= 'x') then
01108 do n=1,n_aero
01109 write(nchar,'(i3.3)') n
01110 write(vname_in,'(a,a)') 'aerosnossl', trim(nchar)
01111 call define_hist_field(n_aerosn1(n,:),vname_in,"kg/kg", &
01112 tstr, tcstr,"snow ssl aerosol mass","none", c1, c0, &
01113 ns1, f_aero)
01114 write(vname_in,'(a,a)') 'aerosnoint', trim(nchar)
01115 call define_hist_field(n_aerosn2(n,:),vname_in,"kg/kg", &
01116 tstr, tcstr,"snow int aerosol mass","none", c1, c0, &
01117 ns1, f_aero)
01118 write(vname_in,'(a,a)') 'aeroicessl', trim(nchar)
01119 call define_hist_field(n_aeroic1(n,:),vname_in,"kg/kg", &
01120 tstr, tcstr,"ice ssl aerosol mass","none", c1, c0, &
01121 ns1, f_aero)
01122 write(vname_in,'(a,a)') 'aeroiceint', trim(nchar)
01123 call define_hist_field(n_aeroic2(n,:),vname_in,"kg/kg", &
01124 tstr, tcstr,"ice int aerosol mass","none", c1, c0, &
01125 ns1, f_aero)
01126 enddo
01127 endif
01128
01129 if (f_faero_atm(1:1) /= 'x') then
01130 do n=1,n_aero
01131 write(nchar,'(i3.3)') n
01132 write(vname_in,'(a,a)') 'faero_atm', trim(nchar)
01133 call define_hist_field(n_faero_atm(n,:),vname_in,"kg/m^2 s", &
01134 tstr, tcstr,"aerosol deposition rate","none", c1, c0, &
01135 ns1, f_faero_atm)
01136 enddo
01137 endif
01138
01139 if (f_faero_ocn(1:1) /= 'x') then
01140 do n=1,n_aero
01141 write(nchar,'(i3.3)') n
01142 write(vname_in,'(a,a)') 'faero_ocn', trim(nchar)
01143 call define_hist_field(n_faero_ocn(n,:),vname_in,"kg/m^2 s", &
01144 tstr, tcstr,"aerosol flux to ocean","none", c1, c0, &
01145 ns1, f_faero_ocn)
01146 enddo
01147 endif
01148
01149
01150 if (f_apond(1:1) /= 'x') then
01151 call define_hist_field(n_apond,"apond","%",tstr, tcstr, &
01152 "melt pond concentration", &
01153 "none", c100, c0, &
01154 ns1, f_apond)
01155 endif
01156
01157 if (f_apondn(1:1) /= 'x') then
01158 do n=1,ncat_hist
01159 write(nchar,'(i3.3)') n
01160 write(vname_in,'(a,a)') 'apond', trim(nchar)
01161 stmp = 'melt pond concentration, category '
01162 write(vdesc_in,'(a,2x,a)') trim(stmp), trim(nchar)
01163 call define_hist_field(n_apondn(n,:),vname_in,"%", &
01164 tstr, tcstr, vdesc_in, "none", c100, c0, &
01165 ns1, f_apondn)
01166 enddo
01167 endif
01168
01169 enddo
01170
01171 allocate(aa(nx_block,ny_block,num_avail_hist_fields,max_blocks))
01172
01173
01174
01175
01176
01177 igrd=.true.
01178
01179 igrd(n_tmask ) = f_tmask
01180 igrd(n_tarea ) = f_tarea
01181 igrd(n_uarea ) = f_uarea
01182 igrd(n_dxt ) = f_dxt
01183 igrd(n_dyt ) = f_dyt
01184 igrd(n_dxu ) = f_dxu
01185 igrd(n_dyu ) = f_dyu
01186 igrd(n_HTN ) = f_HTN
01187 igrd(n_HTE ) = f_HTE
01188 igrd(n_ANGLE ) = f_ANGLE
01189 igrd(n_ANGLET ) = f_ANGLET
01190
01191 ntmp(:) = 0
01192 if (my_task == master_task) then
01193 write(nu_diag,*) ' '
01194 write(nu_diag,*) 'The following variables will be ', &
01195 'written to the history tape: '
01196 write(nu_diag,101) 'description','units','variable','frequency','x'
01197 do n=1,num_avail_hist_fields
01198 if (avail_hist_fields(n)%vhistfreq_n /= 0) &
01199 write(nu_diag,100) avail_hist_fields(n)%vdesc, &
01200 avail_hist_fields(n)%vunit, avail_hist_fields(n)%vname, &
01201 avail_hist_fields(n)%vhistfreq,avail_hist_fields(n)%vhistfreq_n
01202 do ns = 1, nstreams
01203 if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) &
01204 ntmp(ns)=ntmp(ns)+1
01205 enddo
01206 enddo
01207 write(nu_diag,*) ' '
01208 endif
01209 100 format (1x,a40,2x,a12,2x,a12,1x,a1,2x,i6)
01210 101 format (2x,a19,10x,a12,9x,a12,2x,a,3x,a1)
01211
01212 call broadcast_array(ntmp, master_task)
01213 do ns = 1, nstreams
01214 if (ntmp(ns)==0) histfreq_n(ns) = 0
01215 enddo
01216
01217
01218
01219
01220 aa(:,:,:,:) = c0
01221 avgct(:) = c0
01222
01223 if (restart .and. yday >= c2) then
01224
01225 mlt_onset = 999._dbl_kind
01226 frz_onset = 999._dbl_kind
01227 else
01228 mlt_onset = c0
01229 frz_onset = c0
01230 endif
01231
01232 end subroutine init_hist
01233
01234
01235
01236
01237
01238
01239
01240
01241
01242 subroutine ice_write_hist (dt)
01243
01244
01245
01246
01247
01248
01249
01250
01251
01252
01253
01254 use ice_blocks
01255 use ice_domain
01256 use ice_calendar, only: new_year, secday, yday, write_history, &
01257 write_ic, time, histfreq, nstreams
01258 use ice_state
01259 use ice_constants
01260 use ice_flux
01261 use ice_dyn_evp
01262 use ice_itd, only: ilyr1, slyr1
01263 use ice_timers
01264
01265
01266
01267 real (kind=dbl_kind), intent(in) ::
01268 dt
01269
01270
01271 integer (kind=int_kind) ::
01272 i,j,k,n,nct,ns ,
01273 iblk ,
01274 ilo,ihi,jlo,jhi ,
01275 nstrm
01276
01277 real (kind=dbl_kind) ::
01278 ravgct ,
01279 ravgctz
01280
01281 type (block) ::
01282 this_block
01283
01284 real (kind=dbl_kind) ::
01285 worka(nx_block,ny_block),
01286 workb(nx_block,ny_block)
01287
01288 real (kind=dbl_kind) ::
01289 ai ,
01290 ain ,
01291 hs ,
01292 qs ,
01293 qi ,
01294 uee, vnn ,
01295 hiee, hinn
01296
01297
01298
01299
01300
01301 do ns=1,nstreams
01302 if (.not. hist_avg .or. histfreq(ns) == '1') then
01303 do k=1,num_avail_hist_fields
01304 if (avail_hist_fields(k)%vhistfreq == histfreq(ns)) &
01305 aa(:,:,k,:) = c0
01306 enddo
01307 avgct(ns) = c1
01308 else
01309 avgct(ns) = avgct(ns) + c1
01310 if (avgct(ns) == c1) time_beg(ns) = (time-dt)/int(secday)
01311 endif
01312 enddo
01313
01314
01315
01316
01317
01318
01319 do iblk = 1, nblocks
01320
01321 workb(:,:) = aice_init(:,:,iblk)
01322
01323 this_block = get_block(blocks_ice(iblk),iblk)
01324 ilo = this_block%ilo
01325 ihi = this_block%ihi
01326 jlo = this_block%jlo
01327 jhi = this_block%jhi
01328
01329 if (f_hi(1:1) /= 'x') &
01330 call accum_hist_field(n_hi, iblk, vice(:,:,iblk))
01331 if (f_hs(1:1) /= 'x') &
01332 call accum_hist_field(n_hs, iblk, vsno(:,:,iblk))
01333 if (f_Tsfc(1:1) /= 'x') &
01334 call accum_hist_field(n_Tsfc, iblk, trcr(:,:,nt_Tsfc,iblk))
01335 if (f_aice(1:1) /= 'x') &
01336 call accum_hist_field(n_aice, iblk, aice(:,:,iblk))
01337 if (f_uvel(1:1) /= 'x') &
01338 call accum_hist_field(n_uvel, iblk, uvel(:,:,iblk))
01339 if (f_vvel(1:1) /= 'x') &
01340 call accum_hist_field(n_vvel, iblk, vvel(:,:,iblk))
01341 if (f_transix(1:1) /= 'x') then
01342 worka(:,:) = c0
01343 do j = jlo, jhi
01344 do i = ilo, ihi
01345 uee = p5*(uvel(i,j,iblk)+uvel(i,j-1,iblk))
01346 hiee = p5*(vice(i,j,iblk)+vice(i+1,j,iblk))
01347 worka(i,j) = uee*hiee*HTE(i,j,iblk)*rhoi
01348 enddo
01349 enddo
01350 call accum_hist_field(n_transix, iblk, worka(:,:))
01351 endif
01352 if (f_transiy(1:1) /= 'x') then
01353 worka(:,:) = c0
01354 do j = jlo, jhi
01355 do i = ilo, ihi
01356 vnn = p5*(vvel(i,j,iblk)+vvel(i-1,j,iblk))
01357 hinn = p5*(vice(i,j,iblk)+vice(i,j+1,iblk))
01358 worka(i,j) = vnn*hinn*HTN(i,j,iblk)*rhoi
01359 enddo
01360 enddo
01361 call accum_hist_field(n_transiy, iblk, worka(:,:))
01362 endif
01363
01364 if (f_fswdn(1:1) /= 'x') &
01365 call accum_hist_field(n_fswdn, iblk, fsw(:,:,iblk))
01366 if (f_fswup(1:1) /= 'x') &
01367 call accum_hist_field(n_fswup, iblk, &
01368 (fsw(:,:,iblk)-fswabs(:,:,iblk)*workb(:,:)))
01369 if (f_flwdn(1:1) /= 'x') &
01370 call accum_hist_field(n_flwdn, iblk, flw(:,:,iblk))
01371 if (f_snow(1:1) /= 'x') &
01372 call accum_hist_field(n_snow, iblk, fsnow(:,:,iblk))
01373 if (f_snow_ai(1:1) /= 'x') &
01374 call accum_hist_field(n_snow_ai, iblk, fsnow(:,:,iblk)*workb(:,:))
01375 if (f_rain(1:1) /= 'x') &
01376 call accum_hist_field(n_rain, iblk, frain(:,:,iblk))
01377 if (f_rain_ai(1:1) /= 'x') &
01378 call accum_hist_field(n_rain_ai, iblk, frain(:,:,iblk)*workb(:,:))
01379
01380 if (f_sst(1:1) /= 'x') &
01381 call accum_hist_field(n_sst, iblk, sst(:,:,iblk))
01382 if (f_sss(1:1) /= 'x') &
01383 call accum_hist_field(n_sss, iblk, sss(:,:,iblk))
01384 if (f_uocn(1:1) /= 'x') &
01385 call accum_hist_field(n_uocn, iblk, uocn(:,:,iblk))
01386 if (f_vocn(1:1) /= 'x') &
01387 call accum_hist_field(n_vocn, iblk, vocn(:,:,iblk))
01388 if (f_frzmlt(1:1) /= 'x') &
01389 call accum_hist_field(n_frzmlt, iblk, frzmlt(:,:,iblk))
01390
01391 if (f_fswfac(1:1) /= 'x') &
01392 call accum_hist_field(n_fswfac, iblk, fswfac(:,:,iblk))
01393 if (f_fswabs(1:1) /= 'x') &
01394 call accum_hist_field(n_fswabs, iblk, fswabs(:,:,iblk))
01395 if (f_fswabs_ai(1:1) /= 'x') &
01396 call accum_hist_field(n_fswabs_ai,iblk,fswabs(:,:,iblk)*workb(:,:))
01397
01398 if (f_alvdr(1:1) /= 'x') &
01399 call accum_hist_field(n_alvdr, iblk, alvdr(:,:,iblk))
01400 if (f_alidr(1:1) /= 'x') &
01401 call accum_hist_field(n_alidr, iblk, alidr(:,:,iblk))
01402
01403 if (f_albice (1:1) /= 'x') &
01404 call accum_hist_field(n_albice, iblk, albice(:,:,iblk))
01405 if (f_albsno (1:1) /= 'x') &
01406 call accum_hist_field(n_albsno, iblk, albsno(:,:,iblk))
01407 if (f_albpnd (1:1) /= 'x') &
01408 call accum_hist_field(n_albpnd, iblk, albpnd(:,:,iblk))
01409 if (f_coszen (1:1) /= 'x') &
01410 call accum_hist_field(n_coszen, iblk, coszen(:,:,iblk))
01411
01412 if (f_flat(1:1) /= 'x') &
01413 call accum_hist_field(n_flat, iblk, flat(:,:,iblk))
01414 if (f_flat_ai(1:1) /= 'x') &
01415 call accum_hist_field(n_flat_ai,iblk, flat(:,:,iblk)*workb(:,:))
01416 if (f_fsens(1:1) /= 'x') &
01417 call accum_hist_field(n_fsens, iblk, fsens(:,:,iblk))
01418 if (f_fsens_ai(1:1) /= 'x') &
01419 call accum_hist_field(n_fsens_ai,iblk, fsens(:,:,iblk)*workb(:,:))
01420 if (f_flwup(1:1) /= 'x') &
01421 call accum_hist_field(n_flwup, iblk, flwout(:,:,iblk))
01422 if (f_flwup_ai(1:1) /= 'x') &
01423 call accum_hist_field(n_flwup_ai,iblk, flwout(:,:,iblk)*workb(:,:))
01424 if (f_evap(1:1) /= 'x') &
01425 call accum_hist_field(n_evap, iblk, evap(:,:,iblk))
01426 if (f_evap_ai(1:1) /= 'x') &
01427 call accum_hist_field(n_evap_ai,iblk, evap(:,:,iblk)*workb(:,:))
01428
01429 if (f_Tair(1:1) /= 'x') &
01430 call accum_hist_field(n_Tair, iblk, Tair(:,:,iblk))
01431 if (f_Tref(1:1) /= 'x') &
01432 call accum_hist_field(n_Tref, iblk, Tref(:,:,iblk))
01433 if (f_Qref(1:1) /= 'x') &
01434 call accum_hist_field(n_Qref, iblk, Qref(:,:,iblk))
01435 if (f_congel(1:1) /= 'x') &
01436 call accum_hist_field(n_congel, iblk, congel(:,:,iblk))
01437 if (f_frazil(1:1) /= 'x') &
01438 call accum_hist_field(n_frazil, iblk, frazil(:,:,iblk))
01439 if (f_snoice(1:1) /= 'x') &
01440 call accum_hist_field(n_snoice, iblk, snoice(:,:,iblk))
01441 if (f_meltt(1:1) /= 'x') &
01442 call accum_hist_field(n_meltt, iblk, meltt(:,:,iblk))
01443 if (f_meltb(1:1) /= 'x') &
01444 call accum_hist_field(n_meltb, iblk, meltb(:,:,iblk))
01445 if (f_meltl(1:1) /= 'x') &
01446 call accum_hist_field(n_meltl, iblk, meltl(:,:,iblk))
01447 if (f_melts(1:1) /= 'x') &
01448 call accum_hist_field(n_melts, iblk, melts(:,:,iblk))
01449
01450 if (f_fresh(1:1) /= 'x') &
01451 call accum_hist_field(n_fresh, iblk, fresh(:,:,iblk))
01452 if (f_fresh_ai(1:1) /= 'x') &
01453 call accum_hist_field(n_fresh_ai,iblk, fresh_gbm(:,:,iblk))
01454 if (f_fsalt(1:1) /= 'x') &
01455 call accum_hist_field(n_fsalt, iblk, fsalt(:,:,iblk))
01456 if (f_fsalt_ai(1:1) /= 'x') &
01457 call accum_hist_field(n_fsalt_ai,iblk, fsalt_gbm(:,:,iblk))
01458 if (f_fhocn(1:1) /= 'x') &
01459 call accum_hist_field(n_fhocn, iblk, fhocn(:,:,iblk))
01460 if (f_fhocn_ai(1:1) /= 'x') &
01461 call accum_hist_field(n_fhocn_ai,iblk, fhocn_gbm(:,:,iblk))
01462 if (f_fswthru(1:1) /= 'x') &
01463 call accum_hist_field(n_fswthru, iblk, fswthru(:,:,iblk))
01464 if (f_fswthru_ai(1:1) /= 'x') &
01465 call accum_hist_field(n_fswthru_ai,iblk, fswthru_gbm(:,:,iblk))
01466
01467 if (f_strairx(1:1) /= 'x') &
01468 call accum_hist_field(n_strairx, iblk, strairx(:,:,iblk))
01469 if (f_strairy(1:1) /= 'x') &
01470 call accum_hist_field(n_strairy, iblk, strairy(:,:,iblk))
01471 if (f_strtltx(1:1) /= 'x') &
01472 call accum_hist_field(n_strtltx, iblk, strtltx(:,:,iblk))
01473 if (f_strtlty(1:1) /= 'x') &
01474 call accum_hist_field(n_strtlty, iblk, strtlty(:,:,iblk))
01475 if (f_strcorx(1:1) /= 'x') &
01476 call accum_hist_field(n_strcorx, iblk, fm(:,:,iblk)*vvel(:,:,iblk))
01477 if (f_strcory(1:1) /= 'x') &
01478 call accum_hist_field(n_strcory, iblk,-fm(:,:,iblk)*uvel(:,:,iblk))
01479 if (f_strocnx(1:1) /= 'x') &
01480 call accum_hist_field(n_strocnx, iblk, strocnx(:,:,iblk))
01481 if (f_strocny(1:1) /= 'x') &
01482 call accum_hist_field(n_strocny, iblk, strocny(:,:,iblk))
01483 if (f_strintx(1:1) /= 'x') &
01484 call accum_hist_field(n_strintx, iblk, strintx(:,:,iblk))
01485 if (f_strinty(1:1) /= 'x') &
01486 call accum_hist_field(n_strinty, iblk, strinty(:,:,iblk))
01487 if (f_strength(1:1) /= 'x') &
01488 call accum_hist_field(n_strength, iblk, strength(:,:,iblk))
01489
01490
01491
01492
01493
01494
01495 if (f_divu(1:1) /= 'x') &
01496 call accum_hist_field(n_divu, iblk, divu(:,:,iblk))
01497 if (f_shear(1:1) /= 'x') &
01498 call accum_hist_field(n_shear, iblk, shear(:,:,iblk))
01499
01500
01501
01502
01503
01504 if (f_dvidtt(1:1) /= 'x') &
01505 call accum_hist_field(n_dvidtt, iblk, dvidtt(:,:,iblk))
01506 if (f_dvidtd(1:1) /= 'x') &
01507 call accum_hist_field(n_dvidtd, iblk, dvidtd(:,:,iblk))
01508 if (f_daidtt(1:1) /= 'x') &
01509 call accum_hist_field(n_daidtt, iblk, daidtt(:,:,iblk))
01510 if (f_daidtd(1:1) /= 'x') &
01511 call accum_hist_field(n_daidtd, iblk, daidtd(:,:,iblk))
01512
01513 if (f_opening(1:1) /= 'x') &
01514 call accum_hist_field(n_opening, iblk, opening(:,:,iblk))
01515 if (f_dardg1dt(1:1) /= 'x') &
01516 call accum_hist_field(n_dardg1dt, iblk, dardg1dt(:,:,iblk))
01517 if (f_dardg2dt(1:1) /= 'x') &
01518 call accum_hist_field(n_dardg2dt, iblk, dardg2dt(:,:,iblk))
01519 if (f_dvirdgdt(1:1) /= 'x') &
01520 call accum_hist_field(n_dvirdgdt, iblk, dvirdgdt(:,:,iblk))
01521
01522 if (f_fsurf_ai(1:1) /= 'x') &
01523 call accum_hist_field(n_fsurf_ai, iblk, fsurf(:,:,iblk)*workb(:,:))
01524 if (f_fcondtop_ai(1:1) /= 'x') &
01525 call accum_hist_field(n_fcondtop_ai, iblk, &
01526 fcondtop(:,:,iblk)*workb(:,:))
01527
01528 if (f_icepresent(1:1) /= 'x') then
01529 worka(:,:) = c0
01530 do j = jlo, jhi
01531 do i = ilo, ihi
01532 if (aice(i,j,iblk) > puny) worka(i,j) = c1
01533 enddo
01534 enddo
01535 call accum_hist_field(n_icepresent, iblk, worka(:,:))
01536 endif
01537
01538 nct = min(ncat, ncat_hist)
01539 do n=1,nct
01540 workb(:,:) = aicen_init(:,:,n,iblk)
01541 if (f_aicen(1:1) /= 'x') &
01542 call accum_hist_field(n_aicen(n,:), iblk, aicen(:,:,n,iblk))
01543 if (f_vicen(1:1) /= 'x') &
01544 call accum_hist_field(n_vicen(n,:), iblk, vicen(:,:,n,iblk))
01545 if (f_apondn(1:1) /= 'x') &
01546 call accum_hist_field(n_apondn(n,:), iblk, apondn(:,:,n,iblk))
01547 if (f_fsurfn_ai(1:1) /= 'x') &
01548 call accum_hist_field(n_fsurfn_ai(n,:), iblk, &
01549 fsurfn(:,:,n,iblk)*workb(:,:))
01550 if (f_fcondtopn_ai(1:1) /= 'x') &
01551 call accum_hist_field(n_fcondtopn_ai(n,:), iblk, &
01552 fcondtopn(:,:,n,iblk)*workb(:,:))
01553 if (f_flatn_ai(1:1) /= 'x') &
01554 call accum_hist_field(n_flatn_ai(n,:), iblk, &
01555 flatn(:,:,n,iblk)*workb(:,:))
01556
01557
01558
01559 if (f_fmelttn_ai(1:1) /= 'x') &
01560 call accum_hist_field(n_fmelttn_ai(n,:), iblk, &
01561 max(fsurfn(:,:,n,iblk) - fcondtopn(:,:,n,iblk),c0)*workb(:,:))
01562 enddo
01563 if (f_fs(1:1) /= 'x') then
01564 worka(:,:) = c0
01565 nct = min(ncat, ncat_hist)
01566 do n=1,nct
01567 workb(:,:) = aicen_init(:,:,n,iblk)
01568 do j = jlo, jhi
01569 do i = ilo, ihi
01570 hs = c0
01571 if (workb(i,j) > puny) &
01572 hs = vsnon(i,j,n,iblk)/workb(i,j)
01573 if (hs >= hsmin) &
01574 worka(i,j) = worka(i,j) + workb(i,j)*min( hs/hs0, c1 )
01575 enddo
01576 enddo
01577 enddo
01578 call accum_hist_field(n_fs, iblk, worka(:,:))
01579 endif
01580
01581 if (f_qi(1:1) /= 'x') then
01582 worka(:,:) = c0
01583 nct = min(ncat, ncat_hist)
01584 do n=1,nct
01585 workb(:,:) = aicen_init(:,:,n,iblk)
01586 do k=1,nilyr
01587 do j = jlo, jhi
01588 do i = ilo, ihi
01589 qi = c0
01590 if (vicen(i,j,n,iblk) > puny) &
01591 qi = eicen(i,j,ilyr1(n)+k-1,iblk)*tarea(i,j,iblk)
01592 worka(i,j) = worka(i,j) + workb(i,j)*qi
01593 enddo
01594 enddo
01595 enddo
01596 enddo
01597 call accum_hist_field(n_qi, iblk, worka(:,:))
01598 endif
01599
01600 if (f_qs(1:1) /= 'x') then
01601 worka(:,:) = c0
01602 nct = min(ncat, ncat_hist)
01603 do n=1,nct
01604 workb(:,:) = aicen_init(:,:,n,iblk)
01605 do k=1,nslyr
01606 do j = jlo, jhi
01607 do i = ilo, ihi
01608 qs = c0
01609 if (vsnon(i,j,n,iblk) > puny) &
01610 qs = esnon(i,j,slyr1(n)+k-1,iblk)*tarea(i,j,iblk)
01611 worka(i,j) = worka(i,j) + workb(i,j)*qs
01612 enddo
01613 enddo
01614 enddo
01615 enddo
01616 call accum_hist_field(n_qs, iblk, worka(:,:))
01617 endif
01618
01619
01620 if (f_apond(1:1) /= 'x') then
01621 do ns=1, nstreams
01622 if (n_apond(ns) /= 0) then
01623 worka(:,:) = c0
01624 do j = jlo, jhi
01625 do i = ilo, ihi
01626 if (tmask(i,j,iblk)) then
01627 do n=1,nct
01628 worka(i,j) = worka(i,j) + aa(i,j,n_apondn(n,ns),iblk)
01629 enddo
01630 endif
01631 enddo
01632 enddo
01633 aa(:,:,n_apond(ns),iblk) = worka(:,:)
01634 endif
01635 enddo
01636 endif
01637
01638 if (f_fmeltt_ai(1:1) /= 'x') then
01639 do ns=1, nstreams
01640 if (n_fmeltt_ai(ns) /= 0) then
01641 worka(:,:) = c0
01642 do j = jlo, jhi
01643 do i = ilo, ihi
01644 if (tmask(i,j,iblk)) then
01645 do n=1,nct
01646 worka(i,j) = worka(i,j) + aa(i,j,n_fmelttn_ai(n,ns),iblk)
01647 enddo
01648 endif
01649 enddo
01650 enddo
01651 aa(:,:,n_fmeltt_ai(ns),iblk) = worka(:,:)
01652 endif
01653 enddo
01654 endif
01655
01656
01657
01658 if (f_faero_atm(1:1) /= 'x') then
01659 do n=1,n_aero
01660 call accum_hist_field(n_faero_atm(n,:), iblk, &
01661 faero(:,:,n,iblk))
01662 enddo
01663 endif
01664
01665 if (f_faero_ocn(1:1) /= 'x') then
01666 do n=1,n_aero
01667 call accum_hist_field(n_faero_ocn(n,:), iblk, &
01668 fsoot(:,:,n,iblk))
01669 enddo
01670 endif
01671
01672 if (f_aero(1:1) /= 'x') then
01673 do n=1,n_aero
01674 call accum_hist_field(n_aerosn1(n,:), iblk, &
01675 trcr(:,:,nt_aero +4*(n-1),iblk)/rhos)
01676 call accum_hist_field(n_aerosn2(n,:), iblk, &
01677 trcr(:,:,nt_aero+1+4*(n-1),iblk)/rhos)
01678 call accum_hist_field(n_aeroic1(n,:), iblk, &
01679 trcr(:,:,nt_aero+2+4*(n-1),iblk)/rhoi)
01680 call accum_hist_field(n_aeroic2(n,:), iblk, &
01681 trcr(:,:,nt_aero+3+4*(n-1),iblk)/rhoi)
01682 enddo
01683 endif
01684
01685 enddo
01686
01687
01688
01689
01690
01691
01692 nstrm = nstreams
01693 if (write_ic) nstrm = 1
01694
01695 do ns = 1, nstrm
01696
01697 if (write_history(ns) .or. write_ic) then
01698
01699
01700
01701
01702
01703 ravgct = c1/avgct(ns)
01704
01705 do iblk = 1, nblocks
01706 this_block = get_block(blocks_ice(iblk),iblk)
01707 ilo = this_block%ilo
01708 ihi = this_block%ihi
01709 jlo = this_block%jlo
01710 jhi = this_block%jhi
01711
01712 do k = 1, num_avail_hist_fields
01713
01714 if (avail_hist_fields(k)%vhistfreq == histfreq(ns)) then
01715
01716 do j = jlo, jhi
01717 do i = ilo, ihi
01718 if (.not. tmask(i,j,iblk)) then
01719 aa(i,j,k,iblk) = spval
01720 else
01721 aa(i,j,k,iblk) &
01722 = avail_hist_fields(k)%cona*aa(i,j,k,iblk) &
01723 * ravgct + avail_hist_fields(k)%conb
01724 endif
01725 enddo
01726 enddo
01727
01728
01729 if (avail_hist_fields(k)%vname == 'albice') then
01730 do j = jlo, jhi
01731 do i = ilo, ihi
01732 if (tmask(i,j,iblk)) then
01733 ravgctz = c0
01734 if (albcnt(i,j,iblk,ns) > puny) &
01735 ravgctz = c1/albcnt(i,j,iblk,ns)
01736 if (n_albice(ns) /= 0) aa(i,j,n_albice(ns),iblk) = &
01737 aa(i,j,n_albice(ns),iblk)*avgct(ns)*ravgctz
01738 if (n_albsno(ns) /= 0) aa(i,j,n_albsno(ns),iblk) = &
01739 aa(i,j,n_albsno(ns),iblk)*avgct(ns)*ravgctz
01740 if (n_albpnd(ns) /= 0) aa(i,j,n_albpnd(ns),iblk) = &
01741 aa(i,j,n_albpnd(ns),iblk)*avgct(ns)*ravgctz
01742 endif
01743 enddo
01744 enddo
01745 endif
01746
01747 endif
01748 enddo
01749
01750
01751
01752
01753
01754
01755
01756 call principal_stress (nx_block, ny_block, &
01757 stressp_1 (:,:,iblk), &
01758 stressm_1 (:,:,iblk), &
01759 stress12_1(:,:,iblk), &
01760 prs_sig (:,:,iblk), &
01761 sig1 (:,:,iblk), &
01762 sig2 (:,:,iblk))
01763
01764 do j = jlo, jhi
01765 do i = ilo, ihi
01766 if (.not. tmask(i,j,iblk)) then
01767
01768
01769 if (n_sig1 (ns) /= 0) aa(i,j,n_sig1(ns),iblk ) = spval
01770 if (n_sig2 (ns) /= 0) aa(i,j,n_sig2(ns),iblk ) = spval
01771 if (n_mlt_onset(ns) /= 0) aa(i,j,n_mlt_onset(ns),iblk) = spval
01772 if (n_frz_onset(ns) /= 0) aa(i,j,n_frz_onset(ns),iblk) = spval
01773 if (n_hisnap (ns) /= 0) aa(i,j,n_hisnap(ns),iblk) = spval
01774 if (n_aisnap (ns) /= 0) aa(i,j,n_aisnap(ns),iblk) = spval
01775 if (n_trsig (ns) /= 0) aa(i,j,n_trsig(ns),iblk ) = spval
01776 if (n_iage (ns) /= 0) aa(i,j,n_iage(ns),iblk ) = spval
01777 if (n_FY (ns) /= 0) aa(i,j,n_FY(ns),iblk ) = spval
01778 else
01779
01780
01781
01782
01783 if (n_sig1 (ns) /= 0) aa(i,j,n_sig1(ns),iblk) = &
01784 sig1 (i,j,iblk)*avail_hist_fields(n_sig1(ns))%cona
01785 if (n_sig2 (ns) /= 0) aa(i,j,n_sig2(ns),iblk) = &
01786 sig2 (i,j,iblk)*avail_hist_fields(n_sig2(ns))%cona
01787 if (n_mlt_onset(ns) /= 0) aa(i,j,n_mlt_onset(ns),iblk) = &
01788 mlt_onset(i,j,iblk)
01789 if (n_frz_onset(ns) /= 0) aa(i,j,n_frz_onset(ns),iblk) = &
01790 frz_onset(i,j,iblk)
01791 if (n_hisnap (ns) /= 0) aa(i,j,n_hisnap(ns),iblk) = &
01792 vice(i,j,iblk)
01793 if (n_aisnap (ns) /= 0) aa(i,j,n_aisnap(ns),iblk) = &
01794 aice(i,j,iblk)
01795 if (n_trsig (ns) /= 0) aa(i,j,n_trsig(ns),iblk ) = &
01796 p25*(stressp_1(i,j,iblk) &
01797 + stressp_2(i,j,iblk) &
01798 + stressp_3(i,j,iblk) &
01799 + stressp_4(i,j,iblk))
01800 if (n_iage (ns) /= 0) aa(i,j,n_iage(ns),iblk) = &
01801 trcr(i,j,nt_iage,iblk)*avail_hist_fields(n_iage(ns))%cona
01802 if (n_FY (ns) /= 0) aa(i,j,n_FY(ns),iblk) = &
01803 trcr(i,j,nt_FY,iblk)*avail_hist_fields(n_FY(ns))%cona
01804 endif
01805 enddo
01806 enddo
01807
01808 enddo
01809
01810
01811 time_end(ns) = time/int(secday)
01812
01813
01814
01815
01816
01817 call ice_timer_start(timer_readwrite)
01818
01819 if (trim(history_format) == 'nc') then
01820 call icecdf(ns)
01821 else
01822 call icebin(ns)
01823 endif
01824
01825 call ice_timer_stop(timer_readwrite)
01826
01827
01828
01829
01830 if (write_ic) then
01831 aa(:,:,:,:) = c0
01832 avgct(:) = c0
01833 write_ic = .false.
01834 else
01835 avgct(ns) = c0
01836 endif
01837
01838 do k=1,num_avail_hist_fields
01839 if (avail_hist_fields(k)%vhistfreq == histfreq(ns)) aa(:,:,k,:) = c0
01840 enddo
01841
01842 endif
01843 enddo
01844
01845
01846 do iblk = 1, nblocks
01847 this_block = get_block(blocks_ice(iblk),iblk)
01848 ilo = this_block%ilo
01849 ihi = this_block%ihi
01850 jlo = this_block%jlo
01851 jhi = this_block%jhi
01852
01853 if (new_year) then
01854
01855 do j=jlo,jhi
01856 do i=ilo,ihi
01857
01858 if (lmask_n(i,j,iblk)) mlt_onset(i,j,iblk) = c0
01859
01860 if (lmask_s(i,j,iblk)) frz_onset(i,j,iblk) = c0
01861 enddo
01862 enddo
01863 endif
01864
01865 if ((yday >= 181._dbl_kind) .and. &
01866 (yday < 181._dbl_kind+dt/secday)) then
01867
01868 do j=jlo,jhi
01869 do i=ilo,ihi
01870
01871
01872 if (lmask_s(i,j,iblk)) mlt_onset(i,j,iblk) = c0
01873
01874
01875 if (lmask_n(i,j,iblk)) frz_onset(i,j,iblk) = c0
01876 enddo
01877 enddo
01878
01879 endif
01880 enddo
01881
01882 write_ic = .false.
01883
01884 end subroutine ice_write_hist
01885
01886
01887
01888 subroutine define_hist_field(id, vname, vunit, vcoord, vcellmeas, &
01889 vdesc, vcomment, cona, conb, &
01890 ns1, vhistfreq)
01891
01892
01893
01894
01895
01896
01897
01898
01899
01900
01901 use ice_exit
01902 use ice_calendar, only: histfreq, histfreq_n, nstreams
01903
01904
01905
01906 integer (int_kind), dimension(max_nstrm), intent(out) ::
01907 id
01908
01909
01910
01911
01912 character (len=*), intent(in) ::
01913 vname ,
01914 vunit ,
01915 vcoord ,
01916 vcellmeas ,
01917 vdesc ,
01918 vcomment
01919
01920 real (kind=dbl_kind), intent(in) ::
01921 cona ,
01922 conb
01923
01924 character (len=*), intent(in) ::
01925 vhistfreq
01926
01927 integer (kind=int_kind), intent(in) ::
01928 ns1
01929
01930 integer (kind=int_kind) ::
01931 ns ,
01932 lenf
01933
01934 character (len=40) :: stmp
01935
01936 lenf = len(trim(vhistfreq))
01937 if (ns1 == 1) id(:) = 0
01938
01939 do ns = 1, nstreams
01940 if (vhistfreq(ns1:ns1) == histfreq(ns)) then
01941
01942 num_avail_hist_fields = num_avail_hist_fields + 1
01943
01944 if (num_avail_hist_fields > max_avail_hist_fields) &
01945 call abort_ice("Need to increase max_avail_hist_fields")
01946
01947 id(ns) = num_avail_hist_fields
01948
01949 stmp = vname
01950 if (lenf > 1 .and. ns1 > 1) &
01951 write(stmp,'(a,a1,a1)') trim(stmp),'_',vhistfreq(ns1:ns1)
01952
01953 avail_hist_fields(id(ns))%vname = trim(stmp)
01954 avail_hist_fields(id(ns))%vunit = trim(vunit)
01955 avail_hist_fields(id(ns))%vcoord = trim(vcoord)
01956 avail_hist_fields(id(ns))%vcellmeas = trim(vcellmeas)
01957 avail_hist_fields(id(ns))%vdesc = trim(vdesc)
01958 avail_hist_fields(id(ns))%vcomment = trim(vcomment)
01959 avail_hist_fields(id(ns))%cona = cona
01960 avail_hist_fields(id(ns))%conb = conb
01961 avail_hist_fields(id(ns))%vhistfreq = vhistfreq(ns1:ns1)
01962 avail_hist_fields(id(ns))%vhistfreq_n = histfreq_n(ns)
01963
01964 endif
01965 enddo
01966
01967 end subroutine define_hist_field
01968
01969
01970
01971 subroutine accum_hist_field(id, iblk, field_accum)
01972
01973
01974
01975
01976
01977
01978
01979 use ice_domain
01980 use ice_grid, only: tmask
01981 use ice_calendar, only: nstreams
01982
01983
01984
01985 integer (int_kind), dimension(max_nstrm), intent(in) ::
01986 id
01987
01988
01989 integer (kind=int_kind), intent(in) :: iblk
01990
01991 real (kind=dbl_kind), intent(in) :: field_accum(nx_block,ny_block)
01992
01993 type (block) ::
01994 this_block
01995
01996 integer (kind=int_kind) :: i,j, ilo, ihi, jlo, jhi, ns, idns
01997
01998
01999
02000
02001
02002 do ns = 1, nstreams
02003 idns = id(ns)
02004 if (idns > 0) then
02005
02006 this_block = get_block(blocks_ice(iblk),iblk)
02007 ilo = this_block%ilo
02008 ihi = this_block%ihi
02009 jlo = this_block%jlo
02010 jhi = this_block%jhi
02011
02012 do j = jlo, jhi
02013 do i = ilo, ihi
02014 if (tmask(i,j,iblk)) then
02015 aa(i,j,idns, iblk) = aa(i,j,idns, iblk) + field_accum(i,j)
02016 endif
02017 enddo
02018 enddo
02019
02020 endif
02021 enddo
02022
02023 end subroutine accum_hist_field
02024
02025
02026
02027 end module ice_history
02028
02029