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 parse_file(file_idx,nk,slimits1) implicit none c Purpose: to locate line indexes in the given LBL data file where each 1 1/cm spectral interval begins and ends include 'max.inc' include 'formats.inc' character*(Nchar_mx) datafile,results_file integer nl integer ik,nk,j double precision slimits1(1:Nti_mx,2) integer slimits2(1:Nti_mx,2) integer index_tmp,isotope_tmp double precision nu0_tmp integer lminf(1:Nti_mx) integer lmaxf(1:Nti_mx) integer line,af,nfile,file_idx,ifile integer nlines(1:nfiles_mx) character*(Nchar_mx) filename(1:nfiles_mx) integer fformat(1:nfiles_mx) integer ios character*(Nchar_mx) lbl_file c label integer strlen character*(Nchar_mx) label label='subroutine parse_file' lbl_file='./optimizations/LBL_files.txt' open(12,file=lbl_file(1:strlen(lbl_file)) & ,status='old',iostat=ios) if (ios.ne.0) then ! list_file not found call error(label) write(*,*) 'File could not be found:' write(*,*) lbl_file(1:strlen(lbl_file)) stop else read(12,*) nfile do ifile=1,nfile read(12,37) fformat(ifile),nlines(ifile),filename(ifile) enddo close(12) endif ! ios datafile=filename(file_idx) nl=nlines(file_idx) do ik=1,nk lminf(ik)=0 lmaxf(ik)=0 enddo c Debug c write(*,*) datafile(1:strlen(datafile)) c Debug open(14,file=datafile(1:strlen(datafile))) do line=1,nl read(14,80) index_tmp ! molecule number [unitless] & ,isotope_tmp ! isotopologue number [unitless] & ,nu0_tmp ! vacuum wavenumber [cm¯¹] do ik=1,nk if ((nu0_tmp.ge.slimits1(ik,1)) & .and.(nu0_tmp.lt.slimits1(ik,2)) & .and.(lminf(ik).eq.0)) then slimits2(ik,1)=line lminf(ik)=1 endif if ((nu0_tmp.ge.slimits1(ik,2)).and.(lminf(ik).eq.1) & .and.(lmaxf(ik).eq.0)) & then slimits2(ik,2)=line-1 lmaxf(ik)=1 endif enddo ! ik if (line.eq.1) then do ik=1,nk if ((lminf(ik).eq.0).and.(lmaxf(ik).eq.0) & .and.(nu0_tmp.gt.slimits1(ik,2))) then slimits2(ik,1)=-1 slimits2(ik,2)=-1 endif enddo ! ik endif if (line.eq.nl) then do ik=1,nk if ((lminf(ik).eq.1).and.(lmaxf(ik).eq.0) & .and.(nu0_tmp.lt.slimits1(ik,2))) then slimits2(ik,2)=nl lmaxf(ik)=1 endif if ((lminf(ik).eq.0).and.(lmaxf(ik).eq.0) & .and.(nu0_tmp.lt.slimits1(ik,1))) then slimits2(ik,1)=-1 slimits2(ik,2)=-1 endif enddo endif af=1 do ik=1,nk if ((lminf(ik).eq.0).or.(lmaxf(ik).eq.0)) then af=0 goto 111 endif enddo ! ik 111 continue if (af.eq.1) then goto 222 endif enddo ! line 222 continue close(14) c Intervals in which no lines have been found, even if their spectral limits c are within the LBL database do ik=1,nk if ((lminf(ik).eq.0).and.(lmaxf(ik).eq.0)) then slimits2(ik,1)=-1 slimits2(ik,2)=-1 endif enddo call slimits_filename(file_idx,results_file) open(15,file=results_file(1:strlen(results_file))) write(15,*) nk do ik=1,nk write(15,*) (slimits1(ik,j),j=1,2),(slimits2(ik,j),j=1,2) enddo close(15) return end subroutine parse_file_ori(datafile,nl,nk,slimits1,results_file) implicit none include 'max.inc' include 'formats.inc' c Purpose: to locate line indexes in the given LBL data file where each 1 1/cm spectral interval begins and ends character*(Nchar_mx) datafile,results_file integer strlen,nl integer ik,nk,j double precision slimits1(1:Nti_mx,2) integer slimits2(1:Nti_mx,2) double precision numin,numax integer index_tmp,isotope_tmp double precision nu0_tmp integer lminf,lmaxf,line,fline do ik=1,nk numin=slimits1(ik,1) numax=slimits1(ik,2) c Debug write(*,*) numin,numax,ik,nk c Debug open(14,file=datafile(1:strlen(datafile))) if (ik.gt.2) then do line=1,slimits2(ik-2,2) read(14,*) enddo fline=slimits2(ik-2,2)+1 else fline=1 endif do line=fline,nl read(14,80) index_tmp ! molecule number [unitless] & ,isotope_tmp ! isotopologue number [unitless] & ,nu0_tmp ! vacuum wavenumber [cm¯¹] if ((nu0_tmp.ge.numin).and.(lminf.eq.0)) then slimits2(ik,1)=line lminf=1 endif if ((nu0_tmp.ge.numax).and.(lminf.eq.1).and.(lmaxf.eq.0)) & then slimits2(ik,2)=line lmaxf=1 endif if ((line.eq.nl).and.(nu0_tmp.lt.numax) & .and.(lminf.eq.1).and.(lmaxf.eq.0)) then slimits2(ik,2)=nl lmaxf=1 endif if ((lminf.eq.1).and.(lmaxf.eq.1)) then goto 222 endif if ((line.eq.1).and.(lminf.eq.0).and.(lmaxf.eq.0) & .and.(nu0_tmp.gt.numax)) then slimits2(ik,1)=-1 slimits2(ik,2)=-1 goto 222 endif if ((line.eq.nl).and.(lminf.eq.0).and.(lmaxf.eq.0) & .and.(nu0_tmp.lt.numin)) then slimits2(ik,1)=-1 slimits2(ik,2)=-1 goto 222 endif enddo ! line 222 continue close(14) enddo ! ik open(15,file=results_file(1:strlen(results_file))) write(15,*) nk do ik=1,nk write(15,*) (slimits1(ik,j),j=1,2),(slimits2(ik,j),j=1,2) enddo close(15) return end subroutine calc_slimits1(Nb,nu_low,nu_hi,large_db,nk,slimits1) implicit none include 'max.inc' integer Nb,nk,band logical large_db double precision nu_low(1:Nbmx) double precision nu_hi(1:Nbmx) double precision slimits1(1:Nti_mx,2) double precision numin,numax,intw c Debug integer ik,j c Debug c label integer strlen character*(Nchar_mx) label label='subroutine calc_slimits1' c Debug c write(*,*) 'large_db=',large_db c Debug nk=0 do band=1,Nb c Debug c write(*,*) 'band=',band,' /',Nb,' [',nu_low(band), c & '-',nu_hi(band),']' c Debug numin=nu_low(band) c Debug c write(*,*) 'numin=',numin,' nu_hi(',band,')=',nu_hi(band) c Debug call choose_intw(nu_low(band),nu_hi(band),large_db,intw) numax=numin+intw c Debug c write(*,*) 'intw=',intw,' numax=',numax c Debug nk=nk+1 slimits1(nk,1)=numin slimits1(nk,2)=numax do while (numax.lt.nu_hi(band)) nk=nk+1 c Tests if (nk.gt.Nti_mx) then call error(label) write(*,*) 'Nti_mx has been exceeded' stop endif c Tests numin=numax call choose_intw(nu_low(band),nu_hi(band),large_db,intw) numax=numin+intw slimits1(nk,1)=numin slimits1(nk,2)=numax c Debug c write(*,*) 'slimits1(',nk,')=',(slimits1(nk,j),j=1,2) c Debug enddo enddo c Debug c do ik=1,nk c write(*,*) 'slimits1(',ik,')=',(slimits1(ik,j),j=1,2) c enddo c Debug return end