  subroutine rdmyhre_hno3(nopr,iout,indir,outdir,lstrn,iset,&
  lset,llset,nlines,wcm,wavelength,rndat,ridat,iwave)

! *****
! Read in the data from the netCDF file
! Will use netCDF f90 routines
  use netcdf

! param1.inc gives the maximum size of indices and extinction vectors
  include 'param1.inc'

! ****
! See rdinoutdir.f90
  character(len=60) :: indir,outdir
  integer :: nopr,iout,iset

! See main.f90
  integer :: nlines,iwave
  real :: wcm(nwavemax),wavelength(nwavemax)
  real :: rndat(nwavemax),ridat(nwavemax)
  character(len=100) :: lstrn,lset,filnm

! Used here
! Notice that the exact number of elements in rnval,rival works
  integer,parameter :: ndat=6051,mtempmax=6,mcases=3
  integer :: i,j1,j2,k1,k2,mn,nc
  integer :: ncases,ntext,n70
  integer :: start(1),count(1),stride(1)
  integer :: start2(2),count2(2),stride2(2)
  integer :: start3(3),count3(3),stride3(3)
  real :: a1,a2
  real :: wghthno3
! Note that the ndat,mtemp is opposite to the ncdump dump
  real :: numlinesr(mcases),numlinesi(mcases)
  real :: temps(mcases,mtempmax)
  real :: ntemps(mcases)
  real :: whno3(mcases)
  real :: wcmr(ndat,mcases),wavelengthr(ndat,mcases)
  real :: wcmi(ndat,mcases),wavelengthi(ndat,mcases)
  real :: rn(ndat,mcases,mtempmax)
  real :: ri(ndat,mcases,mtempmax)
  character(len=30) :: fil
  character (len=16) :: llset

! ***
! Read in the data from the netCDF file
   lset=lstrn

! ***
! The input netCDF file
!     '123456789012345678901234567890'
  fil='myhre_hno3.nc                 '
   noprf=1
  call getfilnm(noprf,iout,fil,indir,filnm)
  if (noprf .eq. 1) then
!  stop
  end if

! ***
  write(0,fmt=160)
  160 format(2x,"Specify wghthno3 value (e.g. 30, 54, 64) ")
  read(5,*) wghthno3

! The compound string
!       '1234567890123456'
  llset=' Myhre HNO3     '

  write(0,fmt=170) wghthno3
  170 format(2x," wghthno3 ",f10.4)

! ***
! Open the netCDF file
  istat=nf90_open(filnm,nf90_nowrite,id1)
   call handle(iout,istat)

  if (nopr .eq. 1) then
   write(iout,fmt=100) id1,filnm
   100 format(/,2x,"rdmyhre_hno3: id1,filnm ",i8,2x,a100)
  end if

! ***
! Get dimension ids (work with variables part of ncdump)
  istat=nf90_inq_varid(id1,'ncases',idncases)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'ntempmax',idntempmax)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'nlinesr',idnlinesr)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'nlinesi',idnlinesi)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'ntext',idntext)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'n70',idn70)
   call handle(iout,istat)

  if (nopr .eq. 1) then 
   write(iout,fmt=180) idncases,idntempmax,idnlinesr
   180 format(2x,"idncases,idntempmax,idnlinesr ",/,2x,3(1x,i4))
   write(iout,fmt=190) idnlinesi,idntext,idn70
   190 format(2x,"idnlinesi,idntext,idn70 ",/,2x,3(1x,i4))
  end if

! *****
! Get the dimensions
  istat=nf90_get_var(id1,idncases,ncases)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idntempmax,ntempmax)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idnlinesr,nlinesr)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idnlinesi,nlinesi)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idntext,ntext)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idn70,n70)
   call handle(iout,istat)

  if (nopr .eq. 1) then
   write(iout,fmt=200) ncases,ntempmax,nlinesr,nlinesi,ntext,n70
   200 format(/,2x,"ncases,ntempmax,nlinesr,nlinesi,ntext,n70 ",/,&
   2x,6(1x,i4))
  end if

! *****
! Get variable ids
  istat=nf90_inq_varid(id1,'numlinesr',idnumlinesr)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'numlinesi',idnumlinesi)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'temps',idtemps)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'ntemps',idntemps)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'whno3',idwhno3)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'wcmr',idwcmr)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'wavelengthr',idwavelengthr)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'rn',idrn)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'wcmi',idwcmi)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'wavelengthi',idwavelengthi)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'ri',idri)
   call handle(iout,istat)

  if (nopr .eq. 1) then
   write(iout,fmt=210) idnumlinesr,idnumlinesi,idtemps,idntemps,idwhno3
   210 format(/,2x,"idnumlinesr,idnumlinesi,idtemps,idntemps,idwhno3, ",/,&
   2x,5(1x,i4))
   write(iout,fmt=220) idwcmr,idwavelengthr,idrn
   220 format(/,2x,"idwcmr,idwavelengthr,idrn ",/,2x,3(1x,i4))
   write(iout,fmt=230) idwcmi,idwavelengthi,idri
   230 format(/,2x,"idwcmi,idwavelengthi,idri ",/,2x,3(1x,i4))
  end if

! ****
! Get the data

! The number of lines for real indices
   start(1)=1
   count(1)=ncases
   stride(1)=1
  istat=nf90_get_var(id1,idnumlinesr,numlinesr,start,count,stride)

! The number of lines for imaginary indices
   start(1)=1
   count(1)=ncases
   stride(1)=1
  istat=nf90_get_var(id1,idnumlinesi,numlinesi,start,count,stride)

! ****
! The temperatures
   start2(1)=1
   start2(2)=1
   count2(1)=mcases
   count2(2)=mtempmax
   stride2(1)=1
   stride2(2)=1
  istat=nf90_get_var(id1,idtemps,temps,start2,count2,stride2)

! The number of temperatures for each case
   start(1)=1
   count(1)=ncases
   stride(1)=1
  istat=nf90_get_var(id1,idntemps,ntemps,start,count,stride)

! The whno3 values for each case
   start(1)=1
   count(1)=ncases
   stride(1)=1
  istat=nf90_get_var(id1,idwhno3,whno3,start,count,stride)

! ****
! The wcmr values
   start2(1)=1
   start2(2)=1
   count2(1)=ndat
   count2(2)=mcases
   stride2(1)=1
   stride2(2)=1
  istat=nf90_get_var(id1,idwcmr,wcmr,start2,count2,stride2)

! The wavelengthr values
   start2(1)=1
   start2(2)=1
   count2(1)=ndat
   count2(2)=mcases
   stride2(1)=1
   stride2(2)=1
  istat=nf90_get_var(id1,idwavelengthr,wavelengthr,start2,count2,stride2)

! Read in the real indices
! Note that count3(1) and count3(2) etc are opposite to the ncdump dump
   start3(1)=1
   start3(2)=1
   start3(3)=1
   count3(1)=ndat
   count3(2)=mcases
   count3(3)=mtempmax
   stride3(1)=1
   stride3(2)=1
   stride3(3)=1
  istat=nf90_get_var(id1,idrn,rn,start3,count3,stride3)
   call handle(iout,istat)

! ****
! The wcmi values
   start2(1)=1
   start2(2)=1
   count2(1)=ndat
   count2(2)=mcases
   stride2(1)=1
   stride2(2)=1
  istat=nf90_get_var(id1,idwcmi,wcmi,start2,count2,stride2)
   call handle(iout,istat)

! The wavelengthi values
   start2(1)=1
   start2(2)=1
   count2(1)=ndat
   count2(2)=mcases
   stride2(1)=1
   stride2(2)=1
  istat=nf90_get_var(id1,idwavelengthi,wavelengthi,start2,count2,stride2)
   call handle(iout,istat)

! Read in the imaginary indices
! Note that count3(1) and count3(2) etc are opposite to the ncdump dump
   start3(1)=1
   start3(2)=1
   start3(3)=1
   count3(1)=ndat
   count3(2)=mcases
   count3(3)=mtempmax
   stride3(1)=1
   stride3(2)=1
   stride3(3)=1
  istat=nf90_get_var(id1,idri,ri,start3,count3,stride3)
   call handle(iout,istat)

! ****
! Exit define mode
  istat=nf90_enddef(id1)

! Close the data set
  istat=nf90_close(id1)

! *****
! Find the indices for the interpolation

! The whno3 values
  if (nopr .eq. 1) then
   write(iout,fmt=320) 
   320 format(/,2x,"rdmyhre_hno3: i,whno3(i)")
   do i=1,mcases
    write(iout,fmt=330) i,whno3(i)
    330 format(2x,i4,2x,f10.4)
   end do
  end if

! Get case (whno3 value)
    j1=-99
   if (wghthno3 .le. whno3(1)) then
    j1=1
    go to 400
   end if
   if (wghthno3 .ge. whno3(3)) then
    j1=3
    go to 400
   end if
  do i=1,mcases-1
    i2=i+1
   if ((wghthno3 .ge. whno3(i)) .and. (wghthno3 .lt. whno3(i2))) then
    j1=i
    go to 400
   end if
  end do
  400 mn=0

  if (j1 .lt. 0) then
   write(iout,fmt=350) j1
   350 format(/,2x,"rdmyhre_hno3: j1 lt 0 ",i4," will stop")
   stop
  end if

  numtemp=ntemps(j1)

  if (nopr .eq. 1) then
   write(iout,fmt=410) wghthno3,whno3(j1),j1,numtemp
   410 format(/,2x,"rdmyhre_hno3: wghthno3,whno3(j1),j1,numtemp ",/,&
   2x,2(1x,f10.4),2(1x,i4))
  end if

! *****
! Specify the temperature
  if (nopr .eq. 1) then
   write(iout,fmt=420) 
   420 format(/,2x,"rdmyhre_hno3: i,temps(j1,i)")
   do i=1,numtemp
    write(iout,fmt=430) i,temps(j1,i)
    430 format(2x,i4,2x,f10.4)
   end do
  end if

  write(0,fmt=425)
  425 format(2x,"Specify temp from these values ")
  do i=1,numtemp
   write(0,fmt=435) temps(j1,i)
   435 format(2x,f10.4)
  end do

  write(0,fmt=437)
  437 format(2x,"Specify temp ")
  read(5,*) temp

   k1=-99
  if (temp .le. temps(j1,numtemp)) then
   k1=numtemp
   go to 439
  end if
  if (temp .ge. temps(j1,1)) then
   k1=1
   go to 439
  end if
  do i=1,numtemp-1
    i2=i+1
   if ((temp .le. temps(j1,i)) .and. (temp .gt. temps(j1,i2))) then 
    k1=i
    go to 439
   end if
  end do
  439 mn=0

  if (nopr .eq. 1) then
   write(iout,fmt=440) temp,k1,temps(j1,k1)
   440 format(/,2x,"rdmyhre_hno3: temp,k1,temps(j1,k1)  ",/,&
   2x,1(1x,f10.4),1(1x,i4),1(1x,f10.4))
  end if

! *****
! Put values into the output arrays
! rnval(nlinesr,ncases,ntempmax)
! rival(nlinesr,ncases,ntempmax)
   nlines=nlinesr
  do i=1,nlines
   wcm(i)=wcmr(i,j1)
   wavelength(i)=wavelengthr(i,j1)
   a1=rn(i,j1,k1)
   a2=ri(i,j1,k1)
   if (a1 .le. 0.0) then 
    a1=0.0
   endif
   if (a2 .le. 0.0) then 
    a2=0.0
   endif
   rndat(i)=a1
   ridat(i)=a2
  end do

! *****
  if (nopr .eq. 1) then
    write(iout,fmt=800) temp,wghthno3
    800 format(/,2x,"rdmyhre_hno3: temp,wghthno3 ",/,2x,2(1x,f10.4),/,&
    2x,"rdmyhre_hno3: i,wcm,wavelength,rndat,ridat")
   do i=1,nlines
    write(iout,fmt=850) i,wcm(i),wavelength(i),rndat(i),ridat(i)
    850 format(2x,i4,2(1x,f12.4),1p,2(1x,e10.3))
   end do
  end if

! ******
  return
  end
