text:
c
c
	subroutine read_surface_data(i4time,btime,n_obs_g,n_obs_b,time,
     &    wmoid,stations,provider,wx,reptype,autostntype,lat,lon,elev,
     &    t,td,rh,dd,ff,ddg,ffg,alt,stnp,mslp,delpch,delp,vis,solar,
     &    sfct,sfcm,pcp1,pcp3,pcp6,pcp24,snow,kkk_s,max24t,min24t,t_ea,
     &    td_ea,rh_ea,dd_ea,ff_ea,alt_ea,p_ea,vis_ea,solar_ea,sfct_ea,
     &    sfcm_ea,pcp_ea,snow_ea,store_cldamt,store_cldht,
     &    maxsta,jstatus)
c
c*****************************************************************************
c
cdoc    Routine to read the newer LAPS LSO surface data file with the expanded
cdoc    variable list. The data is passed back to the calling routine in 1-d
cdoc    arrays.
c
c	Changes:
c		P. Stamus  04-13-98  Original version (from write_surface_LS2).
c                          05-01-98  Added soil moisture variables.
c               J. Edwards 09-16-98  moved all format definitions to
c                                    src/include/lso_formats.inc
c                                    changed 909 definition to allow for
c                                    missing data
c
c*****************************************************************************
c
	real lat(maxsta), lon(maxsta), elev(maxsta)
	real t(maxsta), t_ea(maxsta), max24t(maxsta), min24t(maxsta)
	real td(maxsta), td_ea(maxsta), rh(maxsta), rh_ea(maxsta)
	real dd(maxsta), ddg(maxsta), dd_ea(maxsta)
	real ff(maxsta), ffg(maxsta), ff_ea(maxsta)
	real alt(maxsta), alt_ea(maxsta), delp(maxsta)
	real stnp(maxsta), mslp(maxsta), p_ea(maxsta)
	real vis(maxsta), vis_ea(maxsta)
	real solar(maxsta), solar_ea(maxsta)
	real sfct(maxsta), sfct_ea(maxsta)
	real sfcm(maxsta), sfcm_ea(maxsta)
	real pcp1(maxsta), pcp3(maxsta), pcp6(maxsta), pcp24(maxsta)
	real snow(maxsta), snow_ea(maxsta), pcp_ea(maxsta)
	real store_cldht(maxsta,5)
c
	integer i4time, wmoid(maxsta), jstatus
	integer time(maxsta), delpch(maxsta), kkk_s(maxsta)
c
	character filetime*9, infile*256, btime*24
	character stations(maxsta)*20, provider(maxsta)*11
        character reptype(maxsta)*6, autostntype(maxsta)*6
	character wx(maxsta)*25, store_cldamt(maxsta,5)*4
c
c
c.....  Blank out the character arrays.
c
	jstatus = 0
	do i=1,maxsta
	   stations(i) = '                    '
	   provider(i) = '           '
	   reptype(i)  = '      '
	   autostntype(i)  = '      '
	   wx(i) = '                         '
	   do j=1,5
	      store_cldamt(i,j) = '    '
	   enddo !j
	enddo !i
c
c.....	Get the file.
c
	call make_fnam_lp(i4time, filetime, istatus)
	call get_directory('lso', infile, len)
	infile = infile(1:len) // filetime(1:9) // '.lso'
c
	open(11,file=infile,status='old',err=999)
c
c.....	Read the header.
c
	read(11,900) btime,		! time
     &               n_obs_g,		! # of obs in the laps grid
     &               n_obs_b		! # of obs in the box
c
c.....	Read the station data.
c
        if(n_obs_b.eq.0) then
         print*, 'no obs in file ',infile,' returning'
         jstatus=-1
         return
        endif
	do k=1,n_obs_b
c
	   read(11,901)  stations(k),            !station id
     &                   wmoid(k),               !WMO id number
     &                   provider(k),            !data provider
     &                   lat(k), lon(k), elev(k),!lat/lon (deg), elev (meters MSL)
     &                   time(k)		 !obs time
c
	  read(11,903)   reptype(k),             !station report type
     &                   autostntype(k),         !station type (manual/auto)
     &                   wx(k)                   !present weather
c
	  read(11,905)   t(k), t_ea(k),          !temp, temp expected accuracy (Deg F)
     &                   td(k), td_ea(k),        !dew point, dew point exp. accuracy (Deg F)
     &                   rh(k), rh_ea(k)         !Rel hum, rh expected accuracy (%)
c
	  read(11,907)   dd(k), ff(k),           !wind dir (deg), wind speed (knots)
     &                   ddg(k), ffg(k),         !wind gust dir (deg), wind gust speed (knots)
     &                   dd_ea(k), ff_ea(k)      !dir expected accuracy (deg), spd exp accuracy (knots)
c
	  read(11,909)   alt(k),                 !altimeter (mb)
     &                   stnp(k),                !station pressure (mb)
     &                   mslp(k),                !MSL pressure (mb)
     &                   delpch(k),              !3-h press change character (FMH-1 Manual, Sec 12.7.2)
     &                   delp(k),                !3-h pressure change (mb)
     &                   p_ea(k), alt_ea(k)      !pressure exp accuracy, alt exp accuracy
c
	  read(11,911)   vis(k), vis_ea(k),      !visibility, vis exp accuracy (miles)
     &                   solar(k), solar_ea(k),  !solar, solar exp accuracy (watts/meter**2)
     &                   sfct(k), sfct_ea(k),    !soil/water temp, soil/water temp exp accuracy (Deg F)
     &                   sfcm(k), sfcm_ea(k)     !soil moist, soil moist temp exp accuracy (%)
c
	  read(11,913)   pcp1(k),                !1-h precipitation (inches)
     &                   pcp3(k),                !3-h precipitation (inches)
     &                   pcp6(k),                !6-h precipitation (inches)
     &                   pcp24(k),               !24-h precipitation
     &                   snow(k),                !snow depth (inches)
     &                   pcp_ea(k), snow_ea(k)   !precip and snow exp accuracy (inches)
c
	  read(11,915)  kkk_s(k),                !num cld layers
     &                  max24t(k),               !24-h max temperature (Deg F)
     &                  min24t(k)                !24-h min temperature (Deg F)
c
c.....	Read the cloud data if we have any.
c
	  if(kkk_s(k) .gt. 0) then
	    do ii=1,kkk_s(k)
  	      read(11,917) store_cldamt(k,ii), store_cldht(k,ii) ! layer cloud amount and height
                                                                 ! (meters MSL)
	    enddo !ii
	  endif
c
	enddo !k
c
!       endfile(11)
	close(11)
c
c..... End of data reading.  Let's go home...
c
	jstatus = 1
	return
999     continue
        print *,'Could not open file: ',infile
	print *,' ++ ERROR opening LSO file in READ_SURFACE_DATA ++'
        jstatus = -1
	return
        include 'lso_formats.inc'
	end
c
c
	subroutine read_sfc_metadata(i4time,n_obs_g,n_obs_b,time
     1              ,wmoid,stations,provider,reptype,autostntype
     1              ,lat,lon,elev,maxsta,jstatus)
c
c*****************************************************************************
c
cdoc    Routine to read the LAPS LSO surface data file and return the station
cdoc    metadata.   The data is passed back to the calling routine in 1-d
cdoc    arrays.
c
c	Changes:
c		P. Stamus  04-13-98  Original version (from write_surface_LS2).
c
c               S. Albers  2000      Expanded for OSSE use
c
c*****************************************************************************
c
	real lat(maxsta), lon(maxsta), elev(maxsta)
c
	integer i4time, wmoid(maxsta), jstatus
        integer time(maxsta)                 ! integer representation of HHMM
c
	character filetime*9, infile*256, btime*24
	character stations(maxsta)*20, provider(maxsta)*11
        character reptype(maxsta)*6, autostntype(maxsta)*6
	character dum*132
c
	jstatus = 0
c
c.....	Get the file.
c
	call make_fnam_lp(i4time, filetime, istatus)
	call get_directory('lso', infile, len)
	infile = infile(1:len) // filetime(1:9) // '.lso'
c
	open(11,file=infile,status='old',err=999)
c
c.....	Read the header.
c
	read(11,900) btime,		! time
     &               n_obs_g,		! # of obs in the laps grid
     &               n_obs_b		! # of obs in the box
c
c.....	Read the station data.
c
	do k=1,n_obs_b
c
	   read(11,901)  stations(k),              !station id
     &                   wmoid(k),                 !WMO id number
     &                   provider(k),              !data provider
     &                   lat(k), lon(k), elev(k),  !lat, lon, elev
     &                   time(k)                   !obs time
c
	  read(11,903)   reptype(k),               !station report type
     &                   autostntype(k)            !station type (manual/auto)
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,915)  kkk_s,                     !num cld layers
     &                  dummy,                     !24-h max temperature
     &                  dummy                      !24-h min temperature
c
c.....	Read the cloud data if we have any.
c
	  if(kkk_s .gt. 0) then
	    do ii=1,kkk_s
	      read(11,919) dum
	    enddo !ii
	  endif
c
	enddo !k
c
!       endfile(11)
	close(11)
c
c..... End of data reading.  Let's go home...
c
	jstatus = 1
	return
999     continue
	print *,' ++ ERROR opening LSO file in READ_SFC_METADATA ++'
        jstatus = -1
	return
        include 'lso_formats.inc'
	end
c
c
	subroutine read_sfc_state(i4time,ext,btime,n_obs_g,n_obs_b,
     &    stations,provider,lat,lon,elev,t,td,rh,dd,ff,alt,stnp,mslp,
     &    maxsta,jstatus)
c
c*****************************************************************************
c
cdoc    Routine to read the LAPS LSO surface data file and return the state
cdoc    variables (wind, temp, dewpt, altimeter & pressure).   The data is
cdoc    passed back to the calling routine in 1-d arrays.
c
c	Changes:
c		P. Stamus  04-13-98  Original version (from write_surface_LS2).
c               J. Edwards 09-16-98  moved all format definitions to
c                                    src/include/lso_formats.inc
c                                    changed 909 definition to allow for
c                                    missing data
c
c*****************************************************************************
c
	real lat(maxsta), lon(maxsta), elev(maxsta)
	real t(maxsta)
	real td(maxsta)
	real rh(maxsta)
	real dd(maxsta)
	real ff(maxsta)
	real alt(maxsta)
	real stnp(maxsta), mslp(maxsta)
c
	integer i4time, jstatus
c
	character filetime*9, infile*256, btime*24
	character stations(maxsta)*20, provider(maxsta)*11
	character dum*132
        character ext*(*)
c
	jstatus = 0
c
c.....	Get the file.
c
	call make_fnam_lp(i4time, filetime, istatus)
	call get_directory('lso', infile, len)
        call s_len(ext,len_ext)
	infile = infile(1:len) // filetime(1:9) // '.' // ext(1:len_ext)
c
	open(11,file=infile,status='old',err=999)
c
c.....	Read the header.
c
	read(11,900) btime,		! time
     &               n_obs_g,		! # of obs in the laps grid
     &               n_obs_b		! # of obs in the box
c
c.....	Read the station data.
c
	do k=1,n_obs_b
c
	   read(11,901)  stations(k),              !station id
     &                   idummy,                   !WMO id number
     &                   provider(k),              !data provider
     &                   lat(k), lon(k), elev(k),  !lat, lon, elev
     &                   idummy                    !obs time
c
	  read(11,919) dum
c
	  read(11,905)   t(k), dummy,              !temp, temp expected accuracy
     &                   td(k), dummy,             !dew point, dew point exp. accuracy
     &                   rh(k), dummy              !Rel hum, rh expected accuracy
c
	  read(11,907)   dd(k), ff(k),             !wind dir, wind speed
     &                   dummy, dummy,             !wind gust dir, wind gust speed
     &                   dummy, dummy              !dir expected accuracy, spd exp accuracy
c
	  read(11,909)   alt(k),                   !altimeter
     &                   stnp(k),                  !station pressure
     &                   mslp(k),                  !MSL pressure
     &                   idummy,                   !3-h press change character
     &                   dummy,                    !3-h pressure change
     &                   dummy, dummy              !pressure exp accuracy, alt exp accuracy
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,915)  kkk_s,                     !num cld layers
     &                  dummy,                     !24-h max temperature
     &                  dummy                      !24-h min temperature
c
c.....	Read the cloud data if we have any.
c
	  if(kkk_s .gt. 0) then
	    do ii=1,kkk_s
  	      read(11,919) dum                     !layer cloud amount and height
	    enddo !ii
	  endif
c
	enddo !k
c
!       endfile(11)
	close(11)
c
c..... End of data reading.  Let's go home...
c
	jstatus = 1
	return
999     continue
	print *,' ++ ERROR opening ',ext(1:len_ext),
     1          ' file in READ_SFC_STATE ++'
        jstatus = -1
	return
        include 'lso_formats.inc'
	end
c
c
	subroutine read_sfc_temp(i4time,btime,n_obs_g,n_obs_b,
     &    stations,provider,lat,lon,elev,t,maxsta,jstatus)
c
c*****************************************************************************
c
cdoc    Routine to read the LAPS LSO surface data file and return the
cdoc    temperature.   The data is passed back to the calling routine in
cdoc    1-d arrays.
c
c	Changes:
c		P. Stamus  04-13-98  Original version (from write_surface_LS2).
c
c*****************************************************************************
c
	real lat(maxsta), lon(maxsta), elev(maxsta)
	real t(maxsta)

c
	integer i4time, jstatus
c
	character filetime*9, infile*256, btime*24
	character stations(maxsta)*20, provider(maxsta)*11
	character dum*132
c
	jstatus = 0
c
c.....	Get the file.
c
	call make_fnam_lp(i4time, filetime, istatus)
	call get_directory('lso', infile, len)
	infile = infile(1:len) // filetime(1:9) // '.lso'
c
	open(11,file=infile,status='old',err=999)
c
c.....	Read the header.
c
	read(11,900) btime,		! time
     &               n_obs_g,		! # of obs in the laps grid
     &               n_obs_b		! # of obs in the box
c
c.....	Read the station data.
c
	do k=1,n_obs_b
c
	   read(11,901)  stations(k),              !station id
     &                   idummy,                   !WMO id number
     &                   provider(k),              !data provider
     &                   lat(k), lon(k), elev(k),  !lat, lon, elev
     &                   idummy                    !obs time
c
	  read(11,919) dum
c
	  read(11,905)   t(k), dummy,              !temp, temp expected accuracy
     &                   dummy, dummy,             !dew point, dew point exp. accuracy
     &                   dummy, dummy              !Rel hum, rh expected accuracy
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,915)  kkk_s,                     !num cld layers
     &                  dummy,                     !24-h max temperature
     &                  dummy                      !24-h min temperature
c
c.....	Read the cloud data if we have any.
c
	  if(kkk_s .gt. 0) then
	    do ii=1,kkk_s
  	      read(11,919) dum                     !layer cloud amount and height
	    enddo !ii
	  endif
c
	enddo !k
c
!       endfile(11)
	close(11)
c
c..... End of data reading.  Let's go home...
c
	jstatus = 1
	return
999     continue
	print *,' ++ ERROR opening LSO file in READ_SFC_TEMP ++'
        jstatus = -1
	return
        include 'lso_formats.inc'
	end
c
c
        subroutine read_sfc_wind(i4time,ext,n_obs_g,n_obs_b,obstime
     1                          ,stations,provider,lat,lon,elev
     1                          ,dd,ff,dd_ea,ff_ea,maxsta,jstatus)
c
c*****************************************************************************
c
cdoc    Routine to read the LAPS LSO surface data file and return wind data.
cdoc    The data is passed back to the calling routine in 1-d arrays.
c
c	Changes:
c		P. Stamus  04-13-98  Original version (from write_surface_LS2).
c
c*****************************************************************************
c
	real lat(maxsta), lon(maxsta), elev(maxsta)
	real dd(maxsta)
	real ff(maxsta)
	real dd_ea(maxsta)
	real ff_ea(maxsta)
c
	integer i4time, jstatus
        integer obstime(maxsta)
c
	character filetime*9, infile*256, btime*24
	character stations(maxsta)*20, provider(maxsta)*11
	character dum*132
        character*(*) ext
c
	jstatus = 0
c
c.....	Get the file.
c
	call make_fnam_lp(i4time, filetime, istatus)
	call get_directory('lso', infile, len)
	infile = infile(1:len) // filetime(1:9) // '.' // ext
c
	open(11,file=infile,status='old',err=999)
c
c.....	Read the header.
c
	read(11,900) btime,		! time
     &               n_obs_g,		! # of obs in the laps grid
     &               n_obs_b		! # of obs in the box
c
c.....	Read the station data.
c
	do k=1,n_obs_b
c
	   read(11,901)  stations(k),              !station id
     &                   idummy,                   !WMO id number
     &                   provider(k),              !data provider
     &                   lat(k), lon(k), elev(k),  !lat, lon, elev
     &                   obstime(k)                !obs time
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,907)   dd(k), ff(k),             !wind dir, wind speed (kt)
     &                   dummy, dummy,             !wind gust dir, wind gust speed
     &                   dd_ea(k), ff_ea(k)        !dir expected accuracy, spd exp accuracy
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,915)  kkk_s,                     !num cld layers
     &                  dummy,                     !24-h max temperature
     &                  dummy                      !24-h min temperature
c
c.....	Read the cloud data if we have any.
c
	  if(kkk_s .gt. 0) then
	    do ii=1,kkk_s
  	      read(11,919) dum
	    enddo !ii
	  endif
c
	enddo !k
c
!       endfile(11)
	close(11)
c
c..... End of data reading.  Let's go home...
c
	jstatus = 1
	return
999     continue
	print *,' ++ ERROR opening LSO file in READ_SFC_WIND ++'
        jstatus = -1
	return
        include 'lso_formats.inc'
	end
c
c
	subroutine read_sfc_press(i4time,btime,n_obs_g,n_obs_b,
     &    stations,provider,lat,lon,elev,alt,stnp,mslp,maxsta,jstatus)
c
c*****************************************************************************
c
cdoc    Routine to read the LAPS LSO surface data file and return pressure
cdoc    data.   The data is passed back to the calling routine in 1-d arrays.
c
c	Changes:
c		P. Stamus  04-13-98  Original version (from write_surface_LS2).
c               J. Edwards 09-16-98  moved all format definitions to
c                                    src/include/lso_formats.inc
c                                    changed 909 definition to allow for
c                                    missing data
c
c*****************************************************************************
c
	real lat(maxsta), lon(maxsta), elev(maxsta)
	real alt(maxsta)
	real stnp(maxsta), mslp(maxsta)
c
	integer i4time, jstatus
c
	character filetime*9, infile*256, btime*24
	character stations(maxsta)*20, provider(maxsta)*11
	character dum*132
c
	jstatus = 0
c
c.....	Get the file.
c
	call make_fnam_lp(i4time, filetime, istatus)
	call get_directory('lso', infile, len)
	infile = infile(1:len) // filetime(1:9) // '.lso'
c
	open(11,file=infile,status='old',err=999)
c
c.....	Read the header.
c
	read(11,900) btime,		! time
     &               n_obs_g,		! # of obs in the laps grid
     &               n_obs_b		! # of obs in the box
c
c.....	Read the station data.
c
	do k=1,n_obs_b
c
	   read(11,901)  stations(k),              !station id
     &                   idummy,                   !WMO id number
     &                   provider(k),              !data provider
     &                   lat(k), lon(k), elev(k),  !lat, lon, elev
     &                   idummy                    !obs time
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,909)   alt(k),                   !altimeter
     &                   stnp(k),                  !station pressure
     &                   mslp(k),                  !MSL pressure
     &                   idummy,                   !3-h press change character
     &                   dummy,                    !3-h pressure change
     &                   dummy, dummy              !pressure exp accuracy, alt exp accuracy
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,915)  kkk_s,                     !num cld layers
     &                  dummy,                     !24-h max temperature
     &                  dummy                      !24-h min temperature
c
c.....	Read the cloud data if we have any.
c
	  if(kkk_s .gt. 0) then
	    do ii=1,kkk_s
  	      read(11,919) dum
	    enddo !ii
	  endif
c
	enddo !k
c
!       endfile(11)
	close(11)
c
c..... End of data reading.  Let's go home...
c
	jstatus = 1
	return
999     continue
	print *,' ++ ERROR opening LSO file in READ_SFC_PRESS ++'
        jstatus = -1
	return
        include 'lso_formats.inc'
	end
c
c
	subroutine read_sfc_precip(i4time,btime,n_obs_g,n_obs_b,
     &    stations,provider,lat,lon,elev,pcp1,pcp3,pcp6,pcp24,snow,
     &    maxsta,jstatus)
c
c*****************************************************************************
c
cdoc    Routine to read the LAPS LSO surface data file and return precipitation
cdoc    data.   The data is passed back to the calling routine in 1-d arrays.
c
c	Changes:
c		P. Stamus  04-13-98  Original version (from write_surface_LS2).
c               J. Edwards 09-16-98  moved all format definitions to
c                                    src/include/lso_formats.inc
c                                    changed 909 definition to allow for
c                                    missing data
c
c*****************************************************************************
c
	real lat(maxsta), lon(maxsta), elev(maxsta)
	real pcp1(maxsta), pcp3(maxsta), pcp6(maxsta), pcp24(maxsta)
	real snow(maxsta)
c
	integer i4time, jstatus
c
	character filetime*9, infile*256, btime*24
	character stations(maxsta)*20, provider(maxsta)*11
	character dum*132
c
	jstatus = 0
        n_obs_g = 0
        n_obs_b = 0
c
c.....	Get the file.
c
	call make_fnam_lp(i4time, filetime, istatus)
	call get_directory('lso', infile, len)
	infile = infile(1:len) // filetime(1:9) // '.lso'
c
	open(11,file=infile,status='old',err=999)
c
c.....	Read the header.
c
	read(11,900) btime,		! time
     &               n_obs_g,		! # of obs in the laps grid
     &               n_obs_b		! # of obs in the box
c
c.....	Read the station data.
c
	do k=1,n_obs_b
c
	   read(11,901)  stations(k),              !station id
     &                   idummy,                   !WMO id number
     &                   provider(k),              !data provider
     &                   lat(k), lon(k), elev(k),  !lat, lon, elev
     &                   idummy                    !obs time
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,919) dum
c
	  read(11,913)   pcp1(k),                  !1-h precipitation
     &                   pcp3(k),                  !3-h precipitation
     &                   pcp6(k),                  !6-h precipitation
     &                   pcp24(k),                 !24-h precipitation
     &                   snow(k),                  !snow depth
     &                   dummy, dummy              !precip and snow exp accuracy
c
	  read(11,915)  kkk_s,                     !num cld layers
     &                  dummy,                     !24-h max temperature
     &                  dummy                      !24-h min temperature
c
c.....	Read the cloud data if we have any.
c
	  if(kkk_s .gt. 0) then
	    do ii=1,kkk_s
  	      read(11,919) dum
	    enddo !ii
	  endif
c
	enddo !k
c
!       endfile(11)
	close(11)
c
c..... End of data reading.  Let's go home...
c
	jstatus = 1
	return
 999    continue
	print *,' ++ ERROR opening LSO file in READ_SFC_PRECIP ++'
        jstatus 	= -1
	return
        include 'lso_formats.inc'
	end
c
c
	subroutine read_surface_dataqc(i4time,btime,n_obs_g,n_obs_b,time,
     &    wmoid,stations,provider,wx,reptype,autostntype,lat,lon,elev,
     &    t,td,rh,dd,ff,ddg,ffg,alt,stnp,mslp,delpch,delp,vis,solar,
     &    sfct,sfcm,pcp1,pcp3,pcp6,pcp24,snow,kkk_s,max24t,min24t,t_ea,
     &    td_ea,rh_ea,dd_ea,ff_ea,alt_ea,p_ea,vis_ea,solar_ea,sfct_ea,
     &    sfcm_ea,pcp_ea,snow_ea,store_cldamt,store_cldht,
     &    maxsta,jstatus)
c
c*****************************************************************************
c
cdoc    Routine to read the LAPS LSO_QC surface data file with the expanded
cdoc    variable list. The data is passed back to the calling routine in 1-d
cdoc    arrays.
c
c	Changes:
c		P. Stamus  12-17-98  Original version (from read_surface_data).
c
c*****************************************************************************
c
	real lat(maxsta), lon(maxsta), elev(maxsta)
	real t(maxsta), t_ea(maxsta), max24t(maxsta), min24t(maxsta)
	real td(maxsta), td_ea(maxsta), rh(maxsta), rh_ea(maxsta)
	real dd(maxsta), ddg(maxsta), dd_ea(maxsta)
	real ff(maxsta), ffg(maxsta), ff_ea(maxsta)
	real alt(maxsta), alt_ea(maxsta), delp(maxsta)
	real stnp(maxsta), mslp(maxsta), p_ea(maxsta)
	real vis(maxsta), vis_ea(maxsta)
	real solar(maxsta), solar_ea(maxsta)
	real sfct(maxsta), sfct_ea(maxsta)
	real sfcm(maxsta), sfcm_ea(maxsta)
	real pcp1(maxsta), pcp3(maxsta), pcp6(maxsta), pcp24(maxsta)
	real snow(maxsta), snow_ea(maxsta), pcp_ea(maxsta)
	real store_cldht(maxsta,5)
c
	integer i4time, wmoid(maxsta), jstatus
	integer time(maxsta), delpch(maxsta), kkk_s(maxsta)
c
	character filetime*9, infile*256, btime*24
	character stations(maxsta)*20, provider(maxsta)*11
        character reptype(maxsta)*6, autostntype(maxsta)*6
	character wx(maxsta)*25, store_cldamt(maxsta,5)*4
c
c
c.....  Blank out the character arrays.
c
	jstatus = 0
	do i=1,maxsta
	   stations(i) = '                    '
	   provider(i) = '           '
	   reptype(i)  = '      '
	   autostntype(i)  = '      '
	   wx(i) = '                         '
	   do j=1,5
	      store_cldamt(i,j) = '    '
	   enddo !j
	enddo !i
c
c.....	Get the file.
c
	call make_fnam_lp(i4time, filetime, istatus)
	call get_directory('lso', infile, len)
	infile = infile(1:len) // filetime(1:9) // '.lso_qc'
c
	open(11,file=infile,status='old',err=999)
c
c.....	Read the header.
c
	read(11,900) btime,		! time
     &               n_obs_g,		! # of obs in the laps grid
     &               n_obs_b		! # of obs in the box
c
c.....	Read the station data.
c
	do k=1,n_obs_b
c
	   read(11,901)  stations(k),              !station id
     &                   wmoid(k),                 !WMO id number
     &                   provider(k),              !data provider
     &                   lat(k), lon(k), elev(k),  !lat, lon, elev
     &                   time(k)		   !obs time
c
	  read(11,903)   reptype(k),               !station report type
     &                   autostntype(k),           !station type (manual/auto)
     &                   wx(k)                     !present weather
c
	  read(11,905)   t(k), t_ea(k),            !temp, temp expected accuracy
     &                   td(k), td_ea(k),          !dew point, dew point exp. accuracy
     &                   rh(k), rh_ea(k)           !Rel hum, rh expected accuracy
c
	  read(11,907)   dd(k), ff(k),             !wind dir, wind speed
     &                   ddg(k), ffg(k),           !wind gust dir, wind gust speed
     &                   dd_ea(k), ff_ea(k)        !dir expected accuracy, spd exp accuracy
c
	  read(11,909)   alt(k),                   !altimeter
     &                   stnp(k),                  !station pressure
     &                   mslp(k),                  !MSL pressure
     &                   delpch(k),                !3-h press change character
     &                   delp(k),                  !3-h pressure change
     &                   p_ea(k), alt_ea(k)        !pressure exp accuracy, alt exp accuracy
c
	  read(11,911)   vis(k), vis_ea(k),        !visibility, vis exp accuracy
     &                   solar(k), solar_ea(k),    !solar, solar exp accuracy
     &                   sfct(k), sfct_ea(k),      !soil/water temp, soil/water temp exp accuracy
     &                   sfcm(k), sfcm_ea(k)       !soil moist, soil moist temp exp accuracy
c
	  read(11,913)   pcp1(k),                  !1-h precipitation
     &                   pcp3(k),                  !3-h precipitation
     &                   pcp6(k),                  !6-h precipitation
     &                   pcp24(k),                 !24-h precipitation
     &                   snow(k),                  !snow depth
     &                   pcp_ea(k), snow_ea(k)     !precip and snow exp accuracy
c
	  read(11,915)  kkk_s(k),                  !num cld layers
     &                  max24t(k),                 !24-h max temperature
     &                  min24t(k)                  !24-h min temperature
c
c.....	Read the cloud data if we have any.
c
	  if(kkk_s(k) .gt. 0) then
	    do ii=1,kkk_s(k)
  	      read(11,917) store_cldamt(k,ii), store_cldht(k,ii) ! layer cloud amount and height
                                                                 ! (meters MSL)
	    enddo !ii
	  endif
c
	enddo !k
c
!       endfile(11)
	close(11)
c
c..... End of data reading.  Let's go home...
c
	jstatus = 1
	return
999     continue
        print *,'Could not open file: ',infile
	print *,' ++ ERROR opening LSO_QC file in READ_SURFACE_DATAQC ++'
        jstatus = -1
	return
        include 'lso_formats.inc'
	end

max_stars_repo_name: maxinye/laps-mirror

max_stars_repo_licenses: Intel Unlicense OLDAP-2.2.1 NetCDF

max_issues_repo_name: longwosion/laps-mirror

max_issues_repo_licenses: Intel NetCDF OLDAP-2.2.1 Unlicense

max_forks_repo_name: longwosion/laps-mirror

max_forks_repo_licenses: Intel NetCDF OLDAP-2.2.1 Unlicense


text:
module to_lower_mod

  implicit none

  public :: to_lower_index

contains

  pure function to_lower_index(c) result(res)
    character(len=1), intent(in) :: c
    character(len=1) :: res

    integer :: i
    character(len=26), parameter :: &
      upper = "ABCDEFGHIJKLMNOPQRSTUVWXYZ", &
      lower = "abcdefghijklmnopqrstuvwxyz"

    res = c
    i = index(upper, c)
    if (i > 0) res = lower(i:i)
  end function

end module

max_stars_repo_name: ivan-pi/fortran-ascii

max_issues_repo_name: ivan-pi/fortran-ascii

max_forks_repo_name: ivan-pi/fortran-ascii


MODULE ldftra
  USE oce
  USE dom_oce
  USE phycst
  USE ldfslp
  USE ldfc1d_c2d
  USE diaptr
  USE in_out_manager
  USE iom
  USE lib_mpp
  USE lbclnk
  IMPLICIT NONE
  PRIVATE
  PUBLIC :: ldf_tra_init
  PUBLIC :: ldf_tra
  PUBLIC :: ldf_eiv_init
  PUBLIC :: ldf_eiv
  PUBLIC :: ldf_eiv_trp
  PUBLIC :: ldf_eiv_dia
  LOGICAL, PUBLIC :: ln_traldf_OFF
  LOGICAL, PUBLIC :: ln_traldf_lap
  LOGICAL, PUBLIC :: ln_traldf_blp
  LOGICAL, PUBLIC :: ln_traldf_lev
  LOGICAL, PUBLIC :: ln_traldf_hor
  LOGICAL, PUBLIC :: ln_traldf_msc
  INTEGER, PUBLIC :: nn_aht_ijk_t
  REAL(KIND = wp), PUBLIC :: rn_Ud
  REAL(KIND = wp), PUBLIC :: rn_Ld
  LOGICAL, PUBLIC :: ln_ldfeiv
  LOGICAL, PUBLIC :: ln_ldfeiv_dia
  INTEGER, PUBLIC :: nn_aei_ijk_t
  REAL(KIND = wp), PUBLIC :: rn_Ue
  REAL(KIND = wp), PUBLIC :: rn_Le
  INTEGER, PARAMETER, PUBLIC :: np_ERROR = - 10
  INTEGER, PARAMETER, PUBLIC :: np_no_ldf = 00
  INTEGER, PARAMETER, PUBLIC :: np_lap = 10, np_blp = 20
  INTEGER, PARAMETER, PUBLIC :: np_lap_i = 11, np_blp_i = 21
  INTEGER, PARAMETER, PUBLIC :: np_lap_it = 12, np_blp_it = 22
  INTEGER, PUBLIC :: nldf_tra = 0
  LOGICAL, PUBLIC :: l_ldftra_time = .FALSE.
  LOGICAL, PUBLIC :: l_ldfeiv_time = .FALSE.
  REAL(KIND = wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:, :, :) :: ahtu, ahtv
  REAL(KIND = wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:, :, :) :: aeiu, aeiv
  REAL(KIND = wp) :: aht0, aei0
  REAL(KIND = wp) :: r1_2 = 0.5_wp
  REAL(KIND = wp) :: r1_4 = 0.25_wp
  REAL(KIND = wp) :: r1_12 = 1._wp / 12._wp
  CONTAINS
  SUBROUTINE ldf_tra_init
    INTEGER :: jk
    INTEGER :: ioptio, ierr, inum, ios, inn
    REAL(KIND = wp) :: zah_max, zUfac
    CHARACTER(LEN = 5) :: cl_Units
    NAMELIST /namtra_ldf/ ln_traldf_OFF, ln_traldf_lap, ln_traldf_blp, ln_traldf_lev, ln_traldf_hor, ln_traldf_triad, ln_traldf_iso, ln_traldf_msc, rn_slpmax, ln_triad_iso, ln_botmix_triad, rn_sw_triad, nn_aht_ijk_t, rn_Ud, rn_Ld
    IF (lwp) THEN
      WRITE(numout, FMT = *)
      WRITE(numout, FMT = *) 'ldf_tra_init : lateral tracer diffusion'
      WRITE(numout, FMT = *) '~~~~~~~~~~~~ '
    END IF
    REWIND(UNIT = numnam_ref)
    READ(numnam_ref, namtra_ldf, IOSTAT = ios, ERR = 901)
901 IF (ios /= 0) CALL ctl_nam(ios, 'namtra_ldf in reference namelist', lwp)
    REWIND(UNIT = numnam_cfg)
    READ(numnam_cfg, namtra_ldf, IOSTAT = ios, ERR = 902)
902 IF (ios > 0) CALL ctl_nam(ios, 'namtra_ldf in configuration namelist', lwp)
    IF (lwm) WRITE(numond, namtra_ldf)
    IF (lwp) THEN
      WRITE(numout, FMT = *) '   Namelist : namtra_ldf --- lateral mixing parameters (type, direction, coefficients)'
      WRITE(numout, FMT = *) '      type :'
      WRITE(numout, FMT = *) '         no explicit diffusion                   ln_traldf_OFF   = ', ln_traldf_OFF
      WRITE(numout, FMT = *) '         laplacian operator                      ln_traldf_lap   = ', ln_traldf_lap
      WRITE(numout, FMT = *) '         bilaplacian operator                    ln_traldf_blp   = ', ln_traldf_blp
      WRITE(numout, FMT = *) '      direction of action :'
      WRITE(numout, FMT = *) '         iso-level                               ln_traldf_lev   = ', ln_traldf_lev
      WRITE(numout, FMT = *) '         horizontal (geopotential)               ln_traldf_hor   = ', ln_traldf_hor
      WRITE(numout, FMT = *) '         iso-neutral Madec operator              ln_traldf_iso   = ', ln_traldf_iso
      WRITE(numout, FMT = *) '         iso-neutral triad operator              ln_traldf_triad = ', ln_traldf_triad
      WRITE(numout, FMT = *) '            use the Method of Stab. Correction   ln_traldf_msc   = ', ln_traldf_msc
      WRITE(numout, FMT = *) '            maximum isoppycnal slope             rn_slpmax       = ', rn_slpmax
      WRITE(numout, FMT = *) '            pure lateral mixing in ML            ln_triad_iso    = ', ln_triad_iso
      WRITE(numout, FMT = *) '            switching triad or not               rn_sw_triad     = ', rn_sw_triad
      WRITE(numout, FMT = *) '            lateral mixing on bottom             ln_botmix_triad = ', ln_botmix_triad
      WRITE(numout, FMT = *) '      coefficients :'
      WRITE(numout, FMT = *) '         type of time-space variation            nn_aht_ijk_t    = ', nn_aht_ijk_t
      WRITE(numout, FMT = *) '            lateral diffusive velocity (if cst)  rn_Ud           = ', rn_Ud, ' m/s'
      WRITE(numout, FMT = *) '            lateral diffusive length   (if cst)  rn_Ld           = ', rn_Ld, ' m'
    END IF
    nldf_tra = np_ERROR
    ioptio = 0
    IF (ln_traldf_off) THEN
      nldf_tra = np_no_ldf
      ioptio = ioptio + 1
    END IF
    IF (ln_traldf_lap) THEN
      ioptio = ioptio + 1
    END IF
    IF (ln_traldf_blp) THEN
      ioptio = ioptio + 1
    END IF
    IF (ioptio /= 1) CALL ctl_stop('tra_ldf_init: use ONE of the 3 operator options (NONE/lap/blp)')
    IF (.NOT. ln_traldf_OFF) THEN
      ioptio = 0
      IF (ln_traldf_lev) ioptio = ioptio + 1
      IF (ln_traldf_hor) ioptio = ioptio + 1
      IF (ln_traldf_iso) ioptio = ioptio + 1
      IF (ln_traldf_triad) ioptio = ioptio + 1
      IF (ioptio /= 1) CALL ctl_stop('tra_ldf_init: use ONE direction (level/hor/iso/triad)')
      ierr = 0
      IF (ln_traldf_lap) THEN
        IF (ln_zco) THEN
          IF (ln_traldf_lev) nldf_tra = np_lap
          IF (ln_traldf_hor) nldf_tra = np_lap
          IF (ln_traldf_iso) nldf_tra = np_lap_i
          IF (ln_traldf_triad) nldf_tra = np_lap_it
        END IF
        IF (ln_zps) THEN
          IF (ln_traldf_lev) ierr = 1
          IF (ln_traldf_hor) nldf_tra = np_lap
          IF (ln_traldf_iso) nldf_tra = np_lap_i
          IF (ln_traldf_triad) nldf_tra = np_lap_it
        END IF
        IF (ln_sco) THEN
          IF (ln_traldf_lev) nldf_tra = np_lap
          IF (ln_traldf_hor) nldf_tra = np_lap_i
          IF (ln_traldf_iso) nldf_tra = np_lap_i
          IF (ln_traldf_triad) nldf_tra = np_lap_it
        END IF
      END IF
      IF (ln_traldf_blp) THEN
        IF (ln_zco) THEN
          IF (ln_traldf_lev) nldf_tra = np_blp
          IF (ln_traldf_hor) nldf_tra = np_blp
          IF (ln_traldf_iso) nldf_tra = np_blp_i
          IF (ln_traldf_triad) nldf_tra = np_blp_it
        END IF
        IF (ln_zps) THEN
          IF (ln_traldf_lev) ierr = 1
          IF (ln_traldf_hor) nldf_tra = np_blp
          IF (ln_traldf_iso) nldf_tra = np_blp_i
          IF (ln_traldf_triad) nldf_tra = np_blp_it
        END IF
        IF (ln_sco) THEN
          IF (ln_traldf_lev) nldf_tra = np_blp
          IF (ln_traldf_hor) nldf_tra = np_blp_it
          IF (ln_traldf_iso) nldf_tra = np_blp_i
          IF (ln_traldf_triad) nldf_tra = np_blp_it
        END IF
      END IF
      IF (ierr == 1) CALL ctl_stop('iso-level in z-partial step, not allowed')
    END IF
    IF (ln_ldfeiv .AND. .NOT. (ln_traldf_iso .OR. ln_traldf_triad)) CALL ctl_stop('ln_ldfeiv=T requires iso-neutral laplacian diffusion')
    IF (ln_isfcav .AND. ln_traldf_triad) CALL ctl_stop(' ice shelf cavity and traldf_triad not tested')
    IF (nldf_tra == np_lap_i .OR. nldf_tra == np_lap_it .OR. nldf_tra == np_blp_i .OR. nldf_tra == np_blp_it) l_ldfslp = .TRUE.
    IF (ln_traldf_blp .AND. (ln_traldf_iso .OR. ln_traldf_triad)) THEN
      IF (.NOT. ln_traldf_msc) CALL ctl_stop('tra_ldf_init: iso-neutral bilaplacian requires ln_traldf_msc=.true.')
    END IF
    IF (lwp) THEN
      WRITE(numout, FMT = *)
      SELECT CASE (nldf_tra)
      CASE (np_no_ldf)
        WRITE(numout, FMT = *) '   ==>>>   NO lateral diffusion'
      CASE (np_lap)
        WRITE(numout, FMT = *) '   ==>>>   laplacian iso-level operator'
      CASE (np_lap_i)
        WRITE(numout, FMT = *) '   ==>>>   Rotated laplacian operator (standard)'
      CASE (np_lap_it)
        WRITE(numout, FMT = *) '   ==>>>   Rotated laplacian operator (triad)'
      CASE (np_blp)
        WRITE(numout, FMT = *) '   ==>>>   bilaplacian iso-level operator'
      CASE (np_blp_i)
        WRITE(numout, FMT = *) '   ==>>>   Rotated bilaplacian operator (standard)'
      CASE (np_blp_it)
        WRITE(numout, FMT = *) '   ==>>>   Rotated bilaplacian operator (triad)'
      END SELECT
      WRITE(numout, FMT = *)
    END IF
    l_ldftra_time = .FALSE.
    IF (ln_traldf_OFF) THEN
      IF (lwp) WRITE(numout, FMT = *) '   ==>>>   No diffusive operator selected. ahtu and ahtv are not allocated'
      RETURN
    ELSE
      ALLOCATE(ahtu(jpi, jpj, jpk), ahtv(jpi, jpj, jpk), STAT = ierr)
      IF (ierr /= 0) CALL ctl_stop('STOP', 'ldf_tra_init: failed to allocate arrays')
      !$ACC KERNELS
      ahtu(:, :, jpk) = 0._wp
      ahtv(:, :, jpk) = 0._wp
      !$ACC END KERNELS
      IF (ln_traldf_lap) THEN
        zufac = r1_2 * rn_ud
        inn = 1
        cl_units = ' m2/s'
      ELSE IF (ln_traldf_blp) THEN
        zufac = r1_12 * rn_ud
        inn = 3
        cl_units = ' m4/s'
      END IF
      aht0 = zUfac * rn_Ld ** inn
      zah_max = zUfac * (ra * rad) ** inn
      SELECT CASE (nn_aht_ijk_t)
      CASE (0)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy diffusivity = constant = ', aht0, cl_Units
        !$ACC KERNELS
        ahtu(:, :, 1 : jpkm1) = aht0
        ahtv(:, :, 1 : jpkm1) = aht0
        !$ACC END KERNELS
      CASE (10)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy diffusivity = F( depth )'
        IF (lwp) WRITE(numout, FMT = *) '           surface eddy diffusivity = constant = ', aht0, cl_Units
        !$ACC KERNELS
        ahtu(:, :, 1) = aht0
        ahtv(:, :, 1) = aht0
        !$ACC END KERNELS
        CALL ldf_c1d('TRA', ahtu(:, :, 1), ahtv(:, :, 1), ahtu, ahtv)
      CASE (- 20)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy diffusivity = F(i,j) read in eddy_diffusivity.nc file'
        CALL iom_open('eddy_diffusivity_2D.nc', inum)
        CALL iom_get(inum, jpdom_data, 'ahtu_2D', ahtu(:, :, 1))
        CALL iom_get(inum, jpdom_data, 'ahtv_2D', ahtv(:, :, 1))
        CALL iom_close(inum)
        !$ACC KERNELS
        DO jk = 2, jpkm1
          ahtu(:, :, jk) = ahtu(:, :, 1)
          ahtv(:, :, jk) = ahtv(:, :, 1)
        END DO
        !$ACC END KERNELS
      CASE (20)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy diffusivity = F( e1, e2 ) or F( e1^3, e2^3 ) (lap or blp case)'
        IF (lwp) WRITE(numout, FMT = *) '           using a fixed diffusive velocity = ', rn_Ud, ' m/s   and   Ld = Max(e1,e2)'
        CALL ldf_c2d('TRA', zUfac, inn, ahtu, ahtv)
      CASE (21)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy diffusivity = F( latitude, longitude, time )'
        IF (lwp) WRITE(numout, FMT = *) '                            = F( growth rate of baroclinic instability )'
        IF (lwp) WRITE(numout, FMT = *) '            min value = 0.2 * aht0 (with aht0= 1/2 rn_Ud*rn_Ld)'
        IF (lwp) WRITE(numout, FMT = *) '            max value =       aei0 (with aei0=1/2 rn_Ue*Le  increased to aht0 within 20N-20S'
        l_ldftra_time = .TRUE.
        IF (ln_traldf_blp) CALL ctl_stop('ldf_tra_init: aht=F( growth rate of baroc. insta .)', '              incompatible with bilaplacian operator')
      CASE (- 30)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy diffusivity = F(i,j,k) read in eddy_diffusivity.nc file'
        CALL iom_open('eddy_diffusivity_3D.nc', inum)
        CALL iom_get(inum, jpdom_data, 'ahtu_3D', ahtu)
        CALL iom_get(inum, jpdom_data, 'ahtv_3D', ahtv)
        CALL iom_close(inum)
      CASE (30)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy diffusivity = F( latitude, longitude, depth )'
        IF (lwp) WRITE(numout, FMT = *) '           using a fixed diffusive velocity = ', rn_Ud, ' m/s   and   Ld = Max(e1,e2)'
        CALL ldf_c2d('TRA', zUfac, inn, ahtu, ahtv)
        CALL ldf_c1d('TRA', ahtu(:, :, 1), ahtv(:, :, 1), ahtu, ahtv)
      CASE (31)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy diffusivity = F( latitude, longitude, depth , time )'
        IF (lwp) WRITE(numout, FMT = *) '           proportional to the velocity : 1/2 |u|e or 1/12 |u|e^3'
        l_ldftra_time = .TRUE.
      CASE DEFAULT
        CALL ctl_stop('ldf_tra_init: wrong choice for nn_aht_ijk_t, the type of space-time variation of aht')
      END SELECT
      IF (.NOT. l_ldftra_time) THEN
        IF (ln_traldf_lap) THEN
          !$ACC KERNELS
          ahtu(:, :, 1 : jpkm1) = ahtu(:, :, 1 : jpkm1) * umask(:, :, 1 : jpkm1)
          ahtv(:, :, 1 : jpkm1) = ahtv(:, :, 1 : jpkm1) * vmask(:, :, 1 : jpkm1)
          !$ACC END KERNELS
        ELSE IF (ln_traldf_blp) THEN
          !$ACC KERNELS
          ahtu(:, :, 1 : jpkm1) = SQRT(ahtu(:, :, 1 : jpkm1)) * umask(:, :, 1 : jpkm1)
          ahtv(:, :, 1 : jpkm1) = SQRT(ahtv(:, :, 1 : jpkm1)) * vmask(:, :, 1 : jpkm1)
          !$ACC END KERNELS
        END IF
      END IF
    END IF
  END SUBROUTINE ldf_tra_init
  SUBROUTINE ldf_tra(kt)
    INTEGER, INTENT(IN) :: kt
    INTEGER :: ji, jj, jk
    REAL(KIND = wp) :: zaht, zahf, zaht_min, zDaht, z1_f20
    IF (ln_ldfeiv .AND. nn_aei_ijk_t == 21) THEN
      CALL ldf_eiv(kt, aei0, aeiu, aeiv)
    END IF
    SELECT CASE (nn_aht_ijk_t)
    CASE (21)
      IF (ln_ldfeiv .AND. nn_aei_ijk_t == 21) THEN
        !$ACC KERNELS
        ahtu(:, :, 1) = aeiu(:, :, 1)
        ahtv(:, :, 1) = aeiv(:, :, 1)
        !$ACC END KERNELS
      ELSE
        CALL ldf_eiv(kt, aht0, ahtu, ahtv)
      END IF
      !$ACC KERNELS
      z1_f20 = 1._wp / (2._wp * omega * SIN(rad * 20._wp))
      zaht_min = 0.2_wp * aht0
      zDaht = aht0 - zaht_min
      DO jj = 1, jpj
        DO ji = 1, jpi
          zaht = (1._wp - MIN(1._wp, ABS(ff_t(ji, jj) * z1_f20))) * zDaht
          zahf = (1._wp - MIN(1._wp, ABS(ff_f(ji, jj) * z1_f20))) * zDaht
          ahtu(ji, jj, 1) = (MAX(zaht_min, ahtu(ji, jj, 1)) + zaht)
          ahtv(ji, jj, 1) = (MAX(zaht_min, ahtv(ji, jj, 1)) + zahf)
        END DO
      END DO
      DO jk = 1, jpkm1
        ahtu(:, :, jk) = ahtu(:, :, 1) * umask(:, :, jk)
        ahtv(:, :, jk) = ahtv(:, :, 1) * vmask(:, :, jk)
      END DO
      !$ACC END KERNELS
    CASE (31)
      IF (ln_traldf_lap) THEN
        !$ACC KERNELS
        DO jk = 1, jpkm1
          ahtu(:, :, jk) = ABS(ub(:, :, jk)) * e1u(:, :) * r1_12
          ahtv(:, :, jk) = ABS(vb(:, :, jk)) * e2v(:, :) * r1_12
        END DO
        !$ACC END KERNELS
      ELSE IF (ln_traldf_blp) THEN
        !$ACC KERNELS
        DO jk = 1, jpkm1
          ahtu(:, :, jk) = SQRT(ABS(ub(:, :, jk)) * e1u(:, :) * r1_12) * e1u(:, :)
          ahtv(:, :, jk) = SQRT(ABS(vb(:, :, jk)) * e2v(:, :) * r1_12) * e2v(:, :)
        END DO
        !$ACC END KERNELS
      END IF
    END SELECT
    CALL iom_put("ahtu_2d", ahtu(:, :, 1))
    CALL iom_put("ahtv_2d", ahtv(:, :, 1))
    CALL iom_put("ahtu_3d", ahtu(:, :, :))
    CALL iom_put("ahtv_3d", ahtv(:, :, :))
    IF (ln_ldfeiv) THEN
      CALL iom_put("aeiu_2d", aeiu(:, :, 1))
      CALL iom_put("aeiv_2d", aeiv(:, :, 1))
      CALL iom_put("aeiu_3d", aeiu(:, :, :))
      CALL iom_put("aeiv_3d", aeiv(:, :, :))
    END IF
  END SUBROUTINE ldf_tra
  SUBROUTINE ldf_eiv_init
    INTEGER :: jk
    INTEGER :: ierr, inum, ios, inn
    REAL(KIND = wp) :: zah_max, zUfac
    NAMELIST /namtra_eiv/ ln_ldfeiv, ln_ldfeiv_dia, nn_aei_ijk_t, rn_Ue, rn_Le
    IF (lwp) THEN
      WRITE(numout, FMT = *)
      WRITE(numout, FMT = *) 'ldf_eiv_init : eddy induced velocity parametrization'
      WRITE(numout, FMT = *) '~~~~~~~~~~~~ '
    END IF
    REWIND(UNIT = numnam_ref)
    READ(numnam_ref, namtra_eiv, IOSTAT = ios, ERR = 901)
901 IF (ios /= 0) CALL ctl_nam(ios, 'namtra_eiv in reference namelist', lwp)
    REWIND(UNIT = numnam_cfg)
    READ(numnam_cfg, namtra_eiv, IOSTAT = ios, ERR = 902)
902 IF (ios > 0) CALL ctl_nam(ios, 'namtra_eiv in configuration namelist', lwp)
    IF (lwm) WRITE(numond, namtra_eiv)
    IF (lwp) THEN
      WRITE(numout, FMT = *) '   Namelist namtra_eiv : '
      WRITE(numout, FMT = *) '      Eddy Induced Velocity (eiv) param.         ln_ldfeiv     = ', ln_ldfeiv
      WRITE(numout, FMT = *) '      eiv streamfunction & velocity diag.        ln_ldfeiv_dia = ', ln_ldfeiv_dia
      WRITE(numout, FMT = *) '      coefficients :'
      WRITE(numout, FMT = *) '         type of time-space variation            nn_aei_ijk_t  = ', nn_aht_ijk_t
      WRITE(numout, FMT = *) '         lateral diffusive velocity (if cst)     rn_Ue         = ', rn_Ue, ' m/s'
      WRITE(numout, FMT = *) '         lateral diffusive length   (if cst)     rn_Le         = ', rn_Le, ' m'
      WRITE(numout, FMT = *)
    END IF
    l_ldfeiv_time = .FALSE.
    IF (.NOT. ln_ldfeiv) THEN
      IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy induced velocity param is NOT used'
      ln_ldfeiv_dia = .FALSE.
    ELSE
      IF (lwp) WRITE(numout, FMT = *) '   ==>>>   use eddy induced velocity parametrization'
      IF (lwp) WRITE(numout, FMT = *)
      IF (ln_traldf_blp) CALL ctl_stop('ldf_eiv_init: eddy induced velocity ONLY with laplacian diffusivity')
      ALLOCATE(aeiu(jpi, jpj, jpk), aeiv(jpi, jpj, jpk), STAT = ierr)
      IF (ierr /= 0) CALL ctl_stop('STOP', 'ldf_eiv: failed to allocate arrays')
      !$ACC KERNELS
      aeiu(:, :, jpk) = 0._wp
      aeiv(:, :, jpk) = 0._wp
      zUfac = r1_2 * rn_Ue
      inn = 1
      aei0 = zUfac * rn_Le ** inn
      zah_max = zUfac * (ra * rad) ** inn
      !$ACC END KERNELS
      SELECT CASE (nn_aei_ijk_t)
      CASE (0)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy induced velocity coef. = constant = ', aei0, ' m2/s'
        !$ACC KERNELS
        aeiu(:, :, 1 : jpkm1) = aei0
        aeiv(:, :, 1 : jpkm1) = aei0
        !$ACC END KERNELS
      CASE (10)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy induced velocity coef. = F( depth )'
        IF (lwp) WRITE(numout, FMT = *) '           surface eddy diffusivity = constant = ', aht0, ' m2/s'
        !$ACC KERNELS
        aeiu(:, :, 1) = aei0
        aeiv(:, :, 1) = aei0
        !$ACC END KERNELS
        CALL ldf_c1d('TRA', aeiu(:, :, 1), aeiv(:, :, 1), aeiu, aeiv)
      CASE (- 20)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy induced velocity coef. = F(i,j) read in eddy_diffusivity_2D.nc file'
        CALL iom_open('eddy_induced_velocity_2D.nc', inum)
        CALL iom_get(inum, jpdom_data, 'aeiu', aeiu(:, :, 1))
        CALL iom_get(inum, jpdom_data, 'aeiv', aeiv(:, :, 1))
        CALL iom_close(inum)
        !$ACC KERNELS
        DO jk = 2, jpkm1
          aeiu(:, :, jk) = aeiu(:, :, 1)
          aeiv(:, :, jk) = aeiv(:, :, 1)
        END DO
        !$ACC END KERNELS
      CASE (20)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy induced velocity coef. = F( e1, e2 )'
        IF (lwp) WRITE(numout, FMT = *) '           using a fixed diffusive velocity = ', rn_Ue, ' m/s   and   Le = Max(e1,e2)'
        CALL ldf_c2d('TRA', zUfac, inn, aeiu, aeiv)
      CASE (21)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy induced velocity coef. = F( latitude, longitude, time )'
        IF (lwp) WRITE(numout, FMT = *) '                                       = F( growth rate of baroclinic instability )'
        IF (lwp) WRITE(numout, FMT = *) '           maximum allowed value: aei0 = ', aei0, ' m2/s'
        l_ldfeiv_time = .TRUE.
      CASE (- 30)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy induced velocity coef. = F(i,j,k) read in eddy_diffusivity_3D.nc file'
        CALL iom_open('eddy_induced_velocity_3D.nc', inum)
        CALL iom_get(inum, jpdom_data, 'aeiu', aeiu)
        CALL iom_get(inum, jpdom_data, 'aeiv', aeiv)
        CALL iom_close(inum)
      CASE (30)
        IF (lwp) WRITE(numout, FMT = *) '   ==>>>   eddy induced velocity coef. = F( latitude, longitude, depth )'
        CALL ldf_c2d('TRA', zUfac, inn, aeiu, aeiv)
        CALL ldf_c1d('TRA', aeiu(:, :, 1), aeiv(:, :, 1), aeiu, aeiv)
      CASE DEFAULT
        CALL ctl_stop('ldf_tra_init: wrong choice for nn_aei_ijk_t, the type of space-time variation of aei')
      END SELECT
      IF (.NOT. l_ldfeiv_time) THEN
        !$ACC KERNELS
        DO jk = 1, jpkm1
          aeiu(:, :, jk) = aeiu(:, :, jk) * umask(:, :, jk)
          ahtv(:, :, jk) = ahtv(:, :, jk) * vmask(:, :, jk)
        END DO
        !$ACC END KERNELS
      END IF
    END IF
  END SUBROUTINE ldf_eiv_init
  SUBROUTINE ldf_eiv(kt, paei0, paeiu, paeiv)
    INTEGER, INTENT(IN   ) :: kt
    REAL(KIND = wp), INTENT(INOUT) :: paei0
    REAL(KIND = wp), DIMENSION(jpi, jpj, jpk), INTENT(INOUT) :: paeiu, paeiv
    INTEGER :: ji, jj, jk
    REAL(KIND = wp) :: zfw, ze3w, zn2, z1_f20, zaht, zaht_min, zzaei
    REAL(KIND = wp), DIMENSION(jpi, jpj) :: zn, zah, zhw, zRo, zaeiw
    !$ACC KERNELS
    zn(:, :) = 0._wp
    zhw(:, :) = 5._wp
    zah(:, :) = 0._wp
    zRo(:, :) = 0._wp
    !$ACC END KERNELS
    IF (ln_traldf_triad) THEN
      !$ACC KERNELS
      DO jk = 1, jpk
        DO jj = 2, jpjm1
          DO ji = 2, jpim1
            zn2 = MAX(rn2b(ji, jj, jk), 0._wp)
            zn(ji, jj) = zn(ji, jj) + SQRT(zn2) * e3w_n(ji, jj, jk)
            ze3w = e3w_n(ji, jj, jk) * tmask(ji, jj, jk)
            zah(ji, jj) = zah(ji, jj) + zn2 * wslp2(ji, jj, jk) * ze3w
            zhw(ji, jj) = zhw(ji, jj) + ze3w
          END DO
        END DO
      END DO
      !$ACC END KERNELS
    ELSE
      !$ACC KERNELS
      DO jk = 1, jpk
        DO jj = 2, jpjm1
          DO ji = 2, jpim1
            zn2 = MAX(rn2b(ji, jj, jk), 0._wp)
            zn(ji, jj) = zn(ji, jj) + SQRT(zn2) * e3w_n(ji, jj, jk)
            ze3w = e3w_n(ji, jj, jk) * tmask(ji, jj, jk)
            zah(ji, jj) = zah(ji, jj) + zn2 * (wslpi(ji, jj, jk) * wslpi(ji, jj, jk) + wslpj(ji, jj, jk) * wslpj(ji, jj, jk)) * ze3w
            zhw(ji, jj) = zhw(ji, jj) + ze3w
          END DO
        END DO
      END DO
      !$ACC END KERNELS
    END IF
    !$ACC KERNELS
    DO jj = 2, jpjm1
      DO ji = 2, jpim1
        zfw = MAX(ABS(2. * omega * SIN(rad * gphit(ji, jj))), 1.E-10)
        zRo(ji, jj) = MAX(2.E3, MIN(.4 * zn(ji, jj) / zfw, 40.E3))
        zaeiw(ji, jj) = zRo(ji, jj) * zRo(ji, jj) * SQRT(zah(ji, jj) / zhw(ji, jj)) * tmask(ji, jj, 1)
      END DO
    END DO
    z1_f20 = 1._wp / (2._wp * omega * SIN(rad * 20._wp))
    DO jj = 2, jpjm1
      DO ji = 2, jpim1
        zzaei = MIN(1._wp, ABS(ff_t(ji, jj) * z1_f20)) * zaeiw(ji, jj)
        zaeiw(ji, jj) = MIN(zzaei, paei0)
      END DO
    END DO
    !$ACC END KERNELS
    CALL lbc_lnk(zaeiw(:, :), 'W', 1.)
    !$ACC KERNELS
    DO jj = 2, jpjm1
      DO ji = 2, jpim1
        paeiu(ji, jj, 1) = 0.5_wp * (zaeiw(ji, jj) + zaeiw(ji + 1, jj)) * umask(ji, jj, 1)
        paeiv(ji, jj, 1) = 0.5_wp * (zaeiw(ji, jj) + zaeiw(ji, jj + 1)) * vmask(ji, jj, 1)
      END DO
    END DO
    !$ACC END KERNELS
    CALL lbc_lnk_multi(paeiu(:, :, 1), 'U', 1., paeiv(:, :, 1), 'V', 1.)
    !$ACC KERNELS
    DO jk = 2, jpkm1
      paeiu(:, :, jk) = paeiu(:, :, 1) * umask(:, :, jk)
      paeiv(:, :, jk) = paeiv(:, :, 1) * vmask(:, :, jk)
    END DO
    !$ACC END KERNELS
  END SUBROUTINE ldf_eiv
  SUBROUTINE ldf_eiv_trp(kt, kit000, pun, pvn, pwn, cdtype)
    INTEGER, INTENT(IN   ) :: kt
    INTEGER, INTENT(IN   ) :: kit000
    CHARACTER(LEN = 3), INTENT(IN   ) :: cdtype
    REAL(KIND = wp), DIMENSION(jpi, jpj, jpk), INTENT(INOUT) :: pun
    REAL(KIND = wp), DIMENSION(jpi, jpj, jpk), INTENT(INOUT) :: pvn
    REAL(KIND = wp), DIMENSION(jpi, jpj, jpk), INTENT(INOUT) :: pwn
    INTEGER :: ji, jj, jk
    REAL(KIND = wp) :: zuwk, zuwk1, zuwi, zuwi1
    REAL(KIND = wp) :: zvwk, zvwk1, zvwj, zvwj1
    REAL(KIND = wp), DIMENSION(jpi, jpj, jpk) :: zpsi_uw, zpsi_vw
    IF (kt == kit000) THEN
      IF (lwp) WRITE(numout, FMT = *)
      IF (lwp) WRITE(numout, FMT = *) 'ldf_eiv_trp : eddy induced advection on ', cdtype, ' :'
      IF (lwp) WRITE(numout, FMT = *) '~~~~~~~~~~~   add to velocity fields the eiv component'
    END IF
    !$ACC KERNELS
    zpsi_uw(:, :, 1) = 0._wp
    zpsi_vw(:, :, 1) = 0._wp
    zpsi_uw(:, :, jpk) = 0._wp
    zpsi_vw(:, :, jpk) = 0._wp
    DO jk = 2, jpkm1
      DO jj = 1, jpjm1
        DO ji = 1, jpim1
          zpsi_uw(ji, jj, jk) = - r1_4 * e2u(ji, jj) * (wslpi(ji, jj, jk) + wslpi(ji + 1, jj, jk)) * (aeiu(ji, jj, jk - 1) + aeiu(ji, jj, jk)) * umask(ji, jj, jk)
          zpsi_vw(ji, jj, jk) = - r1_4 * e1v(ji, jj) * (wslpj(ji, jj, jk) + wslpj(ji, jj + 1, jk)) * (aeiv(ji, jj, jk - 1) + aeiv(ji, jj, jk)) * vmask(ji, jj, jk)
        END DO
      END DO
    END DO
    DO jk = 1, jpkm1
      DO jj = 1, jpjm1
        DO ji = 1, jpim1
          pun(ji, jj, jk) = pun(ji, jj, jk) - (zpsi_uw(ji, jj, jk) - zpsi_uw(ji, jj, jk + 1))
          pvn(ji, jj, jk) = pvn(ji, jj, jk) - (zpsi_vw(ji, jj, jk) - zpsi_vw(ji, jj, jk + 1))
        END DO
      END DO
    END DO
    DO jk = 1, jpkm1
      DO jj = 2, jpjm1
        DO ji = 2, jpim1
          pwn(ji, jj, jk) = pwn(ji, jj, jk) + (zpsi_uw(ji, jj, jk) - zpsi_uw(ji - 1, jj, jk) + zpsi_vw(ji, jj, jk) - zpsi_vw(ji, jj - 1, jk))
        END DO
      END DO
    END DO
    !$ACC END KERNELS
    IF (ln_ldfeiv_dia .AND. cdtype == 'TRA') CALL ldf_eiv_dia(zpsi_uw, zpsi_vw)
  END SUBROUTINE ldf_eiv_trp
  SUBROUTINE ldf_eiv_dia(psi_uw, psi_vw)
    REAL(KIND = wp), DIMENSION(jpi, jpj, jpk), INTENT(INOUT) :: psi_uw, psi_vw
    INTEGER :: ji, jj, jk
    REAL(KIND = wp) :: zztmp
    REAL(KIND = wp), DIMENSION(jpi, jpj) :: zw2d
    REAL(KIND = wp), DIMENSION(jpi, jpj, jpk) :: zw3d
    CALL lbc_lnk_multi(psi_uw, 'U', - 1., psi_vw, 'V', - 1.)
    !$ACC KERNELS
    zw3d(:, :, jpk) = 0._wp
    DO jk = 1, jpkm1
      zw3d(:, :, jk) = (psi_uw(:, :, jk + 1) - psi_uw(:, :, jk)) / (e2u(:, :) * e3u_n(:, :, jk))
    END DO
    !$ACC END KERNELS
    CALL iom_put("uoce_eiv", zw3d)
    !$ACC KERNELS
    DO jk = 1, jpkm1
      zw3d(:, :, jk) = (psi_vw(:, :, jk + 1) - psi_vw(:, :, jk)) / (e1v(:, :) * e3v_n(:, :, jk))
    END DO
    !$ACC END KERNELS
    CALL iom_put("voce_eiv", zw3d)
    !$ACC KERNELS
    DO jk = 1, jpkm1
      DO jj = 2, jpjm1
        DO ji = 2, jpim1
          zw3d(ji, jj, jk) = (psi_vw(ji, jj, jk) - psi_vw(ji, jj - 1, jk) + psi_uw(ji, jj, jk) - psi_uw(ji - 1, jj, jk)) / e1e2t(ji, jj)
        END DO
      END DO
    END DO
    !$ACC END KERNELS
    CALL lbc_lnk(zw3d, 'T', 1.)
    CALL iom_put("woce_eiv", zw3d)
    zztmp = 0.5_wp * rau0 * rcp
    IF (iom_use('ueiv_heattr') .OR. iom_use('ueiv_heattr3d')) THEN
      !$ACC KERNELS
      zw2d(:, :) = 0._wp
      zw3d(:, :, :) = 0._wp
      DO jk = 1, jpkm1
        DO jj = 2, jpjm1
          DO ji = 2, jpim1
            zw3d(ji, jj, jk) = zw3d(ji, jj, jk) + (psi_uw(ji, jj, jk + 1) - psi_uw(ji, jj, jk)) * (tsn(ji, jj, jk, jp_tem) + tsn(ji + 1, jj, jk, jp_tem))
            zw2d(ji, jj) = zw2d(ji, jj) + zw3d(ji, jj, jk)
          END DO
        END DO
      END DO
      !$ACC END KERNELS
      CALL lbc_lnk(zw2d, 'U', - 1.)
      CALL lbc_lnk(zw3d, 'U', - 1.)
      CALL iom_put("ueiv_heattr", zztmp * zw2d)
      CALL iom_put("ueiv_heattr3d", zztmp * zw3d)
    END IF
    !$ACC KERNELS
    zw2d(:, :) = 0._wp
    zw3d(:, :, :) = 0._wp
    DO jk = 1, jpkm1
      DO jj = 2, jpjm1
        DO ji = 2, jpim1
          zw3d(ji, jj, jk) = zw3d(ji, jj, jk) + (psi_vw(ji, jj, jk + 1) - psi_vw(ji, jj, jk)) * (tsn(ji, jj, jk, jp_tem) + tsn(ji, jj + 1, jk, jp_tem))
          zw2d(ji, jj) = zw2d(ji, jj) + zw3d(ji, jj, jk)
        END DO
      END DO
    END DO
    !$ACC END KERNELS
    CALL lbc_lnk(zw2d, 'V', - 1.)
    CALL iom_put("veiv_heattr", zztmp * zw2d)
    CALL iom_put("veiv_heattr", zztmp * zw3d)
    IF (ln_diaptr) CALL dia_ptr_hst(jp_tem, 'eiv', 0.5 * zw3d)
    zztmp = 0.5_wp * 0.5
    IF (iom_use('ueiv_salttr') .OR. iom_use('ueiv_salttr3d')) THEN
      !$ACC KERNELS
      zw2d(:, :) = 0._wp
      zw3d(:, :, :) = 0._wp
      DO jk = 1, jpkm1
        DO jj = 2, jpjm1
          DO ji = 2, jpim1
            zw3d(ji, jj, jk) = zw3d(ji, jj, jk) * (psi_uw(ji, jj, jk + 1) - psi_uw(ji, jj, jk)) * (tsn(ji, jj, jk, jp_sal) + tsn(ji + 1, jj, jk, jp_sal))
            zw2d(ji, jj) = zw2d(ji, jj) + zw3d(ji, jj, jk)
          END DO
        END DO
      END DO
      !$ACC END KERNELS
      CALL lbc_lnk(zw2d, 'U', - 1.)
      CALL lbc_lnk(zw3d, 'U', - 1.)
      CALL iom_put("ueiv_salttr", zztmp * zw2d)
      CALL iom_put("ueiv_salttr3d", zztmp * zw3d)
    END IF
    !$ACC KERNELS
    zw2d(:, :) = 0._wp
    zw3d(:, :, :) = 0._wp
    DO jk = 1, jpkm1
      DO jj = 2, jpjm1
        DO ji = 2, jpim1
          zw3d(ji, jj, jk) = zw3d(ji, jj, jk) + (psi_vw(ji, jj, jk + 1) - psi_vw(ji, jj, jk)) * (tsn(ji, jj, jk, jp_sal) + tsn(ji, jj + 1, jk, jp_sal))
          zw2d(ji, jj) = zw2d(ji, jj) + zw3d(ji, jj, jk)
        END DO
      END DO
    END DO
    !$ACC END KERNELS
    CALL lbc_lnk(zw2d, 'V', - 1.)
    CALL iom_put("veiv_salttr", zztmp * zw2d)
    CALL iom_put("veiv_salttr", zztmp * zw3d)
    IF (ln_diaptr) CALL dia_ptr_hst(jp_sal, 'eiv', 0.5 * zw3d)
  END SUBROUTINE ldf_eiv_dia
END MODULE ldftra


! RUN: %S/test_errors.sh %s %t %flang_fc1
! REQUIRES: shell
! Confirm enforcement of constraints and restrictions in 7.7
! C7107, C7108, C7109

subroutine bozchecks
  ! Type declaration statements
  integer :: f, realpart = B"0101", img = B"1111", resint
  logical :: resbit
  complex :: rescmplx
  real :: dbl, e
  ! C7107
  !ERROR: Invalid digit ('a') in BOZ literal 'b"110a"'
  integer, parameter :: a = B"110A"
  !ERROR: Invalid digit ('2') in BOZ literal 'b"1232"'
  integer, parameter :: b = B"1232"
  !ERROR: BOZ literal 'b"010101010101010101010101011111111111111111111111111111111111111111111111111111111111111111111111111111111111000000000000000000000000000000000000"' too large
  integer, parameter :: b1 = B"010101010101010101010101011111111111111111111&
                              &111111111111111111111111111111111111111111111&
                              &111111111111111111000000000000000000000000000&
                              &000000000"
  ! C7108
  !ERROR: Invalid digit ('8') in BOZ literal 'o"8"'
  integer :: c = O"8"
  !ERROR: Invalid digit ('a') in BOZ literal 'o"a"'
  integer :: d = O"A"

  ! C7109
  !    A) can appear only in data statement
  !    B) Argument to intrinsics listed from 16.9 below
  !       BGE, BGT, BLE, BLT, CMPLX, DBLE, DSHIFTL,
  !       DSHIFTR, IAND, IEOR, INT, IOR, MERGE_BITS, REAL

  ! part A
  data f / Z"AA" / ! OK
  !ERROR: DATA statement value could not be converted to the type 'COMPLEX(4)' of the object 'rescmplx'
  data rescmplx / B"010101" /
  ! part B
  resbit = BGE(B"0101", B"1111")
  resbit = BGT(Z"0101", B"1111")
  resbit = BLE(B"0101", B"1111")
  resbit = BLT(B"0101", B"1111")

  res = CMPLX (realpart, img, 4)
  res = CMPLX (B"0101", B"1111", 4)

  dbl = DBLE(B"1111")
  dbl = DBLE(realpart)

  !ERROR: Typeless (BOZ) not allowed for both 'i=' & 'j=' arguments
  dbl = DSHIFTL(B"0101",B"0101",2)
  !ERROR: Typeless (BOZ) not allowed for both 'i=' & 'j=' arguments
  dbl = DSHIFTR(B"1010",B"1010",2)
  dbl = DSHIFTL(B"0101",5,2) ! OK
  dbl = DSHIFTR(B"1010",5,2) ! OK

  !ERROR: Typeless (BOZ) not allowed for both 'i=' & 'j=' arguments
  resint = IAND(B"0001", B"0011")
  resint = IAND(B"0001", 3)

  !ERROR: Typeless (BOZ) not allowed for both 'i=' & 'j=' arguments
  resint = IEOR(B"0001", B"0011")
  resint = IEOR(B"0001", 3)

  resint = INT(B"1010")

  !ERROR: Typeless (BOZ) not allowed for both 'i=' & 'j=' arguments
  res = IOR(B"0101", B"0011")
  res = IOR(B"0101", 3)

  res = MERGE_BITS(13,3,11)
  res = MERGE_BITS(B"1101",3,11)
  !ERROR: Typeless (BOZ) not allowed for both 'i=' & 'j=' arguments
  res = MERGE_BITS(B"1101",B"0011",11)
  !ERROR: Typeless (BOZ) not allowed for both 'i=' & 'j=' arguments
  res = MERGE_BITS(B"1101",B"0011",B"1011")
  res = MERGE_BITS(B"1101",3,B"1011")

  res = REAL(B"1101")
end subroutine


text:
	SUBROUTINE TB_USTN  ( lun, stid, stnnam, istnm, stat, coun,
     +			      slat, slon, selv, ispri, iregn, iinst,
     +			      iret )
C************************************************************************
C* TB_USTN								*
C*									*
C* This subroutine reads the next record from the GEMPAK station table.	*
C*									*
C* TB_USTN  ( LUN, STID, STNNAM, ISTNM, STAT, COUN, SLAT, SLON,		*
C*            SELV, ISPRI, IREGN, IINST, IRET )				*
C*									*
C* Input parameters:							*
C*	LUN		INTEGER		Logical unit number		*
C*									*
C* Output parameters:							*
C*	STID		CHAR*8		Station identifier		*
C*	STNNAM		CHAR*32		Station name			*
C*	ISTNM 		INTEGER		Station number			*
C*	STAT		CHAR*2		State				*
C*	COUN		CHAR*2  	Country				*
C*	SLAT		REAL		Station latitude		*
C*	SLON		REAL		Station longitude		*
C*	SELV		REAL		Station elevation		*
C*	ISPRI		INTEGER		Station priority		*
C*	IREGN		INTEGER		Station region number		*
C*	IINST		INTEGER		Station instrument type		*
C*	IRET		INTEGER		Return code			*
C*				  	   0 = normal return		*
C*				 	  -1 = end of file 		*
C*				 	  -2 = read error		*
C**									*
C* Log:									*
C* M. desJardins/GSFC	 3/84						*
C* M. desJardins/GSFC	10/88	GEMPAK4					*
C* K. Brill/NMC		 8/93	Change for 8-char IDs and ISPRI		*
C* J. Ator/NP12		04/96	Copy from TB_RSTN, remove 4-char ID's,	*
C*				add IREGN and IINST.			*
C************************************************************************
	INCLUDE		'GEMPRM.PRM'
	CHARACTER*(*) 	stid, stnnam, stat, coun
C*
	CHARACTER  	sid*8, sname*32, st*2, cn*2
	CHARACTER*132	buffer
C------------------------------------------------------------------------
	iret = 0
C
C*	Read the record into a buffer.
C
	READ ( lun, 500, IOSTAT = iostat ) buffer
500	FORMAT ( A )
C
C*	Read in next record.
C
	IF ( iostat .eq. 0 ) THEN
	    READ ( buffer, 1000, IOSTAT = iostat ) sid, isn, sname,
     +		   st, cn, lat, lon, ihght, ispr, irgn, iist
1000	    FORMAT (A, 1X, I6, 1X, A, 1X, A, 1X, A, 1X,
     +              I5, 1X, I6, 1X, I5, 1X, I2, I1, I2 )
	END IF
C
C*	Check for error.
C
	IF  ( iostat .ne. 0 )  THEN
	    IF  ( iostat .lt. 0 )  THEN
		iret = -1
	      ELSE
		iret = -2
	    END IF
	    stid   = ' '
	    stnnam = ' '
	    istnm  = 0
	    stat   = ' '
	    coun   = ' '
	    slat   = 0.0
	    slon   = 0.0
	    selv   = 0.0
	    ispri  = 0
	    iregn  = IMISSD
	    iinst  = IMISSD
	  ELSE
C
C*	    Set output variables.
C
	    iret   = 0
	    stid   = sid
	    istnm  = isn
	    stnnam = sname
	    stat   = st
	    coun   = cn
	    slat   = lat / 100.
	    slon   = lon / 100.
	    selv   = ihght
	    ispri  = ispr
	    iregn  = irgn
	    iinst  = iist
	END IF
C*
	RETURN
	END

ext: f

max_stars_repo_name: oxelson/gempak

max_issues_repo_name: oxelson/gempak

max_forks_repo_name: oxelson/gempak


text:
	SUBROUTINE GDXSFM ( ivcord, xgrd, rlvl, nhxs, nvxs,
     +			    vclsfc, iret)
C************************************************************************
C* GDXSFM								*
C*									*
C* This subroutine sets subsurface values to missing in a cross		*
C* section grid.							*
C*									*
C* GDXSFM  ( IVCORD, XGRD, RLVL, NHXS, NVXS, VCLSFC, IRET )		*
C*									*
C* Input parameters:							*
C*									*
C*	IVCORD		  INTEGER	Vertical coordinate type	*
C*					  1 = PRES			*
C*				          2 = THTA			*
C*					  3 = HGHT			*
C*	XGRD (NHXS, NVXS) REAL		Array of cross section values	*
C*      RLVL (NVXS)	  REAL		Vertical levels in grid		*
C*	NHXS		  INTEGER	Number of xsect pts in horiz.	*
C*      NVXS              INTEGER       Number of xsect pts in vert.	*
C*      VCLSFC (NHXS)     REAL 		Vert coord location of sfc	*
C*									*
C* Output parameters:							*
C*	IRET		  INTEGER	Return code			*
C*					  0 = normal return		*
C**									*
C* Log:									*
C* K. F. Brill/GSC       7/89						*
C************************************************************************
	INCLUDE		'GEMPRM.PRM'
C*
	REAL		xgrd (nhxs, nvxs), vclsfc (nhxs), rlvl (nvxs)
	LOGICAL		done
	INCLUDE		'ERMISS.FNC'
C------------------------------------------------------------------------
	iret = 0
C
C* 	Set all values below the surface to missing.
C
	done = .false.
	i = 1
	DO WHILE ( .not. done )
	  test = rlvl ( i )
	  icnt = 0
	  DO ij = 1, nhxs
	    IF ( ERMISS ( vclsfc ( ij ) ) ) THEN
	    ELSE
	      IF ( ivcord .eq. 1 ) THEN
	        IF ( test .gt. vclsfc ( ij ) ) THEN
	          xgrd ( ij, i ) = RMISSD
	          icnt = icnt + 1
	        END IF
	      ELSE
	        IF ( test .lt. vclsfc ( ij ) ) THEN
	 	  xgrd ( ij, i ) = RMISSD
		  icnt = icnt + 1
	        END IF
	      END IF
	    END IF
	  END DO
	  IF ( icnt .eq. 0 ) done = .true.
	  i = i + 1
	  IF ( i .gt. nvxs ) done = .true.
        END DO
C*
	RETURN
	END


text:
subroutine  CO2FLX
  use simsphere_mod, only: ustar, uten, uaf, xlai, rzascr, rst, fco2, co, ci, &
                           frveg, ccan, raf
  implicit none

  real :: CHAX, PEX, RAIR, RROE, RAFCANOPY, RRTOT

!       INCLUDE 'modvars.h'

!      created 11/2/90; MODIFIED 06/17/91

!      CI =  120. suggested value for c4 plants
!      CI =  210  suggested value for c3 plants
!      CO =  330  50 m above the canopy

  CHAX = USTAR**2 / ( UTEN - UAF )

  PEX  = ( XLAI / 2.0 ) + 1.
  RAIR =  1.0 / CHAX + RZASCR
  RROE =  1.83
  RAFCANOPY = RAF * PEX / XLAI

!    RAFCANOPY IS RAF SCALED FOR XLAI TO A CANOPY
!    RROE IS DENSITY OF PURE CARBON DIOXIDE (KG M-3) AT 20 C
!    SUBROUTINE COMPUTES CARBON DIOXIDE FLUXES
!    WE FOLLOW CLOSELY GOUDRIAAN'S BOOK, CHAPTER 3
!      (CROP MICROMETEOROLOGY:  A SIMULATION STUDY; 1977)

!    CI IS INTERNAL CO2 CONCENTRATION IN LEAF
!      (120 VPM FOR C4 PLANT (CORN); 220 FOR C3 PLANT (BEANS))
!      - CO IS AMBIENT CO2 CONCENTRATION ABOVE CANOPY - - 330 VPM
!    MULTIPLY CONCENTRATIONS BY 10-6 TO GET FRACTIONAL AMOUNT

!   RESISTANCES FOR STOMATAL RESISTANCE (RST)
!                   LEAF BOUNDARY LAYER (RAF)
!                   TOP OF CANOPY (1/CHAX) TO 2 METERS
!                   SCREEN LEVEL (2 M) TO 50 METERS (RZASCR)
!   LEAF BOUNDARY LAYER RESISTANCE RAF CORRECTED FOR CO2 MOLECULAR
!                   DIFFUSIVITY  (S M-1): FACTOR OF 1.32
!   STOMATAL RESISTANCE CORRECTED FOR CO2 DIFFUSIVITY (FACTOR OF 1.66)
!          CALCULATE FLUX OF CO2 TO ENVIRONMENT
!          FACTORS 1.32 AND 1.66 ADJUST VALUES FROM H20 TO CO2.

  RRTOT = 1.32 * RAFCANOPY + 1.66*RST + RAIR

  FCO2 = RROE * ( CO - CI )  / RRTOT

  FCO2 = (FCO2 * FRVEG) / 0.044 ! In moles/m2/s
!
!      CALCULATE CO2 FLUX IN CANOPY FCO2; ADJUST FOR
!         PARTIAL VEGETATION FRACTION FRVEG
!      UNITS OF FLUX  (FCO2) ARE KG (CO2) PER METER SQUARED PER SECOND
!         IN CONCENTRATION (CCAM)  PARTS VOLUME PER MILLION (VPM)

!         **** NOTE SOME PAPERS LIST CO2 FLUX IN MOLES PER M2 PER S ****
!
!         *****    ONE MOLE OF CO2 IS 46 GRAMS    *****

  CCAN  =  CO - ( CO - CI ) * RAIR * FRVEG / RRTOT

!      CALCULATE CO2 CONCENTRATION JUST INSIDE CANOPY - - CCAN

  return
end

max_stars_repo_name: derekmorr/simsphere

max_issues_repo_name: derekmorr/simsphere

max_forks_repo_name: derekmorr/simsphere


text:
      ! -------------------------------------------------------------
      !
      ! This file contains various basic utility routines used to work
      ! with the DLR
      !
      ! -------------------------------------------------------------
      !
      ! Copyright (C) 2021 The Simons Foundation
      !
      ! Author: Jason Kaye
      !
      ! -------------------------------------------------------------
      !
      ! libdlr is licensed under the Apache License, Version 2.0 (the
      ! "License"); you may not use this file except in compliance with
      ! the License.  You may obtain a copy of the License at
      !
      !     http://www.apache.org/licenses/LICENSE-2.0
      !
      ! Unless required by applicable law or agreed to in writing,
      ! software distributed under the License is distributed on an "AS
      ! IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either
      ! express or implied.  See the License for the specific language
      ! governing permissions and limitations under the License.
      !
      ! -------------------------------------------------------------

      !> Get equispaced points on [0,1] in relative format
      !!
      !! For the definitions of absolute and relative formats, see the
      !! readme.
      !!
      !! Endpoints of [0,1] are included
      !!
      !! @param[in]   n number of points on [0,1]
      !! @param[out]  t equispaced points on [0,1], in relative format

      subroutine eqpts_rel(n,t)

      implicit none
      integer n
      real *8 t(n)

      integer i

      do i=1,n-1

        if (i.le.n/2) then
          t(i) = (i-1)*1.0d0/(n-1)
        else
          t(i) = -(n-i)*1.0d0/(n-1)
        endif

      enddo

      t(n) = 1.0d0

      end subroutine eqpts_rel

      !> Convert points on [0,1] from relative format to absolute format
      !!
      !! For the definitions of absolute and relative formats, see the
      !! readme.
      !!
      !! Note: converting a point from relative to absolute format will,
      !! in general, result in a loss of relative accuracy in the
      !! location of the point if the point is close to t = 1. For
      !! example, in three-digit arithmetic, the point t = 0.999111 could
      !! be stored as t* = -0.889e-3 in the relative format, but only as
      !! t = 0.999 in the absolute format.
      !!
      !! @param[in]   n     number of points to convert
      !! @param[in]   t     array of points in relative format
      !! @param[out]  tabs  array of points in absolute format

      subroutine rel2abs(n,t,tabs)

      implicit none
      integer n
      real *8 t(n),tabs(n)

      integer i

      do i=1,n

        if (t(i).lt.0.0d0) then
          tabs(i) = t(i)+1.0d0
        else
          tabs(i) = t(i)
        endif

      enddo

      end subroutine rel2abs

      !> Convert a point on [0,1] from absolute format to relative
      !! format
      !!
      !! For the definitions of absolute and relative formats, see the
      !! readme.
      !!
      !! Note: In order to specify points -- for example points at which
      !! to sample or evaluate a DLR -- in absolute format, those points
      !! must first be converted to relative format using this subroutine
      !! before using them as inputs into any other subroutines. Of
      !! course, in order to maintain full relative precision in all
      !! calculations, it is necessary to specify points in relative
      !! format from the beginning, but in most cases at most a mild loss
      !! of accuracy will result from using the absolute format.
      !!
      !! @param[in]   n     number of points to convert
      !! @param[in]   tabs  array of points in absolute format
      !! @param[out]  t     array of points in relative format

      subroutine abs2rel(n,tabs,t)

      implicit none
      integer n
      real *8 t(n),tabs(n)

      integer i

      do i=1,n

        if (t(i).gt.0.5d0.and.t(i).lt.1.0d0) then
          t(i) = tabs(i)-1.0d0
        else
          t(i) = tabs(i)
        endif

      enddo

      end subroutine abs2rel

max_stars_repo_name: Wentzell/libdlr

max_issues_repo_name: Wentzell/libdlr

max_forks_repo_name: Wentzell/libdlr


text:
! { dg-do run }
! PR 25806: Functions returning pointers to arrays
program a
    integer, target :: storage(5)
    integer :: s(3)

    print *, x(3)  ! { dg-output " *1  *2  *3" }

    if (ssum(x(3)) /= 6) call abort()

    s = 0
    s = x(3)
    if (any(s /= (/1, 2, 3/))) call abort()

contains

    function x(n) result(t)
        integer, intent(in) :: n
        integer, pointer :: t(:)
        integer :: i

        t => storage(1:n)
        t = (/ (i, i = 1, n) /)

    end function x

    integer function ssum(a)
        integer, intent(in) :: a(:)

        ssum = sum(a)

    end function ssum

end program a

max_stars_repo_name: maurizioabba/rose

max_issues_repo_name: sujankh/rose-matlab

max_forks_repo_name: sujankh/rose-matlab


text:
      subroutine binplotprep
c******************************************************************************
c     This routine concatenates two syntheses in order to make a smoothed
c     binary synthetic spectrum
c******************************************************************************

      implicit real*8 (a-h,o-z)
      include 'Atmos.com'
      include 'Factor.com'
      include 'Linex.com'
      include 'Pstuff.com'
      include 'Equivs.com'
      include 'Multistar.com'
      real*8 normalization
      integer pshift

c*****open the files
      nf7out = 21
      lscreen = 4
      array = 'RAW SYNTHESIS INPUT FOR PRIMARY'
      nchars = 31
      call infile ('input  ',nf7out,'formatted  ',0,nchars,
     .             f7out,lscreen)
      nf8out = 22
      lscreen = lscreen + 2
      array = 'RAW SYNTHESIS INPUT FOR SECONDARY'
      nchars = 33
      call infile ('input  ',nf8out,'formatted  ',0,nchars,
     .             f8out,lscreen)
      nf9out = 23
      lscreen = lscreen + 2
      array = 'RAW SYNTHESIS OUTPUT FOR COMBINED BINARY'
      nchars = 40
      call infile ('output ',nf9out,'formatted  ',0,nchars,
     .             f9out,lscreen)
      if (plotopt .ne. 0) then
         nf10out = 24
         lscreen = lscreen + 2
         array = 'SMOOTHED SYNTHESIS OUTPUT FOR COMBINED BINARY'
         nchars = 45
         call infile ('output ',nf10out,'formatted  ',0,nchars,
     .                f10out,lscreen)
         nf5out = 26
         lscreen = lscreen + 2
         array = 'POSTSCRIPT PLOT OUTPUT FOR COMBINED BINARY'
         nchars = 22
         call infile ('output ',nf5out,'formatted  ',0,nchars,
     .                f5out,lscreen)
      endif
      if (plotopt .gt. 1) then
         nf6out = 27
         lscreen = lscreen + 2
         array = 'SPECTRUM COMPARISON OUTPUT'
         nchars = 27
         call infile ('output ',nf6out,'formatted  ',0,nchars,
     .                f6out,lscreen)
      endif

c*****compute and output the parameters related to the spectrum combination:
c     the declared velocity difference, the luminosity ratio, and
c     for information, the computed fluxes of the two stars; get the user's
c     agreement
      deltawave = deltaradvel/2.9979d+5*(start+sstop)/2.
      pshift = int(deltawave/step+0.00001)
      write (array,1002) fluxprimary, fluxsecondary
      lscreen = lscreen + 2
      nchars = 55
      call putasci (nchars,lscreen)
      write (array,1003) lumratio
      lscreen = lscreen + 1
      nchars = 49
      call putasci (nchars,lscreen)
      write (array,1001) deltaradvel, deltawave, pshift
      lscreen = lscreen + 1
      nchars = 71
      call putasci (nchars,lscreen)
      write (array,*) 'ARE THESE VALUES OK (y/n)? '
      lscreen = lscreen + 2
      nchars = 27

c*****read back the header information from the individual star raw
c     synthetic spectra
99    read (nf7out,1004,end=100) array
      read (nf8out,1004) chinfo
      if     (array(1:7).eq.'Isotopi') then
         write (nf9out,1004) array
         go to 99
      elseif (array(1:7).eq.'ALL abu') then
         write (nf9out,1010) array(1:59), chinfo(54:59)
         go to 99
      elseif (array(1:7).eq.'Changin') then
         write (nf9out,1011) array(1:37), chinfo(31:37)
         go to 99
      elseif (array(1:7).eq.'element') then
         write (nf9out,1012) array(1:30), chinfo(25:30)
         go to 99
      elseif (array(1:7).eq.'MODEL: ') then
         write (nf9out,1013) array(8:43), chinfo(8:43)
         modbin(1) = array(1:80)
         modbin(2) = chinfo(1:80)
      endif
      read (nf7out,1004) array
      read (nf8out,1004) chinfo
      read (array,*) start, sstop, step
      kount = int((sstop - start + (step/4.0) )/step) + 1
      write (nf9out,1004) array

c*****read back the raw synthetic spectra from the individual stars;
c     shift by the appropriate point number, add the spectra (normalized),
c     and dump to combined raw synthetic spectrum file
      read (nf7out,1006) (y(i),i=1,kount)
      read (nf8out,1006) (z(i),i=1,kount)
      do i=1,kount
         y(i) = 1.0 - y(i)
         z(i) = (1.0 - z(i))/lumratio
      enddo
      normalization = 1.0 + 1.0/lumratio
      if (pshift .gt. 0) then
         do i=1,pshift
            dev(i) = 1.0
         enddo
         do i=pshift+1,kount
            dev(i) = (y(i)+z(i-pshift))/normalization
         enddo
      elseif (pshift .eq. 0) then
         do i=1,kount
            dev(i) = (y(i)+z(i))/normalization
         enddo
      else
         do i=1, kount-pshift
            dev(i) = (y(i)+z(i-pshift))/normalization
         enddo
         do i=kount-pshift+1,kount
            dev(i) = 1.0
         enddo
      endif
      do i=1,kount
         y(i) = 1.0 - dev(i)
      enddo
      write (nf9out,1006) (y(i),i=1,kount)
      go to 99
100   return

c*****format statements
1001  format ('INPUT DELTA VELOCITY = ', f6.1,
     .        ';  DELTA WAVELENGTH = ', f6.2,
     .        ';  POINT SHIFT = ', i4)
1002  format ('PRIMARY FLUX = ', 1pe10.3,
     .        ';  SECONDARY FLUX = ', 1pe10.3)
1003  format ('INPUT LUMINOSITY RATIO = ', 0pf7.2)
1004  format (a80)
1006  format (10f7.4)
1010  format (a59, ',', a6, ' dex')
1011  format (a37, ',', a6, ' dex')
1012  format (a30, ',', a6)
1013  format ('MODEL: ', 2a36)

      end

max_stars_repo_name: soylentdeen/MoogPy

max_issues_repo_name: soylentdeen/MoogPy

max_forks_repo_name: soylentdeen/MoogPy


text:
cdis    Forecast Systems Laboratory
cdis    NOAA/OAR/ERL/FSL
cdis    325 Broadway
cdis    Boulder, CO     80303
cdis
cdis    Forecast Research Division
cdis    Local Analysis and Prediction Branch
cdis    LAPS
cdis
cdis    This software and its documentation are in the public domain and
cdis    are furnished "as is."  The United States government, its
cdis    instrumentalities, officers, employees, and agents make no
cdis    warranty, express or implied, as to the usefulness of the software
cdis    and documentation for any purpose.  They assume no responsibility
cdis    (1) for the use of the software and documentation; or (2) to provide
cdis     technical support to users.
cdis
cdis    Permission to use, copy, modify, and distribute this software is
cdis    hereby granted, provided that the entire disclaimer notice appears
cdis    in all copies.  All modifications to this software must be clearly
cdis    documented, and are solely the responsibility of the agent making
cdis    the modifications.  If significant modifications or enhancements
cdis    are made to this software, the FSL Software Policy Manager
cdis    (softwaremgr@fsl.noaa.gov) should be notified.
cdis
cdis
cdis
cdis
cdis
cdis
cdis
      program gen_vrc_llij_lut
      implicit none
      include 'lapsparms.cmn'
      integer istatus
      Call get_laps_config('nest7grid',IStatus)
      if(IStatus.eq.1)then
         write(6,*)'LAPS Parameters obtained'
      else
          write(6,*)'IStatus = ',IStatus,'Error - Get_LAPS_Config'
          write(6,*)'Terminating LAPS-LVD. ISPAN Satellite remapping'
          stop
      end if
c

      call gen_vrc_llij_lut_sub(nx_l_cmn,ny_l_cmn,c_raddat_type)

      stop
      end
c
c-------------------------------------------------------------
      subroutine gen_vrc_llij_lut_sub(nx_l,ny_l,c_raddat_type)
c
      implicit none

      integer nx_l, ny_l
      character*3 c_raddat_type
      real*4    lat(nx_l,ny_l)
      real*4    lon(nx_l,ny_l)
      real*4    ri(nx_l,ny_l)
      real*4    rj(nx_l,ny_l)
      real*4    grid_spacing
      real*4    data(nx_l,ny_l,2)
      real*4    la1,lo1,la2,lo2,latin,lov,rlat,rlon,dx,dy
      integer   imax,jmax,nelems,nlines
      integer   i,j
      integer   istatus
      integer   len
c
c dimensions for lat/lon
c
      character*125 comment_ll(2)
      character*10 units_ll(2)
      character*3 var_ll(2)
      character*200 dir_static
c
c
c Definitions needed for acquiring LAPS latitude and longitude arrays.
c -------------------------------------------------------------------
      call get_directory('static',dir_static,len)

      var_ll(1) = 'LAT'
      var_ll(2) = 'LON'

      write(6,*)'Get LAPS lat/lon grid'
      call rd_laps_static(dir_static,'nest7grid',nx_l,ny_l, 2,
     &     var_ll, units_ll, comment_ll, data, grid_spacing,
     &     istatus)

      if(istatus.eq.1)then
         write(6,*)'LAPS lat/lon grid obtained'
         write(6,*)
         do j=1,ny_l
            do i=1,nx_l
               lat(i,j)=data(i,j,1)
               lon(i,j)=data(i,j,2)
            end do
         end do
      else
         write(6,*)'Unable to get lat/lon data'
         write(6,*)'lvd process terminating'
         stop
      end if
c
c      do i=1,n_radar_types

      if(c_raddat_type.eq.'wfo')then

         write(6,*)'Gen LUT for Conus lambert (WFO) WSI data '
c        call gen_vrc_wfo_cdf_lut(2,nx_l,ny_l,lat,lon,istatus)
c        call

         call gen_rirj_lam(imax,jmax,lat,lon,nelems,nlines
     &        ,la1,lo1,la2,lo2,latin,lov,rlat,rlon,dx,dy,ri,rj
     &        ,istatus)

         if(istatus.eq.1)then
            write(6,*)'CDF look-up table generated'
         else
            write(6,*)'Error generating look-up-table'
            goto 901
         endif

      elseif(c_raddat_type.eq.'wsi')then

         write(6,*)'Gen LUT for wsi - Cylindrical Equidistant'

         print*,'Using latlon_to_ceij routine'
         call gen_rirj_ce(nx_l,ny_l,lat,lon,istatus)

      endif

      if(istatus.eq.1)then
         write(6,*)'CDF look-up table generated'
      else
         write(6,*)'Error generating look-up-table'
         goto 901
      endif

c      enddo
c
      goto 900

901   write(6,*)'Error - genvrclut.f'

900   stop
      end


text:
Aspen Pest Management is a local Pest Control pest control company offering pest control and Termites termite inspections and treatments in Yolo County Yolo, Sacramento County Sacramento, and eastern Solano Counties.

Look for their coupons in the Little Green Coupon Machine.

20080923 09:42:24 nbsp Weve used their bimonthly (every other month) service for a couple of years.  They come by and spray around the house foundations to control ants.  Works pretty well.  They will also drop off traps for mice and so forth.  I like using a local business rather than some huge corporate chain that is not responsive. Users/RalphFinch

20100422 19:13:30 nbsp These guys are GREAT for termites! They will not make a mess, use pretty ECOcool stuff and are extremely reasonable as far as pricing. Over a few years weve only needed 2 sprays  they give great suggestions for how to use the inexpensive termite stakes and not waste money. Also, they show up on time! Users/KimPruett

20101108 15:06:08 nbsp I called for an inspection and there was no charge.  The guy who came out was super helpful, he gave me lots of advice for do it yourself pest management and then even gave me some free ant traps.  Honest customer service is a great asset in a company and Aspen seems to have lots of it.  Aspen is now the only pest control company I will ever use. Users/DagonJones

20110108 15:44:44 nbsp We completely agree with the comment by DagonJones.  Amazing and honest customer service!  Basically, our inspector told us how to resolve the issue....resulting in zero cost to us.  Wow. Users/JAD

20110802 13:18:21 nbsp The technician who came out to inspect our bed bug problem was really friendly and knowledgeable.  Users/CarlosOverstreet

20130717 14:18:35 nbsp I called for an inspection of our shed from both Terminix and Aspen. Terminix said the damage to our shed floor was from termites; however, David from Aspen told us to simply replace the floor with pressure treated lumber as the damage was from rot not termites. He then inspected the garage and found a few small termite tubules. He said he had just finished a job and it would be no problem to drill a few holes in the concrete, insert the termite poison and fill the concrete. There was no charge at all for his assistance! He was knowledgeable and provided amazing service. I highly recommend this company to anyone in Davis.  Users/ChristyannDarwent

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text:
	PROGRAM HPGL
C *** LAST REVISED ON 19-AUG-1987 11:41:56.19
C *** SOURCE FILE: [DL.GRAPHICS.LONGLIB]HPGL.FOR
C
C	CREATED:  DGL 15-JAN-1986
C
C	THIS PROGRAM PERMITS REPLOTTING OF A PRINTER GRAPHICS FILE
C	PRODUCED BY THE LONGLIB GRAPHICS LIBRARY TO A VT220 TERMINAL
C	EQUIPPED WITH AN HPGL COMPATIBLE PEN PLOTTER.  THIS VERSION
C	DOES NOT USE STRIPPING OF THE PLOT IMAGE NOR ATTEMPT TO MINIMIZE
C	PEN MOTIONS.
C
C	NOTE: SINCE THE PRINTER FILE DOES NOT MAINTAIN A COPY OF
C	WHAT DATA WENT TO THE SCREENS, ANY CTERMS USED IN THE ORIGINAL
C	PROGRAM GENERATING THE FILE WILL NOT BE REPEATED BY THIS PROGRAM.
C
C	THIS PROGRAM IS FORTRAN 77 COMPATIBLE WITH THE EXCEPTIONS:
C		1. TABS (^I) ARE USED TO INDENT LINES.
C		2. A VAX DEPENDENT ROUTINE IS USED TO GET COMMAND LINE.
C		3. INTEGER*2 USED IN REGET ROUTINE.  THIS CAN BE CHANGED
C		   TO INTEGER IF THIS IS ALSO DONE IN LONGLIB PRINTER
C		   HISTORY ROUTINES LIBRARY.
C		4. USES ICHAR TO EXTRACT ASCII VALUE FROM CHARACTER VARIABLE.
C
	COMMON /TERM/ITERM
	CHARACTER*80 NAME
C
C	MACHINE DEPENDENT METHOD OF GETING THE FILE NAME FROM THE
C	RUN COMMAND LINE
C
	IERR=LIB$GET_FOREIGN(NAME,,IFLAG)
C
C	IF FILE NAME WAS NOT PRESENT ON ARGUMENT LINE THEN USE DEFAULT NAME
C
	IF (NAME.EQ.' ') NAME='FOR003.DAT'
C
C	DETERMINE TERMINAL TYPE
C
	WRITE(*,10)
10	FORMAT('$ENTER TERMINAL TYPE: (0=PLAIN VT220, 1=VT100+SELNAR) ')
	READ(*,*) ITERM
C
C	TURN PRINTER PORT ON TERMINAL
C
	CALL PLOTAUTO
C
C	CONVERT LONGLIB META FILE TO HPGL PLOT COMMANDS AND WRITE TO TERMINAL
C
	CALL HPGLPLOT(NAME)
C
C	TURN PRINTER PORT OFF ON TERMINAL
C
	CALL PLOTAOFF
	STOP
	END
C
	SUBROUTINE PLOTAOFF
C
C	DISABLE PRINTER PORT ON TERMINAL
C
	COMMON /TERM/ITERM
	CHARACTER*1 ESC,ETB
	PARAMETER (ESC=CHAR(27))
	PARAMETER (ETB=CHAR(23))
	IF (ITERM.EQ.0) THEN
C	VT220
		WRITE (*,1) ESC
1		FORMAT(X,A1,'[?4i')
	ELSE IF (ITERM.EQ.1) THEN
C	VT100+SELNAR GR100
		WRITE (*,1) ETB,ESC
2		FORMAT('+',2A1,'2')
	ENDIF
	RETURN
	END
C
C
	SUBROUTINE PLOTAUTO
C
C	ENABLE PRINTER PORT ON TERMINAL
C
	COMMON /TERM/ITERM
	CHARACTER*1 ESC
	PARAMETER (ESC=CHAR(27))
	IF (ITERM.EQ.0) THEN
C	VT220
	 	WRITE (*,1) ESC
1		FORMAT(X,A1,'[?5i')
	ELSE IF (ITERM.EQ.1) THEN
C	VT100+SELNAR GR100
		WRITE (*,2) ESC,ESC
2		FORMAT(X,A1,'1',A1,';<')
	ENDIF
	RETURN
	END
C
C
	SUBROUTINE HPGLPLOT(NAME)
C
C	MAIN ROUTINE TO REPLOT A PRINTER GRAPHICS FILE TO AN HPGL COMPATIBLE
C	DEVICE FOR THE LONGLIB GRAPHICS LIB.
C
C	INPUT:
C	NAME	(CHARACTER)	FILE NAME FOR PRINTER HISTORY FILE
C
	CHARACTER*(*) NAME
	CHARACTER*1 ANS
	INTEGER LINET(7)
	DATA LINET/-1,1,4,2,3,6,5/
	LOGICAL PENUP
	DATA PENUP/.TRUE./
C
C	INITIALIZE PLOTTER
C
	WRITE (*,23)
23	FORMAT(' IN;SP1;PU;LT-1;')
C
C	OPEN PRINTER HISTORY FILE
C
	MP=999
	OPEN(UNIT=2,FILE=NAME,FORM='UNFORMATTED',STATUS='OLD',
     $		READONLY,ERR=99)
C
C	TOP OF COMMAND READ LOOP
C
1000	CALL REGET(M1,M2,M3,MP,2)
	IF (M3.EQ.11) GOTO 110
	SY=7200./FLOAT(M1)/7.92
	SX=10000./FLOAT(M2)/11.0
	IF (M1.EQ.0.OR.M2.EQ.0) GOTO 110
  1	CONTINUE
	CALL REGET(M1,M2,M3,MP,2)
	Y1=M1*SY
	X1=M2*SX
C
C	I3 IS THE COMMAND CODE (EQUIVALENT TO 3RD ARGUMENT OF PLOT)
C
	I3=M3
C
C	IF I3=999 THEN END OF FILE ENCOUNTERED
C
	IF (M3.GE.999) GOTO 999
C
C	EXECUTE COMMAND
C
	GOTO (10,20,20,10,10,10,10,10,20,100,110) I3
10	GOTO 1
C
20	CONTINUE
C
C	CALL PLOT(X1,Y1,I3)		! EQUIV LONGLIB COMMAND
C
	IF (PENUP) THEN
		IF (I3.EQ.2) WRITE(*,3000)
3000		FORMAT(' PD;')
	ELSE
		IF (I3.EQ.3) WRITE(*,3001)
3001		FORMAT(' PU;')
	ENDIF
	IF (I3.EQ.2) PENUP=.FALSE.
	IF (I3.EQ.3) PENUP=.TRUE.
	WRITE(*,43) IFIX(X1+.5),IFIX(Y1+.5)
43	FORMAT(X,'PA',I6,',',I6,';')
	GOTO 1
C
100	CONTINUE
C
C	CALL NEWPAGE		! EQUIV LONGLIB COMMAND
C
	WRITE(*,3002)
3002	FORMAT(' PU;SP0;')
C
C	DISABLE TERMINAL PRINTER PORT
C
	CALL PLOTAOFF
C
C	PROMPT TERMINAL FOR PAGE CHANGE
C
	WRITE(*,44)
44	FORMAT('$Change paper, <return> when ready, S to stop ')
	READ (*,3500,END=99) ANS
3500	FORMAT(A1)
C
C	RE-ENABLE TERMINAL PRINTER PORT
C
	CALL PLOTAUTO
C
C	IS REPLY IS ASCII S THEN STOP
C
	IF (ICHAR(ANS).EQ.115.OR.ICHAR(ANS).EQ.83) RETURN
C
	GOTO 1
999	GOTO (1000,1001,1002,1003) I3-999
2000	GOTO 1
C
C	LINE TYPE
C
1001	CONTINUE
	IOLD=IABS(M1)
1550	IOLD=MOD(IOLD,7)+1
	IOLD=LINET(IOLD)
	WRITE (*,1555) IOLD
1555	FORMAT(X,'LT',I2',1.2;')
	GOTO 1
C
C	SET PEN COLOR
C
1002	CONTINUE
C
C	IF (M2.NE.0) CALL PLOT(FLOAT(M2),0.,0)	! EQIV LONGLIB NEW COLOR
C
	IF (M2.NE.0) THEN
		I=M2
		IF (M2.EQ.255) I=1
		IF (M2.EQ.10) I=1
		WRITE(*,3020) I
3020		FORMAT(' PU;SP',I2,';')
	END IF
	GOTO 1
C
C	DETERMINE LINE WIDTH
C
1003	I3=M1*10+IOLD
C
C	CALL NEWPEN(I3)		! EQUIV LONGLIB COMMAND
C
	GOTO 1
C
C	CLOSE FILE AND TERMINATE PLOTTING
C
110	CLOSE(2)
 99	WRITE(*,3030)
3030	FORMAT(' PU;SP0;')
	RETURN
	END
C
	SUBROUTINE REGET(M1,M2,M3,MP,ILU)
C
C	READ DATA FROM LONGLIB PRINTER HISTORY META FILE
C
	INTEGER*2 M(128)
	MP=MP+3
	IF (MP.GT.128) THEN
		READ (ILU,ERR=99) M
		MP=3
	ENDIF
	M3=M(MP)
	M2=M(MP-1)
	M1=M(MP-2)
	IF (M3.EQ.999) GOTO 99
	RETURN
99	M3=11
	RETURN
	END

ext: for


text:
C Copyright(C) 2009-2017 National Technology & Engineering Solutions of
C Sandia, LLC (NTESS).  Under the terms of Contract DE-NA0003525 with
C NTESS, the U.S. Government retains certain rights in this software.
C
C Redistribution and use in source and binary forms, with or without
C modification, are permitted provided that the following conditions are
C met:
C
C     * Redistributions of source code must retain the above copyright
C       notice, this list of conditions and the following disclaimer.
C
C     * Redistributions in binary form must reproduce the above
C       copyright notice, this list of conditions and the following
C       disclaimer in the documentation and/or other materials provided
C       with the distribution.
C     * Neither the name of NTESS nor the names of its
C       contributors may be used to endorse or promote products derived
C       from this software without specific prior written permission.
C
C THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
C "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
C LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
C A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
C OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
C SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
C LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
C DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
C THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
C (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
C OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

C=======================================================================
      INTEGER FUNCTION IXVW (INCEMP, IVWNUM)
C=======================================================================

C   --*** IXVW *** (MESH) Return the index of the Nth defined view
C   --   Written by Amy Gilkey - revised 10/08/87
C   --
C   --IXVW returns the index of IVWNUMth defined view, which may include
C   --empty views.  If there are not IVWNUM defined views, zero is returned.
C   --
C   --Parameters:
C   --   INCEMP - IN - a view with an empty mode are returned iff true
C   --   IVWNUM - IN - the view number
C   --
C   --Common Variables:
C   --   Uses MSHDEF of /MSHOPT/

      include 'mshopt.blk'

      LOGICAL INCEMP

      IVW = 0
      DO 100 IXVW = 1, 4
         IF (MSHDEF(IXVW) .NE. 'NONE') THEN
            IF (INCEMP .OR. (MSHDEF(IXVW) .NE. 'EMPTY')) THEN
               IVW = IVW + 1
               IF (IVW .GE. IVWNUM) GOTO 110
            END IF
         END IF
  100 CONTINUE
      IXVW = 0
  110 CONTINUE

      RETURN
      END

max_stars_repo_name: mathstuf/seacas

max_stars_repo_licenses: BSD-3-Clause NetCDF Zlib MIT

max_issues_repo_name: mathstuf/seacas

max_issues_repo_licenses: BSD-3-Clause NetCDF Zlib MIT

max_forks_repo_name: mathstuf/seacas

max_forks_repo_licenses: BSD-3-Clause NetCDF Zlib MIT


text:
Termites are insects from the order Isoptera. These insects are important economically, since they cause damage to wooden structures. Although not closely related to honeybees bees, wasps and ants ants, members of the order Hymenoptera, these social insects display a similar pattern or social behavior. Some believe that these insects are most closely related to cockroaches, and, to quote a rhyme,

A primitive cockroach knocked on wood
tasted it and found it good
and that is why your aunt May
fell through the parlor floor today.

Termite colonies in the soil or wood are composed of several different castes. Workers generally excavate new tunnels and chambers and do other housekeeping tasks in the colonies. Soldiers use fearsome jaws, or, in some species, use tubes on the heads to shoot fluids at intruders. Some species of termites have members that are specialized to plug entrances with their huge corklike heads. The reproductive castes produce the offspring to expand the hive. In some cases, the worker castes are not adults, but juveniles, and some of them are capable of developing into reproductive individuals.

Termites are only able to digest the cellulose in wood because of protozoans present in the gastric caecum, a pit in the digestive tract of the insects. These protozoans are not transmitted at birth; termites often eat the feces of other termites to obtain them, or place their mouthparts onto the terminal end of the abdomen of other termites to introduce them to their digestive tract.

Within the order of Isoptera, there are two economically important families, the Kalotermatidae and the Rhinotermatidae, which are of most concern to us. Kalotermatidae are the drywood, dampwood and powderpost termites, and these build their colonies exclusively in the wood. Rhinotermatidae, on the other hand, are the subterranean termites which build their nest partly in wood and partly in the soil. This is an important distinction, because fumigating the infested wood works well for Kalotermatidae, but doesnt work as well for Rhinotermatidae. There other families of termites, which include the impressive mound building termites which create elaborate solid structures which dot the landscapes of arid regions. It is said that a bit of the soil used to make these mounds is an effective repellant for other insects.


! { dg-do run }
! { dg-additional-options "-w" }

! subroutine reduction with private and firstprivate variables

program reduction
  integer, parameter    :: n = 100
  integer               :: i, j, vsum, cs, arr(n)

  call redsub_private (cs, n, arr)
  call redsub_bogus (cs, n)
  call redsub_combined (cs, n, arr)

  vsum = 0

  ! Verify the results
  do i = 1, n
     vsum = i
     do j = 1, n
        vsum = vsum + 1;
     end do
     if (vsum .ne. arr(i)) call abort ()
  end do
end program reduction

! This subroutine tests a reduction with an explicit private variable.

subroutine redsub_private(sum, n, arr)
  integer :: sum, n, arr(n)
  integer :: i, j, v

  !$acc parallel copyout (arr)
  !$acc loop gang private (v)
  do j = 1, n
     v = j

     !$acc loop vector reduction (+:v)
     do i = 1, 100
        v = v + 1
     end do

     arr(j) = v
  end do
  !$acc end parallel

  ! verify the results
  do i = 1, 10
     if (arr(i) .ne. 100+i) call abort ()
  end do
end subroutine redsub_private

! Bogus reduction on a firstprivate variable.  The results do
! survive the parallel region.  The goal here is to ensure that gfortran
! doesn't ICE.

subroutine redsub_bogus(sum, n)
  integer :: sum, n, arr(n)
  integer :: i

  !$acc parallel firstprivate(sum)
  !$acc loop gang worker vector reduction (+:sum)
  do i = 1, n
     sum = sum + 1
  end do
  !$acc end parallel
end subroutine redsub_bogus

! This reduction involving a firstprivate variable yields legitimate results.

subroutine redsub_combined(sum, n, arr)
  integer :: sum, n, arr(n)
  integer :: i, j

  !$acc parallel copy (arr) firstprivate(sum)
  !$acc loop gang
  do i = 1, n
     sum = i;

     !$acc loop reduction(+:sum)
     do j = 1, n
        sum = sum + 1
     end do

     arr(i) = sum
  end do
  !$acc end parallel
end subroutine redsub_combined