c     kspectrum (http://www.meso-star.com/en_Products.html) - This file is part of kspectrum
c     Copyright (C) 2008-2015 - Méso-Star - Vincent Eymet
c     
c     This file must be used under the terms of the CeCILL license.
c     This source file is licensed as described in the file COPYING, which
c     you should have received as part of this distribution.  The terms
c     are also available at
c     http://www.cecill.info/licences/Licence_CeCILL_V2.1-en.txt
c     
      subroutine read_sgrid(nu_min,nu_max,Nks,sgrid)
      implicit none
      include 'max.inc'
c
c     Purpose: to read user-defined spectral grid 
c
c     Input:
c       + nu_min, nu_max: limits of the band to read
c
c     Output:
c       + Nks: number of 'nu' values in the [nu_min,numax] interval
c       + sgrid: user-defined grid in the [nu_min,numax] interval
c
c     I/O
      double precision nu_min,nu_max
      integer Nks
      double precision sgrid(1:Nkmx)
c     temp
      integer i,i1,i2
      logical nu1_found,nu2_found
      integer ios
      character*(Nchar_mx) sgrid_file
      character*(Nchar_mx) stringp_file
      integer Nk_total
      double precision tmp1
c     label
      integer strlen
      character*(Nchar_mx) label
      label='subroutine read_sgrid'

      sgrid_file='./data/user_spectral_grid.in'

      stringp_file='./nlines0'
      call get_nlines(sgrid_file,stringp_file,Nk_total)

      open(10,file=sgrid_file(1:strlen(sgrid_file))
     &     ,status='old',iostat=ios)
      if (ios.ne.0) then        ! file not found
         write(*,*) 'Error from ',
     &        label(1:strlen(label)),' :'
         write(*,*) 'File could not be opened:',
     &        sgrid_file(1:strlen(sgrid_file))
         stop
      else
         nu1_found=.false.
         i=0
         do while (.not.nu1_found)
            read(10,*) tmp1
            i=i+1
            if (tmp1.ge.nu_min) then
               nu1_found=.true.
               i1=i
            endif
            if ((.not.nu1_found).and.(i.ge.Nk_total)) then
               call error(label)
               write(*,*) 'nu_min=',nu_min
               write(*,*) 'could not be found in file:'
               write(*,*) sgrid_file(1:strlen(sgrid_file))
               stop
            endif
         enddo ! nu1_found
         Nks=0
         if (nu_min.lt.tmp1) then
            Nks=Nks+1
            sgrid(Nks)=nu_min
         endif
         Nks=Nks+1
         sgrid(Nks)=tmp1
         nu2_found=.false.
         do while (.not.nu2_found)
            read(10,*) tmp1
            i=i+1
            Nks=Nks+1
            if (Nks.gt.Nkmx) then
               call error(label)
               write(*,*) 'Nks=',Nks
               write(*,*) 'while Nkmx=',Nkmx
               stop
            else
               sgrid(Nks)=tmp1
            endif
            if (tmp1.ge.nu_max) then
               nu2_found=.true.
               i2=i
            endif
            if ((.not.nu2_found).and.(i.ge.Nk_total)) then
               call error(label)
               write(*,*) 'nu_max=',nu_max
               write(*,*) 'could not be found in file:'
               write(*,*) sgrid_file(1:strlen(sgrid_file))
               stop
            endif
         enddo ! nu2_found
         if (nu_max.gt.sgrid(Nks)) then
            Nks=Nks+1
            sgrid(Nks)=nu_max
         endif
      endif ! ios
      close(10)

      return
      end