program mkbynam C 11/08: /aurora/d/emery/tie1-83/mkbynam.f Code to read IMF and make hourly C namelist of BYIMF_TIME for Heelis runs implicit none #include ! ! All are allocated ndata except for f107d and f107a, which ! are allocated ndays. ! integer :: ndata,ndays real :: missing integer,dimension(:),target,allocatable :: days real,dimension(:),target,allocatable :: date, | uttime,yfrac,bx,by,bz,swvel,swden,kp,f107d,f107a real :: bynam(24) integer :: namd(24),namh(24),namm(24) integer :: nh ! ! Local: integer :: index,index1,index2,i integer(kind=8) :: iyearday,iyyyy,iddd real :: model_uttime,model_yfrac,model_dfrac,ut integer :: istat,ncid,istep character(len=240) :: char240,dskfile integer :: id_ndata,id_ndays,idv_date,idv_uttime,idv_yfrac, | idv_bx,idv_by,idv_bz,idv_swvel,idv_swden,idv_kp, | idv_f107d,idv_f107a,idv_missing,idv_days real :: fmin,fmax ! ! Acquire mss file: dskfile = ' ' call getfile('$TGCMDATA/imf_2000001-2006365.nc',dskfile) ! write(6,"(/,72('-'))") write(6,"('RDIMF: read IMF data file:')") write(6,"('Opened netcdf imf data file ',a)") trim(dskfile) ! ! Open imf data file: istat = nf_open(dskfile,NF_NOWRITE,ncid) if (istat /= NF_NOERR) then write(char240,"('Error opening imf_ncfile ',a)")trim(imf_ncfile) call handle_ncerr(istat,char240) call shutdown('imf_ncfile') endif write(6,"(/,72('-'))") write(6,"('RDIMF: read IMF data file:')") write(6,"('Opened netcdf imf data file ',a)") trim(dskfile) ! ! Get ndata and ndays dimensions: istat = nf_inq_dimid(ncid,'ndata',id_ndata) istat = nf_inq_dimlen(ncid,id_ndata,ndata) istat = nf_inq_dimid(ncid,'ndays',id_ndays) istat = nf_inq_dimlen(ncid,id_ndays,ndays) write(6,"('rdimf: ndata=',i5,' ndays=',i5)") ndata,ndays ! ! Allocate data arrays and initialize structure pointers: if (istep==1) call alloc_imf ! ! Read vars: ! ! Missing value: istat = nf_inq_varid(ncid,'missing',idv_missing) istat = nf_get_var_double(ncid,idv_missing,missing) ! yyyyddd.frac date: istat = nf_inq_varid(ncid,'date',idv_date) istat = nf_get_var_double(ncid,idv_date,date) write(6,"('rdimf: date min,max=',2f14.4)") | minval(date),maxval(date) ! integer days: istat = nf_inq_varid(ncid,'days',idv_days) istat = nf_get_var_int(ncid,idv_days,days) write(6,"('rdimf: ndays=',i5,' days min,max=',2i10)") | ndays,minval(days),maxval(days) ! uttime: ! istat = nf_inq_varid(ncid,'uttime',idv_uttime) ! istat = nf_get_var_double(ncid,idv_uttime,uttime) ! write(6,"('rdimf: uttime min,max=',2f6.2)") ! | minval(uttime),maxval(uttime) ! yfrac: ! istat = nf_inq_varid(ncid,'yfrac',idv_yfrac) ! istat = nf_get_var_double(ncid,idv_yfrac,yfrac) ! write(6,"('rdimf: yfrac min,max=',2f10.4)") ! | minval(yfrac),maxval(yfrac) ! bx: istat = nf_inq_varid(ncid,'bx',idv_bx) istat = nf_get_var_double(ncid,idv_bx,bx) call fminmaxspv(bx,ndata,fmin,fmax,missing) write(6,"('rdimf: bx min,max=',2f8.2)") fmin,fmax ! by: istat = nf_inq_varid(ncid,'by',idv_by) istat = nf_get_var_double(ncid,idv_by,by) call fminmaxspv(by,ndata,fmin,fmax,missing) write(6,"('rdimf: by min,max=',2f8.2)") fmin,fmax ! bz: istat = nf_inq_varid(ncid,'bz',idv_bz) istat = nf_get_var_double(ncid,idv_bz,bz) call fminmaxspv(bz,ndata,fmin,fmax,missing) write(6,"('rdimf: bz min,max=',2f8.2)") fmin,fmax ! swvel: istat = nf_inq_varid(ncid,'swvel',idv_swvel) istat = nf_get_var_double(ncid,idv_swvel,swvel) call fminmaxspv(swvel,ndata,fmin,fmax,missing) write(6,"('rdimf: swvel min,max=',2f8.2)") fmin,fmax ! swden: istat = nf_inq_varid(ncid,'swden',idv_swden) istat = nf_get_var_double(ncid,idv_swden,swden) call fminmaxspv(swden,ndata,fmin,fmax,missing) write(6,"('rdimf: swden min,max=',2f8.2)") fmin,fmax ! kp: istat = nf_inq_varid(ncid,'kp',idv_kp) istat = nf_get_var_double(ncid,idv_kp,kp) call fminmaxspv(kp,ndata,fmin,fmax,missing) write(6,"('rdimf: kp min,max=',2f8.2)") fmin,fmax ! f107d: istat = nf_inq_varid(ncid,'f107d',idv_f107d) istat = nf_get_var_double(ncid,idv_f107d,f107d) call fminmaxspv(f107d,ndays,fmin,fmax,missing) write(6,"('rdimf: f107d min,max=',2f8.2)") fmin,fmax ! f107a: istat = nf_inq_varid(ncid,'f107a',idv_f107a) istat = nf_get_var_double(ncid,idv_f107a,f107a) call fminmaxspv(f107a,ndays,fmin,fmax,missing) write(6,"('rdimf: f107a min,max=',2f8.2)") fmin,fmax write(6,"(72('-'),/)") ! Calculate model ut, fractional year-day, and fractional year: iyear = 2002 do iday=104,112 nh = 0 do ihr=0,23 nh = nh + 1 namd(nh) = iday namh(nh) = ihr namm(nh) = 30 ! Get on hour to be sure is index1 is previous hr and index2 is present hour model_uttime = ihr+0.0 model_dfrac = iyear*1000+real(iday)+model_uttime/24. model_yfrac=iyear+(iday-1+(model_uttime/24.))/365. ! ! Calculate fractional year yfrac from data date (yyyyddd.dayfrac): do i=1,ndata iyearday = int(date(i),8) ! yyyyddd (long int) iyyyy = iyearday/1000 ! yyyy iddd = iyearday-iyyyy*1000 ! ddd ut = (date(i)-float(iyearday))*24. ! ut (hrs) yfrac(i) = float(iyyyy)+(float(iddd-1)+ut/24.)/365. ! write(6,"('getimf: i=',i6,' ndata=',i6,' date(i)=',f14.4, ! | ' year=',i4,' day=',i3,' ut=',f6.2,' yfrac(i)=',f9.4)") ! | i,ndata,date(i),iyyyy,iddd,ut,yfrac(i) enddo ! ! Check that requested date is available: if (model_yfrac < yfrac(1) .or. | model_yfrac > yfrac(ndata)) then write(6,"(/,'>>> getimf: requested date ',f14.4,' is not ', | 'available from imf_ncfile data file:')") model_yfrac write(6,"(4x,' date(1)=',f14.4, | ' date(ndata)=',f14.4,/)") date(1),date(ndata) call shutdown('getimf') endif ! ! Get index to requested time: index=real_bsearch(yfrac,1,ndata,model_yfrac) if (index == -1) then write(6,"('>>> getimf: error from real_bsearch: ', | 'could not find model_yfrac=',i8)") model_yfrac call shutdown('getimf: IMF data not available') endif index1=index index2=index+1 bynam(nh) = by(index2) enddo ! Print out namelist for hourly By write (6,"(13x,4(i3.3,',',i2.2,',',i2.2,',',f6.2,', '))") | (namd(n),namh(n),namm(n),bynam(n),n=1,24) enddo end