  subroutine rdstagg_carbon(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
  integer,parameter :: ndat=10,ncomp=4,mtempmax=4
  integer :: id1,istat,noprf,mset,itemp,i1,i2
  integer :: idwcm,idwave,idrn,idri,idtk
  integer :: start(1),count(1),stride(1)
  integer :: start2(2),count2(2),stride2(2)
  integer :: npyro,namorphous,npropane,nsoot
  integer :: ntpyro,ntamorphous,ntpropane,ntsoot
  integer :: ntext,n70,ntk
  real :: wcmdat(ndat),wavedat(ndat)
  real :: rnval(ndat,mtempmax),rival(ndat,mtempmax)
  real :: temps(mtempmax)
  character(len=30) :: fil
  character (len=16) :: files(ncomp),llset

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

! ***
! The input netCDF file
!     '123456789012345678901234567890'
  fil='stagg_carbon.nc               '

! Obtain the full pathname of the input netCDF file
   noprf=0
  call getfilnm(noprf,iout,fil,indir,filnm)
  if (noprf .eq. 1) then
   stop
  end if

! ***
! 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,"rdstagg_carbon: id1,filnm ",i8,2x,a100)
  end if

! ***
! Get dimension ids (work with variables part of ncdump)
  istat=nf90_inq_varid(id1,'ntext',idntext)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'n70',idn70)
   call handle(iout,istat)

  istat=nf90_inq_varid(id1,'npyro',idnpyro)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'namorphous',idnamorphous)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'npropane',idnpropane)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'nsoot',idnsoot)
   call handle(iout,istat)

  istat=nf90_inq_varid(id1,'ntpyro',idntpyro)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'ntamorphous',idntamorphous)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'ntpropane',idntpropane)
   call handle(iout,istat)
  istat=nf90_inq_varid(id1,'ntsoot',idntsoot)
   call handle(iout,istat)

! *****
! Get the dimensions
  istat=nf90_get_var(id1,idnpyro,npyro)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idnamorphous,namorphous)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idnpropane,npropane)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idnsoot,nsoot)
   call handle(iout,istat)

  istat=nf90_get_var(id1,idntpyro,ntpyro)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idntamorphous,ntamorphous)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idntpropane,ntpropane)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idntsoot,ntsoot)
   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=110) npyro,namorphous,npropane,nsoot
   110 format(/,2x,"npyro,namorphous,npropane,nsoot ",/,&
   2x,4(1x,i4))
   write(iout,fmt=115) ntpyro,ntamorphous,ntpropane,ntsoot
   115 format(/,2x,"ntpyro,ntamorphous,ntpropane,ntsoot ",/,&
   2x,4(1x,i4))
   write(iout,fmt=120) ntext,n70
   120 format(/,2x,"ntext,n70",/,&
   2x,2(1x,i4))
  end if

! ***
! Obtain the listing of the compounds
! Will ask the user to choose which compound to work with

! The different compositions
!          '1234567890123456'
  files(1)='Pyro            '
  files(2)='Amorphous       '
  files(3)='Propane         '
  files(4)='Soot            '

! *
! Write out the compounds included in the netCDF file
  write(0,fmt=140)
  140 format(2x," Will select mset value for a specific compound")
  write(0,fmt=145)
  145 format(2x," i, compound")

  do i=1,ncomp
   write(0,fmt=150) i,files(i)
   150 format(2x,i3,2x,a16)
  end do

  write(0,fmt=160)
  160 format(2x,"Specify mset (1-4) e.g. mset=4 for Soot ")

  read(5,*) mset
  mset=int(mset)

! The compound string
  llset=files(mset)

  write(0,fmt=170) mset,llset
  170 format(2x," mset,llset ",i4,2x,a16)

! ***
! Obtain the specific indices

   if (mset .eq. 1) then 
     nlines=npyro
     ntk=ntpyro
    istat=nf90_inq_varid(id1,'wcmpyro',idwcm)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'wmicpyro',idwave)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'rnpyro',idrn)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'ripyro',idri)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'tkpyro',idtk)
     call handle(iout,istat)
   end if

   if (mset .eq. 2) then
     nlines=namorphous
     ntk=ntamorphous
    istat=nf90_inq_varid(id1,'wcmamorphous',idwcm)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'wmicamorphous',idwave)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'rnamorphous',idrn)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'riamorphous',idri)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'tkamorphous',idtk)
     call handle(iout,istat)
   end if

   if (mset .eq. 3) then
     nlines=npropane
     ntk=ntpropane
    istat=nf90_inq_varid(id1,'wcmpropane',idwcm)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'wmicpropane',idwave)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'rnpropane',idrn)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'ripropane',idri)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'tkpropane',idtk)
     call handle(iout,istat)
   end if

   if (mset .eq. 4) then 
     nlines=nsoot
     ntk=ntsoot
    istat=nf90_inq_varid(id1,'wcmsoot',idwcm)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'wmicsoot',idwave)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'rnsoot',idrn)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'risoot',idri)
     call handle(iout,istat)
    istat=nf90_inq_varid(id1,'tksoot',idtk)
     call handle(iout,istat)
   end if

  if (nopr .eq. 1) then
   write(0,fmt=190) nlines,idwcm,idwave,idrn,idri,idtk
   write(iout,fmt=190) nlines,idwcm,idwave,idrn,idri,idtk
   190 format(/,2x,"nlines, idwcm,idwave,idrn,idri,idtk ",/,&
   2x,1(1x,i4),2x,5(1x,i4))
  end if 

! ****
! Read in the temperatures
   start(1)=1
   count(1)=ntk
   stride(1)=1
  istat=nf90_get_var(id1,idtk,temps,start,count,stride)

   write(0,fmt=200)
   200 format(/,2x,"The temperatures to choose from")
  do i=1,ntk 
   i2=int(temps(i))
   write(0,fmt=210) i2
   210 format(2x,i4)
  end do
   write(0,fmt=220)
   220 format(/,2x,"Specify the temperature (e.g. 300 C )")
   read(5,*) temp

   itemp=-99
  if (temp .lt. temps(1)) then
   itemp=1
   go to 230 
  end if
  if (temp .ge. temps(ntk)) then
   itemp=ntk
   go to 230 
  end if
  do i=1,ntk-1 
    i2=i+1
   if ((temp .ge. temps(i)) .and. (temp .lt. temps(i2))) then
    itemp=i
    go to 230 
   end if
  end do
  230 mn=0

  if (itemp .le. 0) then
   write(iout,fmt=235) itemp
   235 format(/,2x," itemp ",i4," will stop")
   stop
  end if

  write(iout,fmt=240) itemp,temp,temps(itemp)
  240 format(2x,"itemp,temp,temps(itemp)",/,2x,i4,2(1x,f10.4))

! ****
!  Read in the real and imaginary indices

   start(1)=1
   count(1)=nlines
   stride(1)=1
  istat=nf90_get_var(id1,idwcm,wcmdat,start,count,stride)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idwave,wavedat,start,count,stride)
   call handle(iout,istat)

   start2(1)=1
   start2(2)=1
   count2(1)=ndat
   count2(2)=mtempmax
   stride2(1)=1
   stride2(2)=1
  istat=nf90_get_var(id1,idrn,rnval,start2,count2,stride2)
   call handle(iout,istat)
  istat=nf90_get_var(id1,idri,rival,start2,count2,stride2)
   call handle(iout,istat)

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

! Close the data set
  istat=nf90_close(id1)

! *****
! Put values into the output arrays
  do i=1,nlines
   wcm(i)=wcmdat(i)
   wavelength(i)=wavedat(i)
   rndat(i)=rnval(i,itemp)
   ridat(i)=rival(i,itemp)
  end do

! *****
  if (nopr .eq. 1) then
    write(iout,fmt=800) mset,files(mset)
    800 format(/,2x,"rdstagg_carbon: work with mset, files(mset)",/,&
    2x,i3,2x,a16)
    write(iout,fmt=840)
    840 format(2x,"rdstagg_carbon: 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,f10.4),1p,2(1x,e10.3))
   end do
  end if

! ******
  return
  end


