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_isotopic_composition(iso_comp_file,
     &     Nlev,Nmol,mol_name,mol_index,mol_niso,iso_index,
     &     iso_abundance,iso_Qref,iso_g,iso_mass)
      implicit none
      include 'max.inc'
      include 'formats.inc'
c
c     Purpose: to generate an isotopic composition file compatible with kspectrum
c
c     Input:
c       + iso_comp_file: name of the file to generate
c
c     Output:
c       + Nlev: number of atmospheric levels
c       + Nmol: number of molecules
c       + mol_name: description of each molecule
c       + mol_index: index of each molecule
c       + mol_niso: number of isotopes for each molecule
c       + iso_index: index of each isotope, for each molecule
c       + iso_abundance: abundance of each isotope, for each molecule
c       + iso_Qref: reference partition function of each isotope, for each molecule
c       + iso_g: parameter g of each isotope, for each molecule
c       + iso_mass: molar mass of each isotope, for each molecule (g/mol)
c
c     I/O
      character*(Nchar_mx) iso_comp_file
      integer Nlev
      integer Nmol
      character*10 mol_name(1:Nmol_max)
      integer mol_index(1:Nmol_max)
      integer mol_niso(1:Nmol_max)
      integer iso_index(1:Nmol_max,1:Niso_max)
      double precision iso_abundance(1:Nmax,1:Nmol_max,1:Niso_max)
      double precision iso_Qref(1:Nmol_max,1:Niso_max)
      integer iso_g(1:Nmol_max,1:Niso_max)
      double precision iso_mass(1:Nmol_max,1:Niso_max)
c     temp
      integer i,imol,iso,ilev,ios,ncm,nl
      character*(Nchar_mx) string
      character*(Nchar_mx) line,nlines_file
      integer line_index(1:Nmol_max)
c     label
      integer strlen
      character*(Nchar_mx) label
      label='subroutine read_isotopic_composition'

      open(12,file=iso_comp_file(1:strlen(iso_comp_file)),
     &     status='old',iostat=ios)
      if (ios.ne.0) then        ! file not found
         write(*,*) 'Error from routine read_composition:'
         write(*,*) 'Data file could not be found:'
         write(*,*) iso_comp_file(1:strlen(iso_comp_file))
         stop
      else
         write(*,*) 'Reading isotopic composition from file:'
         write(*,*) iso_comp_file(1:strlen(iso_comp_file))
      endif

      nlines_file='./nlines'
      call get_nlines(iso_comp_file,nlines_file,nl)

      read(12,*)
      read(12,*) Nlev
      read(12,*)
      Nmol=0
      ncm=50
      do i=4,nl
         read(12,10) line
         if (strlen(line).lt.ncm) then
            Nmol=Nmol+1
            line_index(Nmol)=i
         endif
      enddo ! i
      rewind(12)

      imol=0
      do i=1,3
         read(12,*)
      enddo
      do i=4,nl
         if (i.eq.line_index(imol+1)) then
            read(12,10) line
            imol=imol+1
            call get_mnami(line,mol_name(imol),mol_index(imol))
            if (imol-1.ge.1) then
               mol_niso(imol-1)=iso
            endif
            iso=0
         else
            iso=iso+1
            read(12,64) iso_index(imol,iso)
     &           ,iso_Qref(imol,iso),iso_g(imol,iso),iso_mass(imol,iso)
     &           ,(iso_abundance(ilev,imol,iso),ilev=1,Nlev)
         endif
         if (i.eq.nl) then
            mol_niso(Nmol)=iso
         endif
      enddo ! i

      close(12)

      return
      end