subroutine wreadnc(ncid,ntim,ntimr,tim,dt1,time,dt) c c Read results from the *.nc file, and write them into the c the itg common variables c c first read the dimensions and integers, then c call separate subroutine to read rest of file include 'netcdf.inc' include 'itg.par' include 'itg.cmn' include 'itg_ri.cmn' c subroutine arguments integer ncid ! netCDF id number, from earlier nc_open call integer ntim ! on entry, the index of last time point previously used. ! on exit, the index of the last time point now used. integer ntimr ! the same as for ntim, but for a copy of phi(t) at ! one point with potentially higher time resolution ! for frequency calculation. c output subroutine arguments: real tim ! final time when dt is the same for all ky's. real dt1 ! time step when dt is the same for all ky's. real time(mz) ! final time for various ky's. real dt(mz) ! time step for various ky's. c local variables: character filename*80 ! name of netcdf file character yorn*1,varname*30 integer iwrite,ri,i,nt,k,id,l,m,n,n1,nt1,nt2,ntf1,ntf2 integer i2(2),i3(3),i4(4),i5(5) real xxx integer emflag c netCDF variables integer status integer nvars integer ld_dim,md_dim,nd_dim,nspec_dim,kdpass_dim integer time_dim,nd4_dim,md4_dim,id_dim,timef_dim integer ri_dim ! extra dimension for complex vars, 1=real 2=im integer ld_id,ldb_id,kd_id,nmin_id,nmax_id ! variable tags integer nd_id,lin_id,kdpass_id integer nstp_id,ne_id,nfreq_id,ntotal_id,md_id integer nspecies_id,ikx_id,iperiod_id integer nparmom_id,nperpmom_id,nemom_id integer iodd_id,iflr_id,iphi00_id,ifilter_id integer iexp_id,igradon_id,mmin_id,mmax_id integer tim_id,shr_id,qsf_id,epsn_id,eps_id,epse_id integer nueeff_id,alpha_id,etai_id,nuii_id,rmu1_id integer tiovte_id,vy_id,vz_id,rmime_id,etae_id integer dt0_id,dt1_id,x0_id,y0_id,z0_id,xp_id integer density_id,potential_id,u_par_id,t_par_id,t_perp_id integer q_perp_id,q_par_id integer e_density_id,e_p_id,e_r_id,e_t_id,phi_ba_id integer phi_bk_id,e_denk_id,phiav_old_id integer semi_imp_id,beta_e_id ! EM integer n_e_id,u_e_id,tpar_e_id,tperp_e_id,apar_id ! EM integer timo_id,wpfx_id,fluxi_id,fluxe_id,qfluxe_id integer gamx_id,mgamx_id,wenx_id,pcerr_id integer wtif_id,wkif_id,wpif_id,psp_id,wenrk_id integer dketa_id,dkpar_id,dktor_id,dktdp_id,wkups_id integer dktot_id,cx_id,cy_id,cz_id,czc_id,czcn_id integer czn_id,cznn_id,ct_id,ninterv_id integer phi00_id,den00_id,upar00_id,tpar00_id,tperp00_id integer drt_id,phirms_id,wakx_id,waky_id,wxsp_id integer utor_id,upol_id,grtmx_id,utim_id,timf_id integer avgflag_id,date_id,note_id integer rmass_id,charge_id,n_I_id,Ln_id,eta_id integer eta_par_id,tau_id,phisq_id,fluximn_id integer qfluximn_id,fluxemn_id,qfluxemn_id integer icrit_id,ntimf_id,time_id,dt_id,rkperp2_id c setup real-imaginary arrays for some complex variables c (no complex type in netCDF) real phi00_ri(2,nez,nz) equivalence (phi00,phi00_ri) real den00_ri(2,nez,nz) equivalence (den00,den00_ri) real upar00_ri(2,nez,nz) equivalence (upar00,upar00_ri) real tpar00_ri(2,nez,nz) equivalence (tpar00,tpar00_ri) real tperp00_ri(2,nez,nz) equivalence (tperp00,tperp00_ri) PRINT *, 'netCDF version ',NF_INQ_LIBVERS() c open netCDF file c 5/98 file is now opened by calling routine, ncid is id # c filename=runname(1:lrunname)//'.nc' c c 10 STATUS = nf_open(filename, 0, NCID) ! read only c IF (STATUS .NE. NF_NOERR) then c PRINT *, NF_STRERROR(STATUS) c print *, ' ' c print *, 'Enter filename (with extension):' c read *, filename c goto 10 c endif c first determine if data file was produced by the electromagnetic c version of the code by searching for the variable beta_e emflag=0 status=nf_inq_nvars(ncid,nvars) do i=1,nvars status=nf_inq_varname(ncid,i,varname) if (varname(1:6).eq.'beta_e') then emflag=1 print *, 'Reading *.nc file from electromagnetic version' print *, ' ' endif enddo if (emflag.eq.0) print *, * 'Reading *.nc file from electrostatic version' c get dimensions id's status=nf_inq_dimid(ncid,'ld',ld_dim) status=nf_inq_dimid(ncid,'md',md_dim) status=nf_inq_dimid(ncid,'nd',nd_dim) status=nf_inq_dimid(ncid,'nspecies',nspec_dim) status=nf_inq_dimid(ncid,'ri',ri_dim) status=nf_inq_dimid(ncid,'time',time_dim) status=nf_inq_dimid(ncid,'timef',timef_dim) ! fast timescale status=nf_inq_dimid(ncid,'nd4',nd4_dim) status=nf_inq_dimid(ncid,'md4',md4_dim) status=nf_inq_dimid(ncid,'id',id_dim) ! drive term dim status=nf_inq_dimid(ncid,'char30',char30_dim) status=nf_inq_dimid(ncid,'char80',char80_dim) c get variable id's status=nf_inq_varid(ncid,'ld',ld_id) status=nf_inq_varid(ncid,'ldb',ldb_id) status=nf_inq_varid(ncid,'kd',kd_id) status=nf_inq_varid(ncid,'nmin',nmin_id) status=nf_inq_varid(ncid,'nmax',nmax_id) status=nf_inq_varid(ncid,'nd',nd_id) status=nf_inq_varid(ncid,'lin',lin_id) status=nf_inq_varid(ncid,'kdpass',kdpass_id) status=nf_inq_varid(ncid,'nstp',nstp_id) status=nf_inq_varid(ncid,'ne',ne_id) status=nf_inq_varid(ncid,'nfreq',nfreq_id) status=nf_inq_varid(ncid,'ntotal',ntotal_id) status=nf_inq_varid(ncid,'md',md_id) status=nf_inq_varid(ncid,'nspecies',nspecies_id) status=nf_inq_varid(ncid,'ikx',ikx_id) status=nf_inq_varid(ncid,'iperiod',iperiod_id) status=nf_inq_varid(ncid,'nparmom',nparmom_id) status=nf_inq_varid(ncid,'nperpmom',nperpmom_id) status=nf_inq_varid(ncid,'nemom',nemom_id) status=nf_inq_varid(ncid,'iodd',iodd_id) status=nf_inq_varid(ncid,'iflr',iflr_id) status=nf_inq_varid(ncid,'iphi00',iphi00_id) status=nf_inq_varid(ncid,'ifilter',ifilter_id) status=nf_inq_varid(ncid,'iexp',iexp_id) status=nf_inq_varid(ncid,'igradon',igradon_id) status=nf_inq_varid(ncid,'mmin',mmin_id) status=nf_inq_varid(ncid,'mmax',mmax_id) status=nf_inq_varid(ncid,'tim',tim_id) status=nf_inq_varid(ncid,'shr',shr_id) status=nf_inq_varid(ncid,'qsf',qsf_id) status=nf_inq_varid(ncid,'epsn',epsn_id) status=nf_inq_varid(ncid,'eps',eps_id) status=nf_inq_varid(ncid,'epse',epse_id) status=nf_inq_varid(ncid,'nueeff',nueeff_id) status=nf_inq_varid(ncid,'alpha',alpha_id) status=nf_inq_varid(ncid,'etai',etai_id) status=nf_inq_varid(ncid,'nuii',nuii_id) status=nf_inq_varid(ncid,'rmu1',rmu1_id) status=nf_inq_varid(ncid,'tiovte',tiovte_id) status=nf_inq_varid(ncid,'vy',vy_id) status=nf_inq_varid(ncid,'vz',vz_id) status=nf_inq_varid(ncid,'rmime',rmime_id) status=nf_inq_varid(ncid,'etae',etae_id) status=nf_inq_varid(ncid,'dt0',dt0_id) status=nf_inq_varid(ncid,'dt1',dt1_id) status=nf_inq_varid(ncid,'x0',x0_id) status=nf_inq_varid(ncid,'y0',y0_id) status=nf_inq_varid(ncid,'z0',z0_id) status=nf_inq_varid(ncid,'xp',xp_id) status=nf_inq_varid(ncid,'potential',potential_id) status=nf_inq_varid(ncid,'density',density_id) status=nf_inq_varid(ncid,'u_par',u_par_id) status=nf_inq_varid(ncid,'t_par',t_par_id) status=nf_inq_varid(ncid,'t_perp',t_perp_id) status=nf_inq_varid(ncid,'q_par',q_par_id) status=nf_inq_varid(ncid,'q_perp',q_perp_id) status=nf_inq_varid(ncid,'phiav_old',phiav_old_id) status=nf_inq_varid(ncid,'timo',timo_id) status=nf_inq_varid(ncid,'wpfx',wpfx_id) status=nf_inq_varid(ncid,'fluxi',fluxi_id) status=nf_inq_varid(ncid,'fluxe',fluxe_id) status=nf_inq_varid(ncid,'qfluxe',qfluxe_id) status=nf_inq_varid(ncid,'gamx',gamx_id) status=nf_inq_varid(ncid,'mgamx',mgamx_id) status=nf_inq_varid(ncid,'wenx',wenx_id) status=nf_inq_varid(ncid,'pcerr',pcerr_id) status=nf_inq_varid(ncid,'wtif',wtif_id) status=nf_inq_varid(ncid,'wkif',wkif_id) status=nf_inq_varid(ncid,'wpif',wpif_id) status=nf_inq_varid(ncid,'psp',psp_id) status=nf_inq_varid(ncid,'wenrk',wenrk_id) status=nf_inq_varid(ncid,'dketa',dketa_id) status=nf_inq_varid(ncid,'dkpar',dkpar_id) status=nf_inq_varid(ncid,'dktor',dktor_id) status=nf_inq_varid(ncid,'dktdp',dktdp_id) status=nf_inq_varid(ncid,'wkups',wkups_id) status=nf_inq_varid(ncid,'dktot',dktot_id) status=nf_inq_varid(ncid,'cx',cx_id) status=nf_inq_varid(ncid,'cy',cy_id) status=nf_inq_varid(ncid,'cz',cz_id) status=nf_inq_varid(ncid,'czc',czc_id) status=nf_inq_varid(ncid,'czcn',czcn_id) status=nf_inq_varid(ncid,'czn',czn_id) status=nf_inq_varid(ncid,'cznn',cznn_id) status=nf_inq_varid(ncid,'ct',ct_id) status=nf_inq_varid(ncid,'ninterv',ninterv_id) status=nf_inq_varid(ncid,'phi00',phi00_id) status=nf_inq_varid(ncid,'den00',den00_id) status=nf_inq_varid(ncid,'upar00',upar00_id) status=nf_inq_varid(ncid,'tpar00',tpar00_id) status=nf_inq_varid(ncid,'tperp00',tperp00_id) status=nf_inq_varid(ncid,'drt',drt_id) status=nf_inq_varid(ncid,'phirms',phirms_id) status=nf_inq_varid(ncid,'wakx',wakx_id) status=nf_inq_varid(ncid,'waky',waky_id) status=nf_inq_varid(ncid,'wxsp',wxsp_id) status=nf_inq_varid(ncid,'utor',utor_id) status=nf_inq_varid(ncid,'upol',upol_id) status=nf_inq_varid(ncid,'grtmx',grtmx_id) status=nf_inq_varid(ncid,'utim',utim_id) status=nf_inq_varid(ncid,'timf',timf_id) status=nf_inq_varid(ncid,'avgflag',avgflag_id) status=nf_inq_varid(ncid,'date',date_id) status=nf_inq_varid(ncid,'note',note_id) status=nf_inq_varid(ncid,'rmass',rmass_id) status=nf_inq_varid(ncid,'charge',charge_id) status=nf_inq_varid(ncid,'n_I',n_I_id) status=nf_inq_varid(ncid,'Ln',Ln_id) status=nf_inq_varid(ncid,'eta',eta_id) status=nf_inq_varid(ncid,'eta_par',eta_par_id) status=nf_inq_varid(ncid,'tau',tau_id) status=nf_inq_varid(ncid,'phisq',phisq_id) status=nf_inq_varid(ncid,'fluximn',fluximn_id) status=nf_inq_varid(ncid,'qfluximn',qfluximn_id) status=nf_inq_varid(ncid,'fluxemn',fluxemn_id) status=nf_inq_varid(ncid,'qfluxemn',qfluxemn_id) status=nf_inq_varid(ncid,'icrit',icrit_id) status=nf_inq_varid(ncid,'ntimf',ntimf_id) status=nf_inq_varid(ncid,'time',time_id) status=nf_inq_varid(ncid,'dt',dt_id) status=nf_inq_varid(ncid,'rkperp2',rkperp2_id) c read 0 and 1 dimensional variables status=nf_get_var_int(ncid,ld_id,ld) status=nf_get_var_int(ncid,ldb_id,ldb) status=nf_get_var_int(ncid,kd_id,kd) status=nf_get_var_int(ncid,nmin_id,nmin) status=nf_get_var_int(ncid,nmax_id,nmax) status=nf_get_var_int(ncid,nd_id,nd) status=nf_get_var_int(ncid,lin_id,lin) status=nf_get_var_int(ncid,kdpass_id,kdpass) status=nf_get_var_int(ncid,nstp_id,nstp) status=nf_get_var_int(ncid,ne_id,ne) status=nf_get_var_int(ncid,nfreq_id,nfreq) status=nf_get_var_int(ncid,ntotal_id,ntotal) status=nf_get_var_int(ncid,md_id,md) status=nf_get_var_int(ncid,nspecies_id,nspecies) status=nf_get_var_int(ncid,ikx_id,ikx) status=nf_get_var_int(ncid,iperiod_id,iperiod) status=nf_get_var_int(ncid,nparmom_id,nparmom) status=nf_get_var_int(ncid,nperpmom_id,nperpmom) status=nf_get_var_int(ncid,nemom_id,nemom) status=nf_get_var_int(ncid,iodd_id,iodd) status=nf_get_var_int(ncid,iflr_id,iflr) status=nf_get_var_int(ncid,iphi00_id,iphi00) status=nf_get_var_int(ncid,ifilter_id,ifilter) status=nf_get_var_int(ncid,iexp_id,iexp) status=nf_get_var_int(ncid,igradon_id,igradon) status=nf_get_var_int(ncid,mmin_id,mmin) status=nf_get_var_int(ncid,mmax_id,mmax) c values for time indices: nt1=ntim+1 nt2=ntim+ne ntf1=ntimr+1 ntf2=ntimr+ntotal status=nf_get_var_double(ncid,tim_id,tim) status=nf_get_var_double(ncid,shr_id,shr) status=nf_get_var_double(ncid,qsf_id,qsf) status=nf_get_var_double(ncid,epsn_id,epsn) status=nf_get_var_double(ncid,eps_id,eps) status=nf_get_var_double(ncid,epse_id,epse) status=nf_get_var_double(ncid,nueeff_id,nueeff) status=nf_get_var_double(ncid,alpha_id,alpha) status=nf_get_var_double(ncid,etai_id,etai) status=nf_get_var_double(ncid,nuii_id,nuii) status=nf_get_var_double(ncid,rmu1_id,rmu1) status=nf_get_var_double(ncid,tiovte_id,tiovte) status=nf_get_var_double(ncid,vy_id,vy) status=nf_get_var_double(ncid,vz_id,vz) status=nf_get_var_double(ncid,rmime_id,rmime) status=nf_get_var_double(ncid,etae_id,etae) status=nf_get_var_double(ncid,dt0_id,dt0) status=nf_get_var_double(ncid,dt1_id,dt1) status=nf_get_var_double(ncid,x0_id,x0) status=nf_get_var_double(ncid,y0_id,y0) status=nf_get_var_double(ncid,z0_id,z0) status=nf_get_var_double(ncid,xp_id,xp) status=nf_get_var_double(ncid,timo_id,timo(nt1:nt2)) status=nf_get_var_double(ncid,fluxe_id,fluxe(nt1:nt2)) status=nf_get_var_double(ncid,qfluxe_id,qfluxe(nt1:nt2)) status=nf_get_var_double(ncid,gamx_id,gamx(nt1:nt2)) status=nf_get_var_double(ncid,wenx_id,wenx(nt1:nt2)) status=nf_get_var_double(ncid,pcerr_id,pcerr(nt1:nt2)) status=nf_get_var_double(ncid,ct_id,ct(nt1:nt2)) status=nf_get_var_int(ncid,ninterv_id,ninterv) status=nf_get_var_double(ncid,phirms_id,phirms(nt1:nt2)) status=nf_get_var_double(ncid,wakx_id,wakx(nt1:nt2)) status=nf_get_var_double(ncid,waky_id,waky(nt1:nt2)) status=nf_get_var_double(ncid,wxsp_id,wxsp(nt1:nt2)) status=nf_get_var_double(ncid,utor_id,utor(nt1:nt2)) status=nf_get_var_double(ncid,upol_id,upol(nt1:nt2)) status=nf_get_var_double(ncid,timf_id,timf(ntf1:ntf2)) status=nf_get_var_double(ncid,avgflag_id,avgflag) status=nf_get_var_text(ncid,date_id,date) status=nf_get_var_text(ncid,note_id,note) status=nf_get_var_double(ncid,rmass_id,rmass) status=nf_get_var_double(ncid,charge_id,charge) status=nf_get_var_double(ncid,n_I_id,n_I) status=nf_get_var_double(ncid,Ln_id,Ln) status=nf_get_var_double(ncid,eta_id,eta) status=nf_get_var_double(ncid,eta_par_id,eta_par) status=nf_get_var_double(ncid,tau_id,tau) status=nf_get_var_int(ncid,icrit_id,icrit) status=nf_get_var_int(ncid,ntimf_id,ntimf) status=nf_get_var_double(ncid,time_id,time) status=nf_get_var_double(ncid,dt_id,dt) c read multi-dimensional arrays n1=nspecies status=nf_get_var_double(ncid,phiav_old_id, * phiav_old_ri(1:2,1:md,1:nd)) status=nf_get_var_double(ncid,potential_id, * potential_ri(1:2,1:ld,1:md,1:nd)) status=nf_get_var_double(ncid,rkperp2_id, * rkperp2(1:ld,1:md,1:nd)) status=nf_get_var_double(ncid,density_id, * density_ri(1:2,1:ld,1:md,1:nd,1:n1)) status=nf_get_var_double(ncid,u_par_id, * u_par_ri(1:2,1:ld,1:md,1:nd,1:n1)) status=nf_get_var_double(ncid,t_par_id, * t_par_ri(1:2,1:ld,1:md,1:nd,1:n1)) status=nf_get_var_double(ncid,t_perp_id, * t_perp_ri(1:2,1:ld,1:md,1:nd,1:n1)) status=nf_get_var_double(ncid,q_par_id, * q_par_ri(1:2,1:ld,1:md,1:nd,1:n1)) status=nf_get_var_double(ncid,q_perp_id, * q_perp_ri(1:2,1:ld,1:md,1:nd,1:n1)) c must be careful with the time variable to allow for calls with c ntim or ntimr nonzero. the netCDF values will be 1...ne, but c they must be read into nt1...nt2. status=nf_get_var_double(ncid,wpfx_id, * wpfx(nt1:nt2,1:nspecies)) status=nf_get_var_double(ncid,fluxi_id, * fluxi(nt1:nt2,1:nspecies)) status=nf_get_var_double(ncid,cx_id, * cx(nt1:nt2,1:4*nd)) status=nf_get_var_double(ncid,cy_id, * cy(nt1:nt2,1:4*md)) status=nf_get_var_double(ncid,cz_id, * cz(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,czc_id, * czc(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,czcn_id, * czcn(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,czn_id, * czn(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,cznn_id, * cznn(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,drt_id, * drt(1:14,nt1:nt2)) status=nf_get_var_double(ncid,phi00_id, * phi00_ri(1:2,nt1:nt2,1:nd)) status=nf_get_var_double(ncid,den00_id, * den00_ri(1:2,nt1:nt2,1:nd)) status=nf_get_var_double(ncid,upar00_id, * upar00_ri(1:2,nt1:nt2,1:nd)) status=nf_get_var_double(ncid,tpar00_id, * tpar00_ri(1:2,nt1:nt2,1:nd)) status=nf_get_var_double(ncid,tperp00_id, * tperp00_ri(1:2,nt1:nt2,1:nd)) status=nf_get_var_double(ncid,mgamx_id, * mgamx(1:md,1:nd,nt1:nt2)) status=nf_get_var_double(ncid,wtif_id, * wtif(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,wkif_id, * wkif(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,wpif_id, * wpif(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,psp_id, * psp(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,wenrk_id, * wenrk(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dketa_id, * dketa(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dkpar_id, * dkpar(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dktor_id, * dktor(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dktdp_id, * dktdp(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,wkups_id, * wkups(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dktot_id, * dktot(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,grtmx_id, * grtmx(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,phisq_id, * phisq(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,fluxemn_id, * fluxemn(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,qfluxemn_id, * qfluxemn(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,fluximn_id, * fluximn(nt1:nt2,1:md,1:nd,1:n1)) status=nf_get_var_double(ncid,qfluximn_id, * qfluximn(nt1:nt2,1:md,1:nd,1:n1)) status=nf_get_var_double(ncid,utim_id, * utim_ri(1:2,ntf1:ntf2,1:md,1:nd)) ! set up bounce-avgd electrons only if they exist ! (avoids zero size dimension and saves disk space) if (epse.gt.0.0) then status=nf_inq_dimid(ncid,'kdpass',kdpass_dim) status=nf_inq_varid(ncid,'e_density',e_density_id) status=nf_inq_varid(ncid,'e_p',e_p_id) status=nf_inq_varid(ncid,'e_r',e_r_id) status=nf_inq_varid(ncid,'e_t',e_t_id) status=nf_inq_varid(ncid,'phi_ba',phi_ba_id) status=nf_inq_varid(ncid,'phi_bk',phi_bk_id) status=nf_inq_varid(ncid,'e_denk',e_denk_id) status=nf_get_var_double(ncid,e_density_id, * e_density_ri(1:2,1:kdpass,1:md,1:nd)) status=nf_get_var_double(ncid,e_p_id, * e_p_ri(1:2,1:kdpass,1:md,1:nd)) status=nf_get_var_double(ncid,e_r_id, * e_r_ri(1:2,1:kdpass,1:md,1:nd)) status=nf_get_var_double(ncid,e_t_id, * e_t_ri(1:2,1:kdpass,1:md,1:nd)) status=nf_get_var_double(ncid,phi_ba_id, * phi_ba_ri(1:2,1:kdpass,1:md,1:nd)) status=nf_get_var_double(ncid,phi_bk_id, * phi_bk_ri(1:2,1:ld,1:md,1:nd)) status=nf_get_var_double(ncid,e_denk_id, * e_denk_ri(1:2,1:ld,1:md,1:nd)) endif c increment ntim and ntimr (for concatenation etc.) ntim=ntim+ne ntimr=ntimr+ntotal c close file (now done by the calling routine) c STATUS = NF_CLOSE(NCID) c IF (STATUS .NE. NF_NOERR) STOP return END cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc SUBROUTINE HANDLE_ERR(STATUS) INCLUDE 'netcdf.inc' INTEGER STATUS IF (STATUS .NE. NF_NOERR) THEN PRINT *, NF_STRERROR(STATUS) c STOP 'Stopped' ENDIF END