anchor.diff:            call readf(enermax)
             !
          case('GAUSSIAN','GAUSS','DOPPL','DOPPLER','RECT','BOX','BIN','STICKS','STICK','GAUS0','DOPP0',&
  -            'LOREN','LORENTZIAN','LORENTZ','MAX','VOIGT','HITRAN','PSEUDO','PSE-ROCCO','PSE-LIU','VOI-QUAD')
  +            'LOREN','LORENTZIAN','LORENTZ','MAX','VOIGT','PSEUDO','PSE-ROCCO','PSE-LIU','VOI-QUAD')
             !
             proftype = trim(w)
             !
  
positive.code:   !
    integer(ik),parameter   :: nfiles_max =1000, max_items = 1000, nspecies_max = 10, nquadmax = 100, filtermax = 100
    !
    integer(ik)   :: GNS=1,npoints=1001,nchar=1,nfiles=1,ipartf=0,verbose=3,ioffset = 10
    real(rk)      :: temp=298.0,partfunc=0,freql=-small_,freqr= 20000.0,thresh=1.0d-90,halfwidth=1e-2,meanmass=1.0,maxtemp=10000.0
    real(rk)      :: voigt_gamma = 0.05, voigt_n = 0.44, offset = -25.0, pressure = 1.0_rk 
    real(rk)      :: enermax = 1e6, abscoef_thresh = 1.0d-50
    integer(ik)   :: nquad = 20      ! Number of quadrature points 
  
    !
    character(len=cl) :: specttype="absorption",enrfilename,intfilename(nfiles_max),proftype="DOPPL",output="output"
    character(4) a_fmt
    character(9) b_fmt
    !
  
negative.code:      !
     endif
     !
     ! open and read broenning files
     !
     do i =1,Nspecies
       !
       if ( trim(species(i)%filename)/="" ) then
         !
         write(ioname, '(a)') 'brodenning file'
         !
         call IOstart(trim(ioname),bunit)
         open(unit=bunit,file=trim(species(i)%filename),action='read',status='old')
         !
         Jmax = 0 
         !
         do
           !
           ! Scan and find Jmax
           read(bunit,*,end=14) gamma_,n_,J_
           !
           Jmax = max(Jmax,J_)
           !
           cycle 
           14  exit
           !
         enddo
         !
         if (Jmax>0) then
           !
           allocate(species(i)%gammaQN(0:Jmax,-1:1),stat=info)
           call ArrayStart('gammaQN',info,size(species(i)%gammaQN),kind(species(i)%gammaQN))
           species(i)%gammaQN(:,:) = species(i)%gamma
           !
           allocate(species(i)%nQN(0:Jmax,-1:1),stat=info)
           call ArrayStart('nQN',info,size(species(i)%nQN),kind(species(i)%nQN))
           species(i)%nQN(:,:) = species(i)%N
           !
          else
           !
           write(out,"('Jmax = 0 in ',a)") trim(species(i)%filename)
           stop 'Jmax = 0 in a broadening file'
           !
         endif
         !
         rewind(bunit)
         !
         ! read in the QN-dependent broadening parameters
         !
         do
           !
           ! model specific read
           !
           select case ( trim(species(i)%model) )
             !
           case ('J','A0')
             !
             read(bunit,*,end=15) gamma_,n_,J
             species(i)%gammaQN(J,:) = gamma_
             species(i)%nQN(J,:) = n_
             !
           case ('JJ','A1')
             !
             read(bunit,*,end=15) gamma_,n_,Jpp,Jp
             !
             if (abs(Jp-Jpp)>1) stop 'Jp-Jpp in broadening file'
             !
             species(i)%gammaQN(Jpp,Jp-Jpp) = gamma_
             species(i)%nQN(Jpp,Jp-Jpp) = n_
             !
           case ('JJKK-X')
             !
             read(bunit,*,end=15) gamma_,n_,Jpp,Jp,Kp,Kpp
             !
             if (abs(Jp-Jpp)>1) stop 'Jp-Jpp in broadening file'
             if (Kp>Jp.or.Kpp>Jpp) stop 'K > J in broadening file'
             !
             !species(i)%gammaQN(Jpp,Jp-Jpp,Kp,Kpp) = gamma_
             !species(i)%nQN(Jpp,Jp-Jpp,Kp,Kpp) = n_
             !
           case default 
             !
             write(out,"('The broadening model',a,' is not implemented yet')") trim(species(i)%model)
             stop 'Illegal broadening model'
             !
           end select 
           !
           cycle 
           15  exit
           !
         enddo
         !
       endif
       !
       close(bunit)
       ! 
     enddo
     !
     ! Intensities
     !
  