00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 module ice_history_write
00012
00013
00014
00015 use ice_kinds_mod
00016 use ice_broadcast
00017 use ice_communicate, only: my_task, master_task, MPI_COMM_ICE
00018 use ice_blocks
00019 use ice_grid
00020 use ice_fileunits
00021 use ice_history_fields
00022
00023
00024
00025 implicit none
00026 save
00027
00028
00029
00030 contains
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040 subroutine icecdf(ns)
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052 use shr_mpi_mod
00053
00054 use ice_gather_scatter
00055 use ice_domain_size
00056 use ice_constants
00057 use ice_calendar, only: time, sec, idate, idate0, nyr, month, &
00058 mday, write_ic, histfreq, histfreq_n, &
00059 year_init, new_year, new_month, new_day, &
00060 dayyr, daymo, days_per_year
00061 use ice_work, only: work_g1, work_gr, work_gr3
00062 use ice_restart, only: lenstr, runid, lcdf64
00063 use ice_domain, only: distrb_info
00064 use ice_itd, only: c_hi_range
00065 use ice_exit
00066 use ice_pio
00067 use pio
00068 use shr_sys_mod, only : shr_sys_flush
00069
00070
00071
00072
00073
00074 integer (kind=int_kind), intent(in) :: ns
00075
00076 integer (kind=int_kind) ::
00077 iblk,ilo,ihi,jlo,jhi,lon,lat
00078 integer (kind=int_kind) :: i,j,n,k,
00079 status,imtid,jmtid,timid,
00080 length,nvertexid,ivertex
00081 integer (kind=int_kind), dimension(2) :: dimid2
00082 integer (kind=int_kind), dimension(3) :: dimid3
00083 integer (kind=int_kind), dimension(3) :: dimid_nverts
00084 real (kind=real_kind) :: ltime
00085 character (char_len) :: title
00086 character (char_len_long) :: ncfile(max_nstrm)
00087
00088 integer (kind=int_kind) :: iyear, imonth, iday
00089 integer (kind=int_kind) :: icategory,ind,i_aice,boundid
00090
00091 character (char_len) :: start_time,current_date,current_time
00092 character (len=16) :: c_aice
00093 character (len=8) :: cdate
00094
00095 type(file_desc_t) :: File
00096 type(io_desc_t) :: iodesc2d, iodesc3d
00097 type(var_desc_t) :: varid
00098 real (kind=dbl_kind) :: tmp1d(1)
00099
00100
00101
00102 INTEGER (kind=int_kind), PARAMETER :: ncoord = 4
00103
00104
00105 INTEGER (kind=int_kind), PARAMETER :: nverts = 4
00106
00107
00108
00109 INTEGER (kind=int_kind), PARAMETER :: nvar_verts = 4
00110
00111 TYPE coord_attributes
00112 character (len=11) :: short_name
00113 character (len=45) :: long_name
00114 character (len=20) :: units
00115 END TYPE coord_attributes
00116
00117 TYPE req_attributes
00118 type (coord_attributes) :: req
00119 character (len=20) :: coordinates
00120 END TYPE req_attributes
00121
00122 TYPE(req_attributes), dimension(nvar) :: var
00123 TYPE(coord_attributes), dimension(ncoord) :: coord_var
00124 TYPE(coord_attributes), dimension(nvar_verts) :: var_nverts
00125 CHARACTER (char_len), dimension(ncoord) :: coord_bounds
00126
00127 real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) ::
00128 work1
00129 real (kind=dbl_kind), dimension(nverts,nx_block,ny_block,max_blocks) ::
00130 work3
00131 character(len=char_len_long) ::
00132 filename
00133
00134 real(kind=dbl_kind) :: mrss, mrss0,msize,msize0
00135
00136
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147 ltime = time/int(secday)
00148
00149 if (my_task == master_task) then
00150 call construct_filename(ncfile(ns),'nc',ns)
00151
00152
00153 if (write_ic) then
00154 ncfile(ns) = trim(incond_dir)//ncfile(ns)
00155 else
00156 ncfile(ns) = trim(history_dir)//ncfile(ns)
00157 endif
00158 filename = ncfile(ns)
00159 end if
00160 call broadcast_scalar(filename, master_task)
00161
00162
00163 File%fh=-1
00164 call ice_pio_init(mode='write', filename=trim(filename), File=File, &
00165 clobber=.true., cdf64=lcdf64)
00166
00167 call ice_pio_initdecomp(iodesc=iodesc2d)
00168 call ice_pio_initdecomp(ndim3=nverts, inner_dim=.true., iodesc=iodesc3d)
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181 ltime = time/int(secday)
00182
00183
00184
00185
00186
00187 if (hist_avg) then
00188 status = pio_def_dim(File,'d2',2, boundid)
00189 endif
00190
00191 status = pio_def_dim(File,'ni',nx_global,imtid)
00192 status = pio_def_dim(File,'nj',ny_global,jmtid)
00193 status = pio_def_dim(File,'time',pio_unlimited,timid)
00194 status = pio_def_dim(File,'nvertices',nverts,nvertexid)
00195
00196
00197
00198
00199
00200 status = pio_def_var(File,'time',pio_real,(/timid/),varid)
00201 status = pio_put_att(File,varid,'long_name','model time')
00202
00203 write(cdate,'(i8.8)') idate0
00204 write(title,'(a,a,a,a,a,a,a)') 'days since ', &
00205 cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00'
00206 status = pio_put_att(File,varid,'units',title)
00207
00208 if (days_per_year == 360) then
00209 status = pio_put_att(File,varid,'calendar','360_day')
00210 else
00211 status = pio_put_att(File,varid,'calendar','noleap')
00212 endif
00213
00214 if (hist_avg) then
00215 status = pio_put_att(File,varid,'bounds','time_bounds')
00216 endif
00217
00218
00219
00220
00221
00222 if (hist_avg) then
00223 dimid2(1) = boundid
00224 dimid2(2) = timid
00225 status = pio_def_var(File,'time_bounds',pio_real,dimid2,varid)
00226 status = pio_put_att(File,varid,'long_name', &
00227 'boundaries for time-averaging interval')
00228 write(cdate,'(i8.8)') idate0
00229 write(title,'(a,a,a,a,a,a,a)') 'days since ', &
00230 cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00'
00231 status = pio_put_att(File,varid,'units',title)
00232 endif
00233
00234
00235
00236
00237
00238 ind = 0
00239 ind = ind + 1
00240 coord_var(ind) = coord_attributes('TLON', &
00241 'T grid center longitude', 'degrees_east')
00242 coord_bounds(ind) = 'lont_bounds'
00243 ind = ind + 1
00244 coord_var(ind) = coord_attributes('TLAT', &
00245 'T grid center latitude', 'degrees_north')
00246 coord_bounds(ind) = 'latt_bounds'
00247 ind = ind + 1
00248 coord_var(ind) = coord_attributes('ULON', &
00249 'U grid center longitude', 'degrees_east')
00250 coord_bounds(ind) = 'lonu_bounds'
00251 ind = ind + 1
00252 coord_var(ind) = coord_attributes('ULAT', &
00253 'U grid center latitude', 'degrees_north')
00254 coord_bounds(ind) = 'latu_bounds'
00255
00256
00257
00258
00259
00260 var(n_tarea)%req = coord_attributes('tarea', &
00261 'area of T grid cells', 'm^2')
00262 var(n_tarea)%coordinates = 'TLON TLAT'
00263 var(n_uarea)%req = coord_attributes('uarea', &
00264 'area of U grid cells', 'm^2')
00265 var(n_uarea)%coordinates = 'ULON ULAT'
00266 var(n_dxt)%req = coord_attributes('dxt', &
00267 'T cell width through middle', 'm')
00268 var(n_dxt)%coordinates = 'TLON TLAT'
00269 var(n_dyt)%req = coord_attributes('dyt', &
00270 'T cell height through middle', 'm')
00271 var(n_dyt)%coordinates = 'TLON TLAT'
00272 var(n_dxu)%req = coord_attributes('dxu', &
00273 'U cell width through middle', 'm')
00274 var(n_dxu)%coordinates = 'ULON ULAT'
00275 var(n_dyu)%req = coord_attributes('dyu', &
00276 'U cell height through middle', 'm')
00277 var(n_dyu)%coordinates = 'ULON ULAT'
00278 var(n_HTN)%req = coord_attributes('HTN', &
00279 'T cell width on North side','m')
00280 var(n_HTN)%coordinates = 'TLON TLAT'
00281 var(n_HTE)%req = coord_attributes('HTE', &
00282 'T cell width on East side', 'm')
00283 var(n_HTE)%coordinates = 'TLON TLAT'
00284 var(n_ANGLE)%req = coord_attributes('ANGLE', &
00285 'angle grid makes with latitude line on U grid', &
00286 'radians')
00287 var(n_ANGLE)%coordinates = 'ULON ULAT'
00288 var(n_ANGLET)%req = coord_attributes('ANGLET', &
00289 'angle grid makes with latitude line on T grid', &
00290 'radians')
00291 var(n_ANGLET)%coordinates = 'TLON TLAT'
00292
00293
00294
00295 var_nverts(n_lont_bnds) = coord_attributes('lont_bounds', &
00296 'longitude boundaries of T cells', 'degrees_east')
00297 var_nverts(n_latt_bnds) = coord_attributes('latt_bounds', &
00298 'latitude boundaries of T cells', 'degrees_north')
00299 var_nverts(n_lonu_bnds) = coord_attributes('lonu_bounds', &
00300 'longitude boundaries of U cells', 'degrees_east')
00301 var_nverts(n_latu_bnds) = coord_attributes('latu_bounds', &
00302 'latitude boundaries of U cells', 'degrees_north')
00303
00304
00305
00306
00307
00308 dimid3(1) = imtid
00309 dimid3(2) = jmtid
00310 dimid3(3) = timid
00311
00312 dimid2(1) = imtid
00313 dimid2(2) = jmtid
00314
00315 do i = 1, ncoord
00316 status = pio_def_var(File, coord_var(i)%short_name, pio_real, &
00317 dimid2, varid)
00318
00319 status = pio_put_att(File,varid,'long_name',coord_var(i)%long_name)
00320 status = pio_put_att(File, varid, 'units', coord_var(i)%units)
00321 status = pio_put_att(File, varid, 'missing_value', spval)
00322 status = pio_put_att(File,varid,'_FillValue',spval)
00323 if (coord_var(i)%short_name == 'ULAT') then
00324 status = pio_put_att(File,varid,'comment', &
00325 'Latitude of NE corner of T grid cell')
00326 endif
00327
00328 if (f_bounds) then
00329 status = pio_put_att(File, varid, 'bounds', coord_bounds(i))
00330 endif
00331 enddo
00332
00333
00334 if (igrd(n_tmask)) then
00335 status = pio_def_var(File,'tmask', pio_real, dimid2, varid)
00336 status = pio_put_att(File, varid, 'long_name', 'ocean grid mask')
00337 status = pio_put_att(File, varid, 'coordinates', 'TLON TLAT')
00338 status = pio_put_att(File, varid, 'comment', '0 = land, 1 = ocean')
00339 endif
00340
00341 do i = 2, nvar
00342 if (igrd(i)) then
00343 status = pio_def_var(File, var(i)%req%short_name, &
00344 pio_real, dimid2, varid)
00345 status = pio_put_att(File, varid, 'long_name', var(i)%req%long_name)
00346 status = pio_put_att(File, varid, 'units', var(i)%req%units)
00347 status = pio_put_att(File, varid, 'coordinates', var(i)%coordinates)
00348 endif
00349 enddo
00350
00351
00352 dimid_nverts(1) = nvertexid
00353 dimid_nverts(2) = imtid
00354 dimid_nverts(3) = jmtid
00355 do i = 1, nvar_verts
00356 if (f_bounds) then
00357 status = pio_def_var(File, var_nverts(i)%short_name, &
00358 pio_real,dimid_nverts, varid)
00359 status = &
00360 pio_put_att(File,varid, 'long_name', var_nverts(i)%long_name)
00361 status = &
00362 pio_put_att(File, varid, 'units', var_nverts(i)%units)
00363 endif
00364 enddo
00365
00366 do n=1,num_avail_hist_fields
00367
00368 if (avail_hist_fields(n)%vhistfreq == histfreq(ns).or.write_ic) then
00369
00370 status = pio_def_var(File, avail_hist_fields(n)%vname, &
00371 pio_real, dimid3, varid)
00372 status = pio_put_att(File,varid,'units', &
00373 avail_hist_fields(n)%vunit)
00374 status = pio_put_att(File,varid, 'long_name', &
00375 avail_hist_fields(n)%vdesc)
00376
00377 status = pio_put_att(File,varid,'coordinates', &
00378 avail_hist_fields(n)%vcoord)
00379 status = pio_put_att(File,varid,'cell_measures', &
00380 avail_hist_fields(n)%vcellmeas)
00381 status = pio_put_att(File,varid,'missing_value',spval)
00382 status = pio_put_att(File,varid,'_FillValue',spval)
00383
00384
00385
00386
00387
00388 c_aice = TRIM(avail_hist_fields(n)%vname)
00389 i_aice = lenstr(c_aice)
00390 if (i_aice > 4 .and. c_aice(1:5) == 'aicen') then
00391 read(c_aice(6:9), '(i3)') icategory
00392 avail_hist_fields(n)%vcomment = &
00393 'Ice range: '//c_hi_range(icategory)
00394 endif
00395 status = pio_put_att(File,varid,'comment', &
00396 avail_hist_fields(n)%vcomment)
00397
00398
00399
00400
00401 if (hist_avg) then
00402 if (TRIM(avail_hist_fields(n)%vname)/='sig1' &
00403 .or.TRIM(avail_hist_fields(n)%vname)/='sig2') then
00404 status = pio_put_att(File,varid,'cell_methods','time: mean')
00405 endif
00406 endif
00407
00408
00409 if (histfreq(ns) == '1' .or. .not. hist_avg &
00410
00411 .or. n==n_sig1(ns) .or. n==n_sig2(ns) &
00412 .or. n==n_trsig(ns) &
00413 .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) &
00414 .or. n==n_hisnap(ns) .or. n==n_aisnap(ns) &
00415 .or. n==n_FY(ns)) then
00416 status = pio_put_att(File,varid,'time_rep','instantaneous')
00417 else
00418 status = pio_put_att(File,varid,'time_rep','averaged')
00419 endif
00420
00421 endif
00422 enddo
00423
00424
00425
00426
00427
00428
00429 #ifdef CCSMCOUPLED
00430 status = pio_put_att(File,pio_global,'title',runid)
00431 #else
00432 title = 'sea ice model output for CICE'
00433 status = pio_put_att(File,pio_global,'title',title)
00434 #endif
00435 title = 'Diagnostic and Prognostic Variables'
00436 status = pio_put_att(File,pio_global,'contents',title)
00437
00438 title = 'sea ice model: Community Ice Code (CICE)'
00439 status = pio_put_att(File,pio_global,'source',title)
00440
00441 write(title,'(a,i3,a)') 'All years have exactly ',int(dayyr),' days'
00442 status = pio_put_att(File,pio_global,'comment',title)
00443
00444 write(title,'(a,i8)') 'File written on model date ',idate
00445 status = pio_put_att(File,pio_global,'comment2',title)
00446
00447 write(title,'(a,i6)') 'seconds elapsed into model date: ',sec
00448 status = pio_put_att(File,pio_global,'comment3',title)
00449
00450 title = 'CF-1.0'
00451 status = &
00452 pio_put_att(File,pio_global,'conventions',title)
00453
00454 call date_and_time(date=current_date, time=current_time)
00455 write(start_time,1000) current_date(1:4), current_date(5:6), &
00456 current_date(7:8), current_time(1:2), &
00457 current_time(3:4), current_time(5:8)
00458 1000 format('This dataset was created on ', &
00459 a,'-',a,'-',a,' at ',a,':',a,':',a)
00460
00461 status = pio_put_att(File,pio_global,'history',start_time)
00462
00463
00464
00465
00466
00467 status = pio_enddef(File)
00468
00469
00470
00471
00472
00473 status = pio_inq_varid(File,'time',varid)
00474 status = pio_put_var(File,varid,ltime)
00475
00476
00477
00478
00479
00480 if (hist_avg) then
00481 status = pio_inq_varid(File,'time_bounds',varid)
00482 time_bounds=(/time_beg(ns),time_end(ns)/)
00483 status = pio_put_var(File,varid,time_bounds)
00484 endif
00485
00486
00487
00488
00489
00490 do i = 1,ncoord
00491 call broadcast_scalar(coord_var(i)%short_name,master_task)
00492
00493 SELECT CASE (coord_var(i)%short_name)
00494 CASE ('TLON')
00495 work1(:,:,:) = tlon(:,:,:)
00496
00497 work1 = work1*rad_to_deg + c360
00498 where (work1 > c360) work1 = work1 - c360
00499 where (work1 < c0 ) work1 = work1 + c360
00500 CASE ('TLAT')
00501 work1(:,:,:) = tlat(:,:,:)
00502 work1 = work1*rad_to_deg
00503 CASE ('ULON')
00504 work1(:,:,:) = ulon(:,:,:)
00505 work1 = work1*rad_to_deg
00506 CASE ('ULAT')
00507 work1(:,:,:) = ulat(:,:,:)
00508 work1 = work1*rad_to_deg
00509 END SELECT
00510
00511 status = pio_inq_varid(File, coord_var(i)%short_name, varid)
00512 call pio_write_darray(File, varid, iodesc2d, &
00513 transfer(work1,tmp1d), status, fillval=spval_dbl)
00514 enddo
00515
00516
00517
00518
00519
00520 if (igrd(n_tmask)) then
00521 work1(:,:,:) = hm(:,:,:)
00522 status = pio_inq_varid(File, 'tmask', varid)
00523 call pio_write_darray(File, varid, iodesc2d, &
00524 transfer(work1,tmp1d), status, fillval=spval_dbl)
00525 endif
00526
00527 do i = 2,nvar
00528 if (igrd(i)) then
00529 call broadcast_scalar(var(i)%req%short_name,master_task)
00530
00531 SELECT CASE (var(i)%req%short_name)
00532 CASE ('tarea')
00533 work1(:,:,:) = tarea(:,:,:)
00534 CASE ('uarea')
00535 work1(:,:,:) = uarea(:,:,:)
00536 CASE ('dxu')
00537 work1(:,:,:) = dxu(:,:,:)
00538 CASE ('dyu')
00539 work1(:,:,:) = dyu(:,:,:)
00540 CASE ('dxt')
00541 work1(:,:,:) = dxt(:,:,:)
00542 CASE ('dyt')
00543 work1(:,:,:) = dyt(:,:,:)
00544 CASE ('HTN')
00545 work1(:,:,:) = HTN(:,:,:)
00546 CASE ('HTE')
00547 work1(:,:,:) = HTE(:,:,:)
00548 CASE ('ANGLE')
00549 work1(:,:,:) = ANGLE(:,:,:)
00550 CASE ('ANGLET')
00551 work1(:,:,:) = ANGLET(:,:,:)
00552 END SELECT
00553
00554 status = pio_inq_varid(File, var(i)%req%short_name, varid)
00555 call pio_write_darray(File, varid, iodesc2d, &
00556 transfer(work1,tmp1d), status, fillval=spval_dbl)
00557 endif
00558 enddo
00559
00560
00561
00562
00563
00564 if (f_bounds) then
00565 work3(:,:,:,:) = c0
00566 do i = 1, nvar_verts
00567 call broadcast_scalar(var_nverts(i)%short_name,master_task)
00568 SELECT CASE (var_nverts(i)%short_name)
00569 CASE ('lont_bounds')
00570 do ivertex = 1, nverts
00571 work3(ivertex,:,:,:) = lont_bounds(ivertex,:,:,:)
00572 enddo
00573 CASE ('latt_bounds')
00574 do ivertex = 1, nverts
00575 work3(ivertex,:,:,:) = latt_bounds(ivertex,:,:,:)
00576 enddo
00577 CASE ('lonu_bounds')
00578 do ivertex = 1, nverts
00579 work3(ivertex,:,:,:) = lonu_bounds(ivertex,:,:,:)
00580 enddo
00581 CASE ('latu_bounds')
00582 do ivertex = 1, nverts
00583 work3(ivertex,:,:,:) = latu_bounds(ivertex,:,:,:)
00584 enddo
00585 END SELECT
00586
00587 status = pio_inq_varid(File, var_nverts(i)%short_name, varid)
00588 call pio_write_darray(File, varid, iodesc3d, &
00589 transfer(work3,tmp1d), status, fillval=spval_dbl)
00590 enddo
00591 endif
00592
00593
00594
00595
00596
00597 do n=1,num_avail_hist_fields
00598 if (avail_hist_fields(n)%vhistfreq == histfreq(ns).or.write_ic) then
00599 status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid)
00600 if (status /= PIO_noerr) call abort_ice( &
00601 'ice: Error getting varid for '//avail_hist_fields(n)%vname)
00602 work1(:,:,:) = aa(:,:,n,:)
00603 call pio_setframe(varid, int(1,kind=PIO_OFFSET))
00604 call pio_write_darray(File, varid, iodesc2d,&
00605 transfer(work1,tmp1d), status, fillval=spval_dbl)
00606 endif
00607 enddo
00608
00609
00610
00611
00612
00613 call pio_closefile(File)
00614 if (my_task == master_task) then
00615 write(nu_diag,*) ' '
00616 write(nu_diag,*) 'Finished writing ',trim(ncfile(ns))
00617 endif
00618
00619
00620
00621
00622
00623
00624
00625
00626
00627
00628
00629
00630
00631
00632 call pio_freedecomp(File,iodesc2d)
00633 call pio_freedecomp(File,iodesc3d)
00634
00635 call ice_pio_finalize
00636
00637
00638
00639
00640
00641
00642
00643
00644
00645
00646
00647
00648 105 format( A, f10.2, A, f10.2, A)
00649
00650 end subroutine icecdf
00651
00652
00653
00654
00655
00656
00657
00658
00659
00660
00661
00662
00663 subroutine icebin(ns)
00664
00665
00666
00667
00668
00669
00670
00671
00672
00673
00674
00675 use ice_gather_scatter
00676 use ice_domain_size
00677 use ice_constants
00678 use ice_restart, only: lenstr, runid
00679 use ice_itd, only: c_hi_range
00680 use ice_calendar, only: write_ic, dayyr, histfreq
00681
00682
00683
00684
00685
00686 integer (kind=int_kind), intent(in) :: ns
00687
00688 integer (kind=int_kind) :: i,j,n,nrec,nbits
00689 character (char_len) :: title
00690 character (char_len_long) :: ncfile(max_nstrm), hdrfile
00691
00692 integer (kind=int_kind) :: icategory,i_aice
00693
00694 character (char_len) :: current_date,current_time
00695 character (len=16) :: c_aice
00696 logical (kind=log_kind) :: dbug
00697
00698 dbug = .false.
00699
00700 if (my_task == master_task) then
00701
00702 call construct_filename(ncfile(ns),'da',ns)
00703
00704
00705 if (write_ic) then
00706 ncfile(ns) = trim(incond_dir)//ncfile(ns)
00707 else
00708 ncfile(ns) = trim(history_dir)//ncfile(ns)
00709 endif
00710 hdrfile = trim(ncfile(ns))//'.hdr'
00711
00712
00713
00714
00715 nbits = 32
00716 call ice_open(nu_history, ncfile(ns), nbits)
00717 open(nu_hdr,file=hdrfile,form='formatted',status='unknown')
00718
00719
00720
00721 title = 'sea ice model: Community Ice Code (CICE)'
00722 write (nu_hdr, 999) 'source',title,' '
00723
00724 write (nu_hdr, 999) 'file name contains model date',trim(ncfile(ns)),' '
00725 #ifdef CCSMCOUPLED
00726 write (nu_hdr, 999) 'runid',runid,' '
00727 #endif
00728 write (nu_hdr, 999) 'calendar','noleap',' '
00729 write (title,'(a,i3,a)') 'All years have exactly ',int(dayyr),' days'
00730 write (nu_hdr, 999) 'comment',title,' '
00731 write (nu_hdr, 999) 'conventions','CICE',' '
00732 write (nu_hdr, 997) 'missing_value',spval
00733 write (nu_hdr, 997) '_FillValue',spval
00734
00735 call date_and_time(date=current_date, time=current_time)
00736 write (nu_hdr,1000) current_date(1:4), current_date(5:6), &
00737 current_date(7:8), current_time(1:2), &
00738 current_time(3:4), current_time(5:8)
00739 write (nu_hdr, * ) ' '
00740 write (nu_hdr, * ) 'Grid size:'
00741 write (nu_hdr, 998) ' ni',nx_global
00742 write (nu_hdr, 998) ' nj',ny_global
00743
00744 write (nu_hdr, * ) 'Grid variables: (left column = nrec)'
00745 nrec = 1
00746 write (nu_hdr, 996) nrec,'tarea','area of T grid cells','m^2'
00747 write (nu_hdr, * ) 'History variables: (left column = nrec)'
00748 endif
00749
00750 call ice_write(nu_history, nrec, tarea, 'rda4', dbug)
00751
00752 do n=1,num_avail_hist_fields
00753
00754 if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) then
00755
00756 nrec = nrec + 1
00757 if (my_task == master_task) then
00758 write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
00759 trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit)
00760
00761
00762 c_aice = TRIM(avail_hist_fields(n)%vname)
00763 i_aice = lenstr(c_aice)
00764 if (i_aice > 4 .and. c_aice(1:5) == 'aicen') then
00765 read(c_aice(6:9), '(i3)') icategory
00766 avail_hist_fields(n)%vcomment = &
00767 'Ice range: '//c_hi_range(icategory)
00768 endif
00769 write (nu_hdr, 995) nrec,trim(avail_hist_fields(n)%vname), &
00770 trim(avail_hist_fields(n)%vcomment)
00771
00772
00773 if (histfreq(ns) == '1' .or. .not. hist_avg &
00774
00775 .or. n==n_sig1(ns) .or. n==n_sig2(ns) &
00776 .or. n==n_trsig(ns) &
00777 .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) &
00778 .or. n==n_hisnap(ns) .or. n==n_aisnap(ns) &
00779 .or. n==n_FY(ns)) then
00780 write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
00781 'time_rep','instantaneous'
00782 else
00783 write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
00784 'time_rep','averaged'
00785 endif
00786 endif
00787
00788 call ice_write(nu_history, nrec, aa(:,:,n,:), 'rda4', dbug)
00789
00790 endif
00791 enddo
00792
00793 995 format(i3,2x,a,' comment: ',a)
00794 996 format(i3,2x,a,': ',a,',',2x,a)
00795 997 format(a,': ',es13.6)
00796 998 format(a,': ',i6)
00797 999 format(a,': ',a,2x,a)
00798 1000 format('This dataset was created on ', &
00799 a,'-',a,'-',a,' at ',a,':',a,':',a)
00800
00801 if (my_task == master_task) then
00802 close (nu_hdr)
00803 close (nu_history)
00804 write (nu_diag,*) ' '
00805 write (nu_diag,*) 'Finished writing ',trim(ncfile(ns))
00806 endif
00807
00808 end subroutine icebin
00809
00810 end module ice_history_write