! subroutine mkncfile(fileout,eof) use readmgcm,only: nlon,nlat,glon,glat,xls,time,frd,sourcefile implicit none #include ! ! Args: character(len=*),intent(in) :: fileout integer,intent(in) :: eof ! ! Local: integer,save :: ncalls=0 integer,save :: id_time,id_lon,id_lat,idv_tn,idv_un, | idv_vn,idv_z,idv_time,idv_xls,idv_lon,idv_lat integer :: istat,ncid,itime,start_3d(3),count_3d(3), | id1(1),id_dims3(3) integer,parameter :: secs_per_sol = 86400 character(len=120) :: char120 character(len=80) :: char80 ! ncalls = ncalls+1 itime = ncalls if (ncalls > 1) goto 100 ! ! If first call, create new netcdf dataset and define metadata. ! istat = nf_create(fileout,NF_CLOBBER,ncid) if (istat /= NF_NOERR) then write(char120,"('Error return from nf_create for netcdf ', | 'file ',a)") trim(fileout) call handle_ncerr(istat,char120) ncid = 0 else write(6,"('Created netcdf file ',a,' ncid=',i8)") | trim(fileout),ncid endif ! ! Define unlimited dimension: istat = nf_def_dim(ncid,"time",NF_UNLIMITED,id_time) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining time dimension') ! ! Define grid dimensions: istat = nf_def_dim(ncid,"lon",nlon,id_lon) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining longitude dimension') istat = nf_def_dim(ncid,"lat",nlat,id_lat) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining latitude dimension') ! ! Define coord arrays lat,lon: id1(1) = id_lon istat = nf_def_var(ncid,"LON",NF_DOUBLE,1,id1,idv_lon) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining longitude dimension variable') write(char80,"('Geographic Longitude (-west, +east)')") istat = nf_put_att_text(ncid,idv_lon,"long_name", | len_trim(char80),trim(char80)) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining long_name of longitude dimension variable') istat = nf_put_att_text(ncid,idv_lon,"units",12,'degrees_east') if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining units of longitude dimension variable') ! id1(1) = id_lat istat = nf_def_var(ncid,"LAT",NF_DOUBLE,1,id1,idv_lat) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining latitude dimension variable') write(char80,"('Geographic Latitude (-south, +north)')") istat = nf_put_att_text(ncid,idv_lat,"long_name", | len_trim(char80),trim(char80)) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining long_name of latitude dimension variable') istat = nf_put_att_text(ncid,idv_lat,"units",13,'degrees_north') if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining units of latitude dimension variable') ! ! Define time(time): istat = nf_def_var(ncid,'TIME',NF_DOUBLE,1,id_time,idv_time) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining TIME variable') write(char80,"('MGCM model time')") istat = nf_put_att_text(ncid,idv_time,"long_name", | len_trim(char80),trim(char80)) istat = nf_put_att_text(ncid,idv_time,"units",4,"SOLS") ! ! Define xls(time): istat = nf_def_var(ncid,'XLS',NF_DOUBLE,1,id_time,idv_xls) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error defining XLS variable') write(char80,"('MGCM aerocentric longitude Ls')") istat = nf_put_att_text(ncid,idv_xls,"long_name", | len_trim(char80),trim(char80)) istat = nf_put_att_text(ncid,idv_xls,"units",7,"DEGREES") write(char80,"('Ls = 90 is Northern Summer Solstice')") istat = nf_put_att_text(ncid,idv_xls,"Ls90", | len_trim(char80),trim(char80)) write(char80,"('Ls = 270 is Southern Summer Solstice')") istat = nf_put_att_text(ncid,idv_xls,"Ls270", | len_trim(char80),trim(char80)) write(char80,"('Ls = 0 and 180 refer to the Equinoxes')") istat = nf_put_att_text(ncid,idv_xls,"Ls180", | len_trim(char80),trim(char80)) ! ! Define the 4 fields z,t,u,v: id_dims3(1) = id_lon id_dims3(2) = id_lat id_dims3(3) = id_time ! ! Define Z: istat = nf_def_var(ncid,'Z',NF_DOUBLE,3,id_dims3,idv_z) if (istat /= NF_NOERR) | call handle_ncerr(istat,'Error defining 3-d field Z.') write(char80,"('GEOPOTENTIAL HEIGHT')") istat = nf_put_att_text(ncid,idv_z,"long_name", | len_trim(char80),trim(char80)) istat = nf_put_att_text(ncid,idv_z,"units",2,"KM") ! ! Define TN: istat = nf_def_var(ncid,'TN',NF_DOUBLE,3,id_dims3,idv_tn) if (istat /= NF_NOERR) | call handle_ncerr(istat,'Error defining 3-d field TN.') write(char80,"('NEUTRAL TEMPERATURE')") istat = nf_put_att_text(ncid,idv_tn,"long_name", | len_trim(char80),trim(char80)) istat = nf_put_att_text(ncid,idv_tn,"units",5,"DEG K") ! ! Define UN: istat = nf_def_var(ncid,'UN',NF_DOUBLE,3,id_dims3,idv_un) if (istat /= NF_NOERR) | call handle_ncerr(istat,'Error defining 3-d field UN.') write(char80,"('NEUTRAL ZONAL WIND (-WEST, +EAST)')") istat = nf_put_att_text(ncid,idv_un,"long_name", | len_trim(char80),trim(char80)) istat = nf_put_att_text(ncid,idv_un,"units",4,"CM/S") ! ! Define VN: istat = nf_def_var(ncid,'VN',NF_DOUBLE,3,id_dims3,idv_vn) if (istat /= NF_NOERR) | call handle_ncerr(istat,'Error defining 3-d field VN.') write(char80,"('NEUTRAL MERIDIONAL WIND (-SOUTH, +NORTH)')") istat = nf_put_att_text(ncid,idv_vn,"long_name", | len_trim(char80),trim(char80)) istat = nf_put_att_text(ncid,idv_vn,"units",4,"CM/S") ! ! Define global attributes: write(char120,"('Data from NASA Ames Mars MGCM C-grid model at', | ' 1.32-microbar level, for import to Mars tgcm (MTGCM)')") istat = nf_put_att_text(ncid,NF_GLOBAL,"Description", | len_trim(char120),trim(char120)) ! ! Seconds per sol: istat = nf_put_att_int(ncid,NF_GLOBAL,"Seconds_per_Sol",NF_INT, | 1,secs_per_sol) ! ! Source file: istat = nf_put_att_text(ncid,NF_GLOBAL,"IEEE_Source_File", | len_trim(sourcefile),trim(sourcefile)) ! ! Person contacts: char120 = ' ' write(char120,"('Dr. Steve Bougher: bougher@engin.umich.edu')") istat = nf_put_att_text(ncid,NF_GLOBAL,"MTGCM_Contact", | len_trim(char120),trim(char120)) ! write(char120,"('Dr. Jim Murphy, murphy@nmsu.edu')") istat = nf_put_att_text(ncid,NF_GLOBAL,"MGCM_Contact", | len_trim(char120),trim(char120)) ! write(char120,"('Ben Foster, foster@hao.ucar.edu')") istat = nf_put_att_text(ncid,NF_GLOBAL,"Software_Engineer", | len_trim(char120),trim(char120)) ! ! Exit define mode: istat = nf_enddef(ncid) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error return from nf_enddef') ! ! Give values to lon,lat grid coord vars: istat = nf_put_var_real(ncid,idv_lon,glon) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error assigning values to longitude array.') ! istat = nf_put_var_real(ncid,idv_lat,glat) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error assigning values to latitude array.') ! 100 continue ! end ncalls==1, begin ncalls > 1 (file has been defined) ! ! Give values to variables dimensioned in time: ! ! TIME: istat = nf_put_var1_real(ncid,idv_time,itime,time) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error giving value to TIME var.') ! ! XLS: istat = nf_put_var1_real(ncid,idv_xls,itime,xls) if (istat /= NF_NOERR) call handle_ncerr(istat, | 'Error giving value to XLS var.') ! ! Set up for 3d fields ! start_3d(1:2) = 1 start_3d(3) = itime count_3d(1) = nlon count_3d(2) = nlat count_3d(3) = 1 ! ! (IMPORTANT NOTE: order of fields in frd is z,t,u,v). ! ! Z: istat = nf_put_vara_real(ncid,idv_z,start_3d,count_3d,frd(:,:,1)) if (istat /= NF_NOERR) then call handle_ncerr(istat,'Error giving values to 3d var Z.') stop 'Write Z' endif ! ! TN: istat = nf_put_vara_real(ncid,idv_tn,start_3d,count_3d,frd(:,:,2)) if (istat /= NF_NOERR) then call handle_ncerr(istat,'Error giving values to 3d var TN.') stop 'Write TN' endif ! ! UN: istat = nf_put_vara_real(ncid,idv_un,start_3d,count_3d,frd(:,:,3)) if (istat /= NF_NOERR) then call handle_ncerr(istat,'Error giving values to 3d var UN.') stop 'Write UN' endif ! ! VN: istat = nf_put_vara_real(ncid,idv_vn,start_3d,count_3d,frd(:,:,4)) if (istat /= NF_NOERR) then call handle_ncerr(istat,'Error giving values to 3d var VN.') stop 'Write VN' endif ! ! Close file (eof was encountered on input binary file): if (eof > 0) then istat = nf_close(ncid,fileout) if (istat /= NF_NOERR) then write(char120,"('Error closing netcdf file ',a)")trim(fileout) call handle_ncerr(istat,char120) else write(6,"('mkncfile: closed netcdf file ',a)") trim(fileout) endif endif end subroutine mkncfile !----------------------------------------------------------------------- subroutine handle_ncerr(istat,msg) implicit none #include "netcdf.inc" ! ! Handle a netcdf lib error: ! integer,intent(in) :: istat character(len=*),intent(in) :: msg ! write(6,"(/72('-'))") write(6,"('>>> Error from netcdf library:')") write(6,"(a)") trim(msg) write(6,"('istat=',i5)") istat write(6,"(a)") nf_strerror(istat) write(6,"(72('-')/)") return end subroutine handle_ncerr