X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Fwham%2Fsrc-NEWSC-NEWCORR%2Freadrtns.F.org;fp=source%2Fwham%2Fsrc-NEWSC-NEWCORR%2Freadrtns.F.org;h=1fa6e46397ccf354b44c9286c9bbd0dedd2c3e2c;hb=d101c97dea752458d76055fdbae49c26fff03c1f;hp=0000000000000000000000000000000000000000;hpb=325eda160c9ad2982501e091ca40606a29043712;p=unres.git

diff --git a/source/wham/src-NEWSC-NEWCORR/readrtns.F.org b/source/wham/src-NEWSC-NEWCORR/readrtns.F.org
new file mode 100644
index 0000000..1fa6e46
--- /dev/null
+++ b/source/wham/src-NEWSC-NEWCORR/readrtns.F.org
@@ -0,0 +1,691 @@
+      subroutine read_general_data(*)
+      implicit none
+      include "DIMENSIONS"
+      include "DIMENSIONS.ZSCOPT"
+      include "DIMENSIONS.FREE"
+      include "COMMON.TORSION"
+      include "COMMON.INTERACT"
+      include "COMMON.IOUNITS"
+      include "COMMON.TIME1"
+      include "COMMON.PROT"
+      include "COMMON.PROTFILES"
+      include "COMMON.CHAIN"
+      include "COMMON.NAMES"
+      include "COMMON.FFIELD"
+      include "COMMON.ENEPS"
+      include "COMMON.WEIGHTS"
+      include "COMMON.FREE"
+      include "COMMON.CONTROL"
+      include "COMMON.ENERGIES"
+      character*800 controlcard
+      integer i,j,k,ii,n_ene_found
+      integer ind,itype1,itype2,itypf,itypsc,itypp
+      integer ilen
+      external ilen
+      character*16 ucase
+      character*16 key
+      external ucase
+
+      call card_concat(controlcard,.true.)
+      call readi(controlcard,"N_ENE",n_ene,max_ene)
+      if (n_ene.gt.max_ene) then
+        write (iout,*) "Error: parameter out of range: N_ENE",n_ene,
+     &    max_ene
+        return1
+      endif
+      call readi(controlcard,"NPARMSET",nparmset,1)
+      if (nparmset.gt.max_parm) then
+        write (iout,*) "Error: parameter out of range: NPARMSET",
+     &    nparmset, Max_Parm
+        return1
+      endif
+      call readi(controlcard,"MAXIT",maxit,5000)
+      call reada(controlcard,"FIMIN",fimin,1.0d-3)
+      call readi(controlcard,"ENSEMBLES",ensembles,0)
+      write (iout,*) "Number of energy parameter sets",nparmset
+      call multreadi(controlcard,"ISAMPL",isampl,nparmset,1)
+      write (iout,*) "MaxSlice",MaxSlice
+      call readi(controlcard,"NSLICE",nslice,1)
+      call flush(iout)
+      if (nslice.gt.MaxSlice) then
+        write (iout,*) "Error: parameter out of range: NSLICE",nslice,
+     &    MaxSlice
+        return1
+      endif
+      write (iout,*) "Frequency of storing conformations",
+     & (isampl(i),i=1,nparmset)
+      write (iout,*) "Maxit",maxit," Fimin",fimin
+      call readi(controlcard,"NQ",nQ,1)
+      if (nQ.gt.MaxQ) then
+        write (iout,*) "Error: parameter out of range: NQ",nq,
+     &    maxq
+        return1
+      endif
+      indpdb=0
+      if (index(controlcard,"CLASSIFY").gt.0) indpdb=1
+      call reada(controlcard,"DELTA",delta,1.0d-2)
+      call readi(controlcard,"EINICHECK",einicheck,2)
+      call reada(controlcard,"DELTRMS",deltrms,5.0d-2)
+      call reada(controlcard,"DELTRGY",deltrgy,5.0d-2)
+      call readi(controlcard,"RESCALE",rescale_mode,1)
+      write (iout,*) "delta",delta
+      write (iout,*) "einicheck",einicheck
+      write (iout,*) "rescale_mode",rescale_mode
+      call flush(iout)
+      bxfile=index(controlcard,"BXFILE").gt.0
+      cxfile=index(controlcard,"CXFILE").gt.0
+      if (nslice .eq. 1 .and. .not.bxfile .and. .not.cxfile)
+     & bxfile=.true.
+      histfile=index(controlcard,"HISTFILE").gt.0
+      entfile=index(controlcard,"ENTFILE").gt.0
+      zscfile=index(controlcard,"ZSCFILE").gt.0
+      return
+      end
+c------------------------------------------------------------------------------
+      subroutine read_efree(iparm,*)
+C
+C Read molecular data
+C
+      implicit none
+      include 'DIMENSIONS'
+      include 'DIMENSIONS.ZSCOPT'
+      include 'DIMENSIONS.COMPAR'
+      include 'DIMENSIONS.FREE'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.TIME1'
+      include 'COMMON.SBRIDGE'
+      include 'COMMON.CONTROL'
+      include 'COMMON.CHAIN'
+      include 'COMMON.HEADER'
+      include 'COMMON.GEO'
+      include 'COMMON.FREE'
+      character*320 controlcard,ucase
+      integer iparm,ib,i,j
+      integer ilen
+      external ilen
+      call card_concat(controlcard,.true.)
+      call readi(controlcard,'NT',nT_h(iparm),1)
+      if (nT_h(iparm).gt.MaxT_h) then
+        write (iout,*)  "Error: parameter out of range: NT",nT_h(iparm),
+     &    MaxT_h
+        return1
+      endif
+      replica(iparm)=index(controlcard,"REPLICA").gt.0
+      umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
+      read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
+      write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
+     &  replica(iparm)," umbrella ",umbrella(iparm),
+     &  " read_iset",read_iset(iparm)
+      call flush(iout)
+      do ib=1,nT_h(iparm)
+        call card_concat(controlcard,.true.)
+        call readi(controlcard,'NR',nR(ib,iparm),1)
+        if (umbrella(iparm)) then
+          nRR(ib,iparm)=1
+        else
+          nRR(ib,iparm)=nR(ib,iparm)
+        endif
+        if (nR(ib,iparm).gt.MaxR) then
+          write (iout,*)  "Error: parameter out of range: NR",
+     &      nR(ib,iparm),MaxR
+        return1
+        endif
+        call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
+        beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
+        call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
+     &    0.0d0)
+        do i=1,nR(ib,iparm)
+          call card_concat(controlcard,.true.)
+          call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
+     &      100.0d0)
+          call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
+     &      0.0d0)
+        enddo
+      enddo
+      do ib=1,nT_h(iparm)
+        write (iout,*) "ib",ib," beta_h",
+     &    1.0d0/(0.001987*beta_h(ib,iparm))
+        write (iout,*) "nR",nR(ib,iparm)
+        write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
+        do i=1,nR(ib,iparm)
+          write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
+     &      "q0",(q0(j,i,ib,iparm),j=1,nQ)
+        enddo
+        call flush(iout)
+      enddo
+      return
+      end
+c-----------------------------------------------------------------------------
+      subroutine read_protein_data(iparm,*)
+      implicit none
+      include "DIMENSIONS"
+      include "DIMENSIONS.ZSCOPT"
+      include "DIMENSIONS.FREE"
+#ifdef MPI
+      include "mpif.h"
+      integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+      include "COMMON.MPI"
+#endif
+      include "COMMON.CHAIN"
+      include "COMMON.IOUNITS"
+      include "COMMON.PROT"
+      include "COMMON.PROTFILES"
+      include "COMMON.NAMES"
+      include "COMMON.FREE"
+      include "COMMON.OBCINKA"
+      character*64 nazwa
+      character*16000 controlcard
+      integer i,ii,ib,iR,iparm,ilen,iroof,nthr
+      external ilen,iroof
+      call flush(iout)
+C Read names of files with conformation data.
+      if (replica(iparm)) then
+        nthr = 1
+      else
+        nthr = nT_h(iparm)
+      endif
+      do ib=1,nthr
+      do ii=1,nRR(ib,iparm)
+      write (iout,*) "Parameter set",iparm," temperature",ib,
+     & " window",ii
+      call card_concat(controlcard,.true.) 
+      write (iout,*) controlcard(:ilen(controlcard))
+      call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
+      call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
+      call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
+      call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
+      call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
+     & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
+      call reada(controlcard,"TIME_START",
+     &  time_start_collect(ii,ib,iparm),0.0d0)
+      call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
+     &  1.0d10)
+      write (iout,*) "rec_start",rec_start(ii,ib,iparm),
+     & " rec_end",rec_end(ii,ib,iparm)
+      write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
+     & " time_end",time_end_collect(ii,ib,iparm)
+      call flush(iout)
+      if (replica(iparm)) then
+        call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
+        write (iout,*) "Number of trajectories",totraj(ii,iparm)
+        call flush(iout)
+      endif
+      if (nfile_bin(ii,ib,iparm).lt.2 
+     &    .and. nfile_asc(ii,ib,iparm).eq.0
+     &    .and. nfile_cx(ii,ib,iparm).eq.0) then
+        write (iout,*) "Error - no action specified!"
+        return1
+      endif
+      if (nfile_bin(ii,ib,iparm).gt.0) then
+        call card_concat(controlcard,.false.)
+        call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
+     &   maxfile_prot,nfile_bin(ii,ib,iparm))
+#ifdef DEBUG
+        write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
+        write(iout,*) (protfiles(i,1,ii,ib,iparm),
+     &    i=1,nfile_bin(ii,ib,iparm))
+#endif
+      endif
+      if (nfile_asc(ii,ib,iparm).gt.0) then
+        call card_concat(controlcard,.false.)
+        call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
+     &   maxfile_prot,nfile_asc(ii,ib,iparm))
+#ifdef DEBUG
+        write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
+        write(iout,*) (protfiles(i,2,ii,ib,iparm),
+     &    i=1,nfile_asc(ii,ib,iparm))
+#endif
+      else if (nfile_cx(ii,ib,iparm).gt.0) then
+        call card_concat(controlcard,.false.)
+        call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
+     &   maxfile_prot,nfile_cx(ii,ib,iparm))
+#ifdef DEBUG
+        write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
+        write(iout,*) (protfiles(i,2,ii,ib,iparm),
+     &    i=1,nfile_cx(ii,ib,iparm))
+#endif
+      endif
+      call flush(iout)
+      enddo
+      enddo
+      return
+      end
+c-------------------------------------------------------------------------------
+      subroutine opentmp(islice,iunit,bprotfile_temp)
+      implicit none
+      include "DIMENSIONS"
+      include "DIMENSIONS.ZSCOPT"
+      include "DIMENSIONS.FREE"
+#ifdef MPI
+      include "mpif.h"
+      integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+      include "COMMON.MPI"
+#endif
+      include "COMMON.IOUNITS"
+      include "COMMON.PROTFILES"
+      include "COMMON.PROT"
+      character*64 bprotfile_temp
+      character*3 liczba
+      character*2 liczba1
+      integer iunit,islice
+      integer ilen,iroof
+      external ilen,iroof
+      logical lerr
+
+      write (liczba1,'(bz,i2.2)') islice
+#ifdef MPI
+      write (liczba,'(bz,i3.3)') me
+      bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
+     &  prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
+      open (iunit,file=bprotfile_temp,status="unknown",
+     &    form="unformatted",access="direct",recl=lenrec)
+#else
+      bprotfile_temp = scratchdir(:ilen(scratchdir))//
+     &  "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
+      open (iunit,file=bprotfile_temp,status="unknown",
+     &    form="unformatted",access="direct",recl=lenrec)
+#endif      
+      return
+      end
+c-------------------------------------------------------------------------------
+      subroutine read_database(*)
+      implicit none
+      include "DIMENSIONS"
+      include "DIMENSIONS.ZSCOPT"
+      include "DIMENSIONS.FREE"
+#ifdef MPI
+      include "mpif.h"
+      integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+      include "COMMON.MPI"
+#endif
+      include "COMMON.CHAIN"
+      include "COMMON.IOUNITS"
+      include "COMMON.PROTFILES"
+      include "COMMON.NAMES"
+      include "COMMON.VAR"
+      include "COMMON.GEO"
+      include "COMMON.ENEPS"
+      include "COMMON.PROT"
+      include "COMMON.INTERACT"
+      include "COMMON.FREE"
+      include "COMMON.SBRIDGE"
+      include "COMMON.OBCINKA"
+      real*4 csingle(3,maxres2)
+      character*64 nazwa,bprotfile_temp
+      character*3 liczba
+      character*2 liczba1
+      integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
+     &  ll(maxslice),mm(maxslice),if
+      integer nrec,nlines,iscor,iunit,islice
+      double precision energ
+      integer ilen,iroof
+      external ilen,iroof
+      double precision rmsdev,energia(0:max_ene),efree,eini,temp
+      double precision prop(maxQ)
+      integer ntot_all(maxslice,0:maxprocs-1)
+      integer iparm,ib,iib,ir,nprop,nthr
+      double precision etot,time
+      integer ixdrf,iret 
+      logical lerr,linit
+
+      lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
+      lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
+      lenrec=lenrec2+8
+      write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
+     &  " lenrec2",lenrec2
+
+      do i=1,nQ
+        prop(i)=0.0d0
+      enddo
+      do islice=1,nslice
+        ll(islice)=0
+        mm(islice)=0
+      enddo
+      write (iout,*) "nparmset",nparmset
+      do iparm=1,nparmset
+
+      if (replica(iparm)) then
+        nthr = 1
+      else
+        nthr = nT_h(iparm)
+      endif
+
+      do ib=1,nthr
+      do iR=1,nRR(ib,iparm)
+
+      write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
+      do islice=1,nslice
+        jj(islice)=0
+        kk(islice)=0
+      enddo
+
+      IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
+c Read conformations from binary DA files (one per batch) and write them to 
+c a binary DA scratchfile.
+        write (liczba,'(bz,i3.3)') me
+        do if=1,nfile_bin(iR,ib,iparm)
+          nazwa=protfiles(if,1,iR,ib,iparm)
+     &     (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
+          open (ientin,file=nazwa,status="old",form="unformatted",
+     &     access="direct",recl=lenrec2,err=1111)
+          ii=0
+          do islice=1,nslice
+            call opentmp(islice,ientout,bprotfile_temp)
+            call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
+     &        mm(islice),iR,ib,iparm)
+            close(ientout)
+          enddo
+          close(ientin)
+        enddo
+      ENDIF ! NFILE_BIN>0
+c
+      IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
+c Read conformations from multiple ASCII int files and write them to a binary
+c DA scratchfile.
+        do if=1,nfile_asc(iR,ib,iparm)
+          nazwa=protfiles(if,2,iR,ib,iparm)
+     &     (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
+          open(unit=ientin,file=nazwa,status='old',err=1111)
+          write(iout,*) "reading ",nazwa(:ilen(nazwa))
+          ii=0
+          call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
+        enddo ! if
+      ENDIF
+      IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
+c Read conformations from cx files and write them to a binary
+c DA scratchfile.
+        do if=1,nfile_cx(iR,ib,iparm)
+          nazwa=protfiles(if,2,iR,ib,iparm)
+     &     (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
+          write(iout,*) "reading ",nazwa(:ilen(nazwa))
+          ii=0
+          print *,"Calling cxread"
+          call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
+     &       *1111)
+          close(ientout)
+          write (iout,*) "exit cxread"
+          call flush(iout)
+        enddo
+      ENDIF
+
+      do islice=1,nslice
+        stot(islice)=stot(islice)+jj(islice)
+      enddo
+
+      enddo
+      enddo
+      write (iout,*) "IPARM",iparm
+      enddo
+
+      if (nslice.eq.1) then
+#ifdef MPI
+        write (liczba,'(bz,i3.3)') me
+        bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
+     &    prefix(:ilen(prefix))//liczba//".xbin.tmp"
+#else
+        bprotfile_temp = scratchdir(:ilen(scratchdir))//
+     &     "/"//prefix(:ilen(prefix))//".xbin.tmp"
+#endif
+        write(iout,*) mm(1)," conformations read",ll(1),
+     &    " conformations written to ", 
+     &    bprotfile_temp(:ilen(bprotfile_temp))
+      else
+        do islice=1,nslice
+          write (liczba1,'(bz,i2.2)') islice
+#ifdef MPI
+          write (liczba,'(bz,i3.3)') me
+          bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
+     &      prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
+#else
+          bprotfile_temp = scratchdir(:ilen(scratchdir))//
+     &       "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
+#endif
+          write(iout,*) mm(islice)," conformations read",ll(islice),
+     &    " conformations written to ", 
+     &    bprotfile_temp(:ilen(bprotfile_temp))
+        enddo
+      endif
+
+#ifdef MPI
+c Check if everyone has the same number of conformations
+      call MPI_Allgather(stot(1),maxslice,MPI_INTEGER,
+     &  ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
+      lerr=.false.
+      do i=0,nprocs-1
+        if (i.ne.me) then
+          do islice=1,nslice
+          if (stot(islice).ne.ntot_all(islice,i)) then
+            write (iout,*) "Number of conformations at processor",i,
+     &       " differs from that at processor",me,
+     &       stot(islice),ntot_all(islice,i)," slice",islice
+            lerr = .true.
+          endif
+          enddo
+        endif
+      enddo 
+      if (lerr) then
+        write (iout,*)
+        write (iout,*) "Numbers of conformations read by processors"
+        write (iout,*)
+        do i=0,nprocs-1
+          write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
+        enddo
+        write (iout,*) "Calculation terminated."
+        call flush(iout)
+        return1
+      endif
+      do islice=1,nslice
+        ntot(islice)=stot(islice)
+      enddo
+      return
+#endif
+ 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
+      call flush(iout)
+      return1
+      end
+c------------------------------------------------------------------------------
+      subroutine card_concat(card,to_upper)
+      implicit none
+      include 'DIMENSIONS.ZSCOPT'
+      include "COMMON.IOUNITS"
+      character*(*) card
+      character*80 karta,ucase
+      logical to_upper
+      integer ilen
+      external ilen
+      read (inp,'(a)') karta
+      if (to_upper) karta=ucase(karta)
+      card=' '
+      do while (karta(80:80).eq.'&')
+        card=card(:ilen(card)+1)//karta(:79)
+        read (inp,'(a)') karta
+        if (to_upper) karta=ucase(karta)
+      enddo
+      card=card(:ilen(card)+1)//karta
+      return
+      end
+c------------------------------------------------------------------------------
+      subroutine readi(rekord,lancuch,wartosc,default)
+      implicit none
+      character*(*) rekord,lancuch
+      integer wartosc,default
+      integer ilen,iread
+      external ilen
+      iread=index(rekord,lancuch(:ilen(lancuch))//"=")
+      if (iread.eq.0) then
+        wartosc=default
+        return
+      endif
+      iread=iread+ilen(lancuch)+1
+      read (rekord(iread:),*) wartosc
+      return
+      end
+c----------------------------------------------------------------------------
+      subroutine reada(rekord,lancuch,wartosc,default)
+      implicit none
+      character*(*) rekord,lancuch
+      character*80 aux
+      double precision wartosc,default
+      integer ilen,iread
+      external ilen
+      iread=index(rekord,lancuch(:ilen(lancuch))//"=")
+      if (iread.eq.0) then
+        wartosc=default
+        return
+      endif
+      iread=iread+ilen(lancuch)+1
+      read (rekord(iread:),*) wartosc
+      return
+      end
+c----------------------------------------------------------------------------
+      subroutine multreadi(rekord,lancuch,tablica,dim,default)
+      implicit none
+      integer dim,i
+      integer tablica(dim),default
+      character*(*) rekord,lancuch
+      character*80 aux
+      integer ilen,iread
+      external ilen
+      do i=1,dim
+        tablica(i)=default
+      enddo
+      iread=index(rekord,lancuch(:ilen(lancuch))//"=")
+      if (iread.eq.0) return
+      iread=iread+ilen(lancuch)+1
+      read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
+   10 return
+      end
+c----------------------------------------------------------------------------
+      subroutine multreada(rekord,lancuch,tablica,dim,default)
+      implicit none
+      integer dim,i
+      double precision tablica(dim),default
+      character*(*) rekord,lancuch
+      character*80 aux
+      integer ilen,iread
+      external ilen
+      do i=1,dim
+        tablica(i)=default
+      enddo
+      iread=index(rekord,lancuch(:ilen(lancuch))//"=")
+      if (iread.eq.0) return
+      iread=iread+ilen(lancuch)+1
+      read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
+   10 return
+      end
+c----------------------------------------------------------------------------
+      subroutine reads(rekord,lancuch,wartosc,default)
+      implicit none
+      character*(*) rekord,lancuch,wartosc,default
+      character*80 aux
+      integer ilen,lenlan,lenrec,iread,ireade
+      external ilen
+      logical iblnk
+      external iblnk
+      lenlan=ilen(lancuch)
+      lenrec=ilen(rekord)
+      iread=index(rekord,lancuch(:lenlan)//"=")
+c      print *,"rekord",rekord," lancuch",lancuch
+c      print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
+      if (iread.eq.0) then
+        wartosc=default
+        return
+      endif
+      iread=iread+lenlan+1
+c      print *,"iread",iread
+c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
+      do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
+        iread=iread+1
+c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
+      enddo
+c      print *,"iread",iread
+      if (iread.gt.lenrec) then
+         wartosc=default
+        return
+      endif
+      ireade=iread+1
+c      print *,"ireade",ireade
+      do while (ireade.lt.lenrec .and.
+     &   .not.iblnk(rekord(ireade:ireade)))
+        ireade=ireade+1
+      enddo
+      wartosc=rekord(iread:ireade)
+      return
+      end
+c----------------------------------------------------------------------------
+      subroutine multreads(rekord,lancuch,tablica,dim,default)
+      implicit none
+      integer dim,i
+      character*(*) rekord,lancuch,tablica(dim),default
+      character*80 aux
+      integer ilen,lenlan,lenrec,iread,ireade
+      external ilen
+      logical iblnk
+      external iblnk
+      do i=1,dim
+        tablica(i)=default
+      enddo
+      lenlan=ilen(lancuch)
+      lenrec=ilen(rekord)
+      iread=index(rekord,lancuch(:lenlan)//"=")
+c      print *,"rekord",rekord," lancuch",lancuch
+c      print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
+      if (iread.eq.0) return
+      iread=iread+lenlan+1
+      do i=1,dim
+c      print *,"iread",iread
+c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
+      do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
+        iread=iread+1
+c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
+      enddo
+c      print *,"iread",iread
+      if (iread.gt.lenrec) return
+      ireade=iread+1
+c      print *,"ireade",ireade
+      do while (ireade.lt.lenrec .and.
+     &   .not.iblnk(rekord(ireade:ireade)))
+        ireade=ireade+1
+      enddo
+      tablica(i)=rekord(iread:ireade)
+      iread=ireade+1
+      enddo
+      end
+c----------------------------------------------------------------------------
+      subroutine split_string(rekord,tablica,dim,nsub)
+      implicit none
+      integer dim,nsub,i,ii,ll,kk
+      character*(*) tablica(dim)
+      character*(*) rekord
+      integer ilen
+      external ilen
+      do i=1,dim
+        tablica(i)=" "
+      enddo
+      ii=1
+      ll = ilen(rekord)
+      nsub=0
+      do i=1,dim
+C Find the start of term name
+        kk = 0
+        do while (ii.le.ll .and. rekord(ii:ii).eq." ") 
+          ii = ii+1
+        enddo
+C Parse the name into TABLICA(i) until blank found
+        do while (ii.le.ll .and. rekord(ii:ii).ne." ") 
+          kk = kk+1
+          tablica(i)(kk:kk)=rekord(ii:ii)
+          ii = ii+1
+        enddo
+        if (kk.gt.0) nsub=nsub+1
+        if (ii.gt.ll) return
+      enddo
+      return
+      end
+c--------------------------------------------------------------------------------
+      integer function iroof(n,m)
+      ii = n/m
+      if (ii*m .lt. n) ii=ii+1
+      iroof = ii
+      return
+      end