20040224: Hyperslab I/O of Varying Size: use of IMAP argument

Ben,

>Date: Tue, 24 Feb 2004 09:47:14 -0600
>From: Ben Howell <address@hidden>
>Organization: University of Wisconsin, Madison/SSEC
>To: Steve Emmerson <address@hidden>
>Subject: Re: 20040220: Hyperslab I/O of Varying Size: use of IMAP argument
> Keywords: 200402181846.i1IIkk2N001457

The above message contained the following:

> Thank you very much for your suggestions.
> I extracted the lines from your message and  put them into my code.
> The latitudes still have only half, i.e. 125 of 250, values followed by 
> 125 zeros.

The code doesn't look obviously wrong.

Are the values being put into the "slat" vector correctly?  What is the
value of the "nlines" variable?

> The longitudes however were written correctly (600 values).
> I have attached my FORTRAN code which contains, toward the end, your 
> suggested
> instructions.  If you still have some patience (and time) let me know if 
> there's any hope.
> 
> Thanks!
> :-)
> ------------------------------------------------------------------------------------------------------------------------
...
>       program MAIN
> 
>       implicit none
>       
>       include 'netcdf.inc'
> 
>       character*31 name 
>       character*80 stamp, ifile, logfile, msg
>       character*256 cbuf
>       integer*4 lentrim, IDCDF, IRET, nlines, nelems, nbands
>       integer*4 imageDate, imageTime, i, j, k, n
>       integer*4 maxbands, maxlines, maxelems, start(3), count(3)
>       
>       
> c ... 23Feb04:      
>       parameter(maxbands=1,maxlines=500, maxelems=1200)
> ccccc      parameter(maxbands=1,maxlines=250, maxelems=600)
>       
>       real*4 xdata(maxelems, maxlines, maxbands)
>       
> c ... Need 2 arrays of brightness temps, one for each output configuration:   
>    
>       real*4 btemp1(maxelems, maxlines)
>       real*4 btemp2(maxlines, maxelems)
>       
>       real*4 dlat(maxelems,maxlines), dlon(maxelems,maxlines)
>       real*4 slat(maxlines), slon(maxelems)
>       integer*4 idlinesDim, idelemsDim, idbandsDim
>       integer*4 idimageDate, idimageTime, idxdata, idlat, idlon
>       integer*4 iyear, idayoy
>       integer*4 nyrday, kyear, kmon, kalday, log
>       character*3 mon
>       real*4 y, ymin, yavg, ymax, hours, hoursit, fmnts
>       logical found
>       integer*4 ihour, mmss, imnts, isecs
> 
> c ... for output file:
>       integer*4 idims(2), ndims, latdim, londim
>       integer*4 idbt, idbch, idrch, idtl, idtlu, idcldy
>       integer*4 iddate, idtime, idate, itime
>       character*80 ofile, comnts
>                         
>       INTEGER*4 IYEREF,IMOREF,IDAREF,IHOREF,IMIREF,ISEREF
> C ... DABTIME RETURNS ELAPSED TIME IN DAYS (DOUBLE PRECISION)
> C ... EDBASE = ELAPSED TIME OF OUR BASE TIME SINCE "EARLY REF."
> C ... EDNOW = ELAPSE TIME (DAYS) FROM "EARLY REF." UNTIL SPECIFIED DATE & 
> TIME.
>       REAL*8 DABTIM, EDREF, EDNOW, DVALU
> 
>       REAL*4 ABSENT
>       INTEGER*4 JDAY, NTHDAY, IOPTDATE, IOPTOUT, NLATS, NLONS
> 
> c ... 17Feb04:
>       REAL*4 SIGN
> 
> C ... 23FEB04:
>       INTEGER*4 STRIDE(2), IMAP(2)
>       DATA STRIDE/1,1/
>                   
>       DATA ABSENT/-9999.0/
>             
> C ... REFERENCE DATE AND TIME:
>       DATA IYEREF/1970/, IMOREF/1/, IDAREF/1/, IHOREF/0/, IMIREF/0/,
>      *ISEREF/0/
> 
> c ---------------------------------------------------------------------       
>      
>       data stamp/
>      *'$Id: sat_globalbt.f,v 1.2(modified) 2004/02/23 $'/
> c ---------------------------------------------------------------------
>             
>       data start/1,1,1/, log/13/, ndims/2/
>       
> c ***
> 
>       cbuf=stamp
>       n=lentrim(cbuf)
>       write(*,*)' '
>       write(*,*)stamp(1:n)
>       write(*,*)' '
> 
>       write(*,*)'Log (output, text) file:'
>       read(*,'(a)')logfile
>                        
>       write(*,*)'Input (MCIDAS 1-band NetCDF)file:'
>       read(*,'(a)')ifile
>       cbuf=ifile
>       n=lentrim(cbuf)
>       inquire(file=ifile, exist=found)
>       if(found)then
>        write(msg,*)'Input file:', cbuf(1:n)
>        call logmsg(log,logfile,msg)
>       else
>          msg='Input file not found:'
>        call logmsg(log, logfile, msg)
>        write(msg,*)'File:',cbuf(1:n)
>        call logmsg(log, logfile, msg)
>          stop 'Input file?'
>       endif
> 
>       write(*,*)'Output (NetCDF) file:'
>       read(*,'(a)')ofile
>       
>       write(*,*)'Choose format for date:'
>       write(*,*)'   1. YYYYDDD'
>       write(*,*)'   2. YYYYMMDD'
>       read(*,*)ioptdate
>             
>       write(*,*)'Choose configuration of output:'
>       write(*,*)'   1. row-wise (latitude, longitude)'
>       write(*,*)'   2. column-wise (longitude, latitude)'
>       read(*,*)ioptout      
>       
>       write(*,*)'Comments(up to 80 characters):'
>       read(*,'(a)')comnts
>       
> c  Make errors non-fatal to allow us to do error handling.
> c  The default is all errors are fatal = call ncpopt(NCVERBOS+NCFATAL)
> c
>       call NCPOPT(0)
> 
> c  Open the netCDF file as read-only NCNOWRIT; get its "ID" IDCDF:
> c
>       cbuf=ifile
>       n=lentrim(cbuf)
>       ifile(n+1:n+1)=char(0)
>       IDCDF = NCOPN(ifile, NCNOWRIT, IRET)
>       if (IRET.ne.0)stop 'Error opening file'
> 
> c ... get dimension IDs from their names:
>       idlinesDim = NCDID(IDCDF, 'lines', IRET)
>       if(IRET.ne.0)stop 'Error, lines dimension'
>       
>       idelemsDim = NCDID(IDCDF, 'elems', IRET)
>       if(IRET.ne.0)stop 'Error, elems dimension'
>       
>       idbandsDim = NCDID(IDCDF, 'bands', IRET)
>       if(IRET.ne.0)stop 'Error, bands dimension'
>       
> c ... Get the size of each dimension
>       call NCDINQ(IDCDF, idlinesDim, name, nlines, IRET)
>       if ( IRET.NE.0 ) stop 'Error getting size of lines dimension'
>       cbuf=name
>       n=lentrim(cbuf)
>       write(msg,*)'dimension name:', name(1:n)
>       call logmsg(log, logfile, msg)
>       write(msg,*)'dimension size:', nlines
>       call logmsg(log, logfile, msg)
>                   
>       call NCDINQ(IDCDF, idelemsDim, name, nelems, IRET)
>       if ( IRET.NE.0 ) stop 'Error getting sizse of elems dimension'
>       cbuf=name
>       n=lentrim(cbuf)
>       write(msg,*)'dimension name:', name(1:n)
>       call logmsg(log, logfile, msg)
>       write(msg,*)'dimension size:', nelems
>       call logmsg(log, logfile, msg)
>       
>       call NCDINQ(IDCDF, idbandsDim, name, nbands, IRET)
>       if ( IRET.NE.0 ) stop 'Error getting size of bands dimension'
>       cbuf=name
>       n=lentrim(cbuf)
>       write(msg,*)'dimension name:', name(1:n)
>       call logmsg(log, logfile, msg)
>       write(msg,*)'dimension size:', nbands
>       call logmsg(log, logfile, msg)
> 
> c ... Check size of "band" dimension:
>       if(nbands.ne.1)then
>          write(msg,*)'nbands:', nbands
>        call logmsg(log, logfile, msg)
>          msg='Wrong band dimension.  Expect 1'
>        call logmsg(log, logfile, msg)
>        stop 'Bands?'
>       endif
>             
> c ... get variable IDs from their names      
>       idimageDate = NCVID(IDCDF, 'imageDate', IRET)
>       if(IRET.ne.0)stop 'Error, imageDate ID'
>       
>       idimageTime = NCVID(IDCDF, 'imageTime', IRET)
>       if(IRET.ne.0)stop 'Error, imageTime ID'
>       
>       idlat = NCVID(IDCDF, 'latitude', IRET)
>       if(IRET.ne.0)stop 'Error, latitude ID'
>       
>       idlon = NCVID(IDCDF, 'longitude', IRET)
>       if(IRET.ne.0)stop 'Error, longitude ID'
> 
>       idxdata = NCVID(IDCDF, 'data', IRET)
>       if(IRET.ne.0)stop 'Error, data ID'
> 
> c ... Set "count" vector
>       count(1)=nelems
>       count(2)=nlines
>       count(3)=nbands
>        
> c  Get data (point values first, then time series)
>       call NCVGT(IDCDF, idimageDate, 1, 1, imageDate, IRET)
>       call NCVGT(IDCDF, idimageTime, 1, 1, imageTime, IRET)
>       
>       write(*,*)'imageDate, imageTime:',imageDate, imageTime
>       write(msg,*)'imageDate, imageTime:',imageDate, imageTime
>       call logmsg( log, logfile, msg)
>       
>       iyear=imageDate/1000
>       idayoy=imageDate-1000*iyear
>       iyear=iyear+1900
>       write(*,*)'year, day-of-year:', iyear, idayoy
>       write(msg,*)'year, day-of-year:', iyear, idayoy
>       call logmsg( log, logfile, msg)
>       
> c ... get calendar month and day
>       nyrday=idayoy + 1000*iyear
>       call calday(nyrday, kyear, kmon, kalday, mon)
>       write(*,*)'Calendar year, month, day:', kyear, kmon, kalday      
>       write(msg,*)'Calendar year, month, day:', kyear, kmon, kalday      
>       call logmsg(log, logfile, msg)
>       
>       if(ioptdate.eq.1)then
> c ...    Save date in YYYYDDD format:      
>          jday=nthday(kyear,kmon,kalday)
>          idate=kyear*1000+jday
>       else
> c ...    Save date in YYYYMMDD format:
>          idate=kyear*10000 + kmon*100 + kalday
>       endif      
>              
>             
> 
> c ... Get hours, minutes and seconds (image time):     
>       ihour=imageTime/10000
>       mmss=imageTime-10000*ihour
>       imnts=mmss/100
>       isecs=mmss-100*imnts
>       write(msg,*)'Image time, hour, minute, second:',
>      *             ihour, imnts, isecs
>       call logmsg(log,logfile,msg)
>       
>       hoursit=float(ihour)+float(imnts)/float(60)
>      *        +float(isecs)/float(3600)
>       write(*,*)'Image time, decimal hours:', hoursit
>       write(msg,*)'Image time, decimal hours:', hoursit
>       call logmsg(log,logfile,msg)
>       
> c ... Get "dlat" ...
>       call NCVGT(IDCDF, idlat, start, count, dlat, IRET)
>       if(IRET.ne.0)stop 'Error getting dlat'
>        
> c ... Get "dlon" ...
>       call NCVGT(IDCDF, idlon, start, count, dlon, IRET)
>       if(IRET.ne.0)stop 'Error getting dlon'
>        
> c ... Store unique latitudes:
>       do j=1,nlines
>        slat(j)=dlat(1,j)
>       enddo
>       nlats=nlines
>       
> c ... Correction of last latitude
>       if(slat(nlines) .lt. -90.0 .or. slat(nlines) .gt. 90.0)then
>          write(msg,*)'slat(nlines), before correction:', slat(nlines)
>        call logmsg(log, logfile, msg)
>          slat(nlines)= -90.0
>          write(msg,*)'lat(nlines), after correction:', slat(nlines)
>        call logmsg(log, logfile, msg)
>       endif                  
>                
> c ... Store unique longitudes:
>       do i=1,nelems
>       
> c ...    17Feb04, input file "GOES1993.nc" appears to have longitudes
> c                 as degrees east, so no conversion is necessary.
> c                 To be safe, check first longitude
>          if(i.eq.1)then
>           if(dlon(1,1).lt.0.0)then
>              sign=1.0
>           else
>              sign= -1.0
>           endif
>        endif               
>          slon(i)=dlon(i,1)*sign
>        
>       enddo
>       nlons=nelems
>                        
> c ... Ingest "data" variable into "xdata" ...       
>       call NCVGT(IDCDF, idxdata, start, count, xdata, IRET)
>       if(IRET.ne.0)stop 'Error getting xdata'
> 
> c  Close the netCDF file
>       call NCCLOS(IDCDF, IRET)
>       if(IRET.ne.0)stop'Error closing file'
> 
>          k=1
>          do i=1,nelems
>           do j=1,nlines
> 
> c ...          18Feb04, check for "missing" data:
>                if(nint(xdata(i,j,k)) .eq. 330)then
>                 xdata(i,j,k)=-9999.0
>              endif
>                   
>                if(ioptout.eq.1)then       
>                   btemp1(i,j)=xdata(i,j,k)
>              else
>                   btemp2(j,i)=xdata(i,j,k)
>              endif
>                 
>           enddo
>        enddo
>        write(*,*)'Stored brightness temp.'
>        
> c ...    Get min, mean, and max for each latitiude(line)
>          msg='Brightness temperature, lat, min, mean, max:'
>        call logmsg(log, logfile, msg)
>        write(*,*)msg
>          do j=1,nlines,50
>           ymin=1.0e+10
>           ymax=-ymin
>           yavg=0.
>           do i=1,nelems
>              if(ioptout.eq.1)then
>                 y=btemp1(i,j)
>              else
>                 y=btemp2(j,i)
>              endif      
>              yavg=yavg+y
>              if(y.lt.ymin)ymin=y
>              if(y.gt.ymax)ymax=y
>           enddo
>           yavg=yavg/float(nelems)
>           write(*,'(3x, f5.1, 3f7.1)')
>      *      slat(j),ymin,yavg,ymax
>           write(msg,'(3x, f5.1, 3f7.1)')
>      *      slat(j),ymin,yavg,ymax
>           call logmsg(log, logfile, msg)
>        enddo
>        
> CCCCCCCCCCC  END OF INPUT SECTION, START OUTPUT SECTION CCCCCCCCC
> 
> c ... Create output (NetCDF) file ...
>       IDCDF=nccre(ofile,NCCLOB,IRET)
>       if(IRET.ne.0)then
>          cbuf=ofile
>        n=lentrim(cbuf)
>          write(*,*)'Unable to open file:', cbuf(1:n)
>          write(msg,*)'Unable to open file:', cbuf(1:n)
>        call logmsg(log, logfile, msg)
>       endif
>       
> C ... DEFINE DIMENSIONS:
> 
>       LATDIM=NCDDEF(IDCDF, 'latitude', nlats,  IRET)
>       LONDIM=NCDDEF(IDCDF, 'longitude', nlons,  IRET)
> 
> C ... DEFINE VARIABLE IDS:      
>       IDDATE=NCVDEF(IDCDF, 'date', NCLONG, 0, 0, IRET)
>       IDTIME=NCVDEF(IDCDF, 'time', NCLONG, 0, 0, IRET)
>       
>       IDLAT=NCVDEF(IDCDF, 'latitude', NCFLOAT, 1, LATDIM, IRET)
>       IDLON=NCVDEF(IDCDF, 'longitude', NCFLOAT, 1, LONDIM, IRET)
>       
>       if(ioptout.eq.1)then 
>          IDIMS(1)=LONDIM
>          IDIMS(2)=LATDIM
>       else     
>          IDIMS(1)=LATDIM
>          IDIMS(2)=LONDIM
>       endif
>       
>       IDBT=NCVDEF(IDCDF, 'brightness_temp',NCFLOAT,NDIMS,IDIMS,IRET)
>       IF(IRET.NE.0)STOP 'IDBT?'
>       
> 
> C ... ASSIGN ATTRIBUTES TO VARIABLES
> 
> c ... Date:
> 
>       if(ioptdate.eq.1)then
>          CALL NCAPTC(IDCDF, IDDATE, 'long_name', NCCHAR, 29,
>      *'   Image date, in YYYYDDD format', IRET)     
>       else
>          CALL NCAPTC(IDCDF, IDDATE, 'long_name', NCCHAR, 30,
>      *'   Image date, in YYYYMMDD format', IRET)     
>       endif
>       
> c ... Time:
>       CALL NCAPTC(IDCDF, IDTIME, 'long_name', NCCHAR, 15,
>      *'Image time, UTC', IRET)     
>       CALL NCAPTC(IDCDF,IDTIME,'units',NCCHAR,31,
>      *'seconds since 1970-1-1 00:00:00',IRET)
>       CALL NCAPT(IDCDF, IDTIME, 'image_time_HHMMSS_UTC', NCLONG, 1,
>      *imageTime, IRET) 
> 
> c ... Latitude:
>       CALL NCAPTC(IDCDF, IDLAT, 'long_name', NCCHAR, 8,
>      *'Latitude', IRET)
>       CALL NCAPTC(IDCDF, IDLAT, 'units', NCCHAR, 13,
>      *'degrees_north', IRET)
>       CALL NCAPTC(IDCDF, IDLAT, 'precision', NCCHAR, 4,
>      *'1E-1', IRET)
>                        
> c ... Longitude:
>       CALL NCAPTC(IDCDF, IDLON, 'long_name', NCCHAR, 9,
>      *'Longitude', IRET)
>       CALL NCAPTC(IDCDF, IDLON, 'units', NCCHAR, 12,
>      *'degrees_east', IRET)
>       CALL NCAPTC(IDCDF, IDLON, 'precision', NCCHAR, 4,
>      *'1E-1', IRET)
>                        
> c ... Brightness temperature
>       CALL NCAPTC(IDCDF, IDBT, 'long_name', NCCHAR, 32,
>      *'Radiative brightness temperature', IRET)
>       CALL NCAPTC(IDCDF, IDBT, 'units', NCCHAR, 14,
>      *'degrees_Kelvin', IRET)
>       CALL NCAPTC(IDCDF, IDBT, 'precision', NCCHAR, 4,
>      *'1E-1', IRET)
>       CALL NCAPT(IDCDF, IDBT, 'missing_value', NCFLOAT, 1,             
>      *ABSENT, IRET)
>                        
> C ... ADD USER COMMENTS TO "Global Attributes" ...
>       CBUF=COMNTS
>       N=LENTRIM(CBUF)
>       IF(N.GT.79)N=79
>       IF(N.LT.1)THEN
>          COMNTS='No comments'
>          N=11
>       ENDIF
>       CALL NCAPTC(IDCDF,NCGLOBAL,'Comments',NCCHAR,N,COMNTS,IRET)
>          
> C ... ADD SOFTWARE VERSION TO "Global Attributes" ...
>       CBUF=STAMP
>       N=LENTRIM(CBUF)
>       IF(N.GT.79)N=79
>       stamp(n+1:n+1)=char(0)
>       CALL NCAPTC(IDCDF,NCGLOBAL,'SoftwareVersion',NCCHAR,N,STAMP,
>      *            IRET)
> 
> C ... LEAVE DEFINE MODE ...
>       CALL NCENDF(IDCDF, IRET)
> 
> C ... COMPUTE ELAPSED DAYS FROM EARLY REFERENCE TO 1 JAN 1970
>       EDREF=DABTIM(IYEREF,IMOREF,IDAREF,IHOREF,IMIREF,ISEREF)
> 
> C ... GET BASE TIME (SECONDS ELAPSED SINCE 1 JAN 1970, 0000Z)
>       EDNOW=DABTIM(KYEAR,KMON,KALDAY,IHOUR,IMNTS,ISECS)
>       DVALU=(EDNOW-EDREF)*DFLOAT(86400)
>       ITIME=DNINT(DVALU)
>       
> 
>       write(*,*)'Writing data to NetCDF file...'
>             
>       START(1)=1
>       COUNT(1)=1
>       
> C ... STORE DATE ...
>       CALL NCVPT(IDCDF, IDDATE, START, COUNT, IDATE, IRET)
>       IF(IRET.EQ.0)THEN
>          WRITE(*,*)'Wrote date:', idate
>       ELSE
>          STOP 'ERROR WRITING DATE'
>       ENDIF    
>        
>        
> C ... STORE TIME ...
>       CALL NCVPT(IDCDF, IDTIME, START, COUNT, ITIME, IRET)
>       IF(IRET.EQ.0)THEN
>          WRITE(*,*)'Wrote image time:', itime
>       ELSE
>          STOP 'ERROR WRITING IMAGE TIME'
>       ENDIF    
> 
> C ... STORE LATITUDES:
> C ... 24FEB04, as suggested by Steve Emmerson(Unidata):
>        
>       
> 
> C ... OUTPUT LATITUDES:
>       START(1)=1
>       COUNT(1)=NLATS
> c     IMAP(1)=1
> 
> c ... temporary:
>       write(*,*)'latitudes: start(1), count(1):', start(1), count(1)
>       pause
>             
> c     IRET=NF_PUT_VARM_REAL(IDCDF,IDLAT,START,COUNT,  ! no
> c    *STRIDE,IMAP,SLAT)                               ! no
>       IRET=NF_PUT_VARA_REAL(IDCDF,IDLAT,START,COUNT,  ! yes
>      *SLAT)                                           ! yes
>  
> C ... OUTPUT LONGITUDES:
>       START(1)=1
>       COUNT(1)=NLONS
> c     IMAP(1)=1                                       ! no
>       
> c ... temporary:
>       write(*,*)'longitudes: start(1), count(1):', start(1), count(1)
>       pause
>             
> c     IRET=NF_PUT_VARM_REAL(IDCDF,IDLON,START,COUNT,  ! no
> c    *STRIDE,IMAP,SLON)                               ! no
>       IRET=NF_PUT_VARA_REAL(IDCDF,IDLON,START,COUNT,  ! yes
>      *SLON)                                           ! yes
>  
> CCCCC CALL NCVPT(IDCDF, IDLAT, START, COUNT, SLAT, IRET)
> 
>       IF(IRET.EQ.0)THEN
>          WRITE(*,*)'Wrote latitudes'
>       ELSE
>          STOP 'ERROR WRITING LATITUDES'
>       ENDIF    
>       
>       START(1)=1
>       COUNT(1)=NLONS
>       IMAP(1)=1
>       
> C ... STORE LONGITUDES:
> 
> C ... 24FEB04, as suggested by Steve Emmerson(Unidata):
>       IRET=NF_PUT_VARA_REAL(IDCDF,IDLON,START,COUNT,SLON)                     
>                       ! yes
> CCCCC CALL NCVPT(IDCDF, IDLON, START, COUNT, SLON, IRET)
> 
>       IF(IRET.EQ.0)THEN
>          WRITE(*,*)'Wrote longitudes'
>       ELSE
>          STOP 'ERROR WRITING LONGITUDES'
>       ENDIF    
>       
> C ... 11FEB04:      
>       IF(IOPTOUT.EQ.1)THEN
>          START(1)=1
>          COUNT(1)=NLONS
>          START(2)=1
>          COUNT(2)=NLATS
>        IMAP(1)=1
>        
> C ...    24FEB04, as suggested by Steve Emmerson(Unidata):
>          IMAP(2)=MAXELEMS       ! element count between 2nd index increments
>        
> C ...    23FEB04:
>          IRET=NF_PUT_VARM_REAL(IDCDF,IDBT,START,COUNT,
>      *   STRIDE,IMAP,BTEMP1)   
> CCCCC    CALL NCVPT(IDCDF, IDBT, START, COUNT, BTEMP1, IRET)
> 
>       ELSE
>          START(1)=1
>          COUNT(1)=NLATS
>          START(2)=1
>          COUNT(2)=NLONS
>        IMAP(1)=1
>        
> C ...    24FEB04, as suggested by Steve Emmerson(Unidata):
>          IMAP(2)=MAXLINES       ! element count between 2nd index increments
>        
> C ...    23FEB04:      
>          IRET=NF_PUT_VARM_REAL(IDCDF,IDBT,START,COUNT,
>      *   STRIDE,IMAP,BTEMP2)   
> CCCCC    CALL NCVPT(IDCDF, IDBT, START, COUNT, BTEMP2, IRET)
> 
>       ENDIF    
>       IF(IRET.EQ.0)THEN
>          WRITE(*,*)'Wrote brightness temps'
>       ELSE
>          STOP 'ERROR WRITING BR.TEMPS'
>       ENDIF    
>       
> 
> C ... CLOSE NETCDF FILE ..
>       CALL NCCLOS(IDCDF,IRET)
>       IF(IRET.NE.0)THEN
>          WRITE(*,*)'Bad status from NCCLOS:', IRET
>          STOP 'NCCLOS error'
>       ENDIF
>        
>       msg='*** End sat_globalbt ***'
>       call logmsg(log, logfile, msg)
>                    
>       stop '*** End sat_globalbt ***'
>       end

Regards,
Steve Emmerson