Adam's changes
[unres.git] / source / wham / src-M-Pawel / readrtns.F
1       subroutine read_general_data(*)
2       implicit none
3       include "DIMENSIONS"
4       include "DIMENSIONS.ZSCOPT"
5       include "DIMENSIONS.FREE"
6       include "COMMON.TORSION"
7       include "COMMON.INTERACT"
8       include "COMMON.IOUNITS"
9       include "COMMON.TIME1"
10       include "COMMON.PROT"
11       include "COMMON.PROTFILES"
12       include "COMMON.CHAIN"
13       include "COMMON.NAMES"
14       include "COMMON.FFIELD"
15       include "COMMON.ENEPS"
16       include "COMMON.WEIGHTS"
17       include "COMMON.FREE"
18       include "COMMON.CONTROL"
19       include "COMMON.ENERGIES"
20       character*800 controlcard
21       integer i,j,k,ii,n_ene_found
22       integer ind,itype1,itype2,itypf,itypsc,itypp
23       integer ilen
24       external ilen
25       character*16 ucase
26       character*16 key
27       external ucase
28
29       call card_concat(controlcard,.true.)
30       call readi(controlcard,"N_ENE",n_ene,max_ene)
31       if (n_ene.gt.max_ene) then
32         write (iout,*) "Error: parameter out of range: N_ENE",n_ene,
33      &    max_ene
34         return1
35       endif
36       call readi(controlcard,"NPARMSET",nparmset,1)
37       separate_parset = index(controlcard,"SEPARATE_PARSET").gt.0
38       call readi(controlcard,"IPARMPRINT",iparmprint,1)
39       write (iout,*) "PARMPRINT",iparmprint
40       if (nparmset.gt.max_parm) then
41         write (iout,*) "Error: parameter out of range: NPARMSET",
42      &    nparmset, Max_Parm
43         return1
44       endif
45       call readi(controlcard,"MAXIT",maxit,5000)
46       call reada(controlcard,"FIMIN",fimin,1.0d-3)
47       call readi(controlcard,"ENSEMBLES",ensembles,0)
48       hamil_rep=index(controlcard,"HAMIL_REP").gt.0
49       write (iout,*) "Number of energy parameter sets",nparmset
50       call multreadi(controlcard,"ISAMPL",isampl,nparmset,1)
51       write (iout,*) "MaxSlice",MaxSlice
52       call readi(controlcard,"NSLICE",nslice,1)
53       call flush(iout)
54       if (nslice.gt.MaxSlice) then
55         write (iout,*) "Error: parameter out of range: NSLICE",nslice,
56      &    MaxSlice
57         return1
58       endif
59       write (iout,*) "Frequency of storing conformations",
60      & (isampl(i),i=1,nparmset)
61       write (iout,*) "Maxit",maxit," Fimin",fimin
62       call readi(controlcard,"NQ",nQ,1)
63       if (nQ.gt.MaxQ) then
64         write (iout,*) "Error: parameter out of range: NQ",nq,
65      &    maxq
66         return1
67       endif
68       indpdb=0
69       if (index(controlcard,"CLASSIFY").gt.0) indpdb=1
70       call reada(controlcard,"DELTA",delta,1.0d-2)
71       call readi(controlcard,"EINICHECK",einicheck,2)
72       call readi(controlcard,"NGRIDT",NGridT,400)
73       call reada(controlcard,"STARTGRIDT",StartGridT,200.0d0)
74       call reada(controlcard,"DELTA_T",Delta_T,1.0d0)
75       call reada(controlcard,"DELTRMS",deltrms,5.0d-2)
76       call reada(controlcard,"DELTRGY",deltrgy,5.0d-2)
77       call readi(controlcard,"RESCALE",rescale_mode,1)
78       check_conf=index(controlcard,"NO_CHECK_CONF").eq.0
79       call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
80       call readi(controlcard,'SYM',symetr,1)
81       write (iout,*) "DISTCHAINMAX",distchainmax
82       write (iout,*) "delta",delta
83       write (iout,*) "einicheck",einicheck
84       write (iout,*) "rescale_mode",rescale_mode
85       call flush(iout)
86       bxfile=index(controlcard,"BXFILE").gt.0
87       cxfile=index(controlcard,"CXFILE").gt.0
88       if (nslice .eq. 1 .and. .not.bxfile .and. .not.cxfile)
89      & bxfile=.true.
90       histfile=index(controlcard,"HISTFILE").gt.0
91       histout=index(controlcard,"HISTOUT").gt.0
92       entfile=index(controlcard,"ENTFILE").gt.0
93       zscfile=index(controlcard,"ZSCFILE").gt.0
94       with_dihed_constr = index(controlcard,"WITH_DIHED_CONSTR").gt.0
95       write (iout,*) "with_dihed_constr ",with_dihed_constr
96       call readi(controlcard,'CONSTR_DIST',constr_dist,0)
97       return
98       end
99 c------------------------------------------------------------------------------
100       subroutine read_efree(*)
101 C
102 C Read molecular data
103 C
104       implicit none
105       include 'DIMENSIONS'
106       include 'DIMENSIONS.ZSCOPT'
107       include 'DIMENSIONS.COMPAR'
108       include 'DIMENSIONS.FREE'
109       include 'COMMON.IOUNITS'
110       include 'COMMON.TIME1'
111       include 'COMMON.SBRIDGE'
112       include 'COMMON.CONTROL'
113       include 'COMMON.CHAIN'
114       include 'COMMON.HEADER'
115       include 'COMMON.GEO'
116       include 'COMMON.FREE'
117       character*320 controlcard,ucase
118       integer iparm,ib,i,j,npars
119       integer ilen
120       external ilen
121      
122       if (hamil_rep) then
123         npars=1
124       else
125         npars=nParmSet
126       endif
127
128       do iparm=1,npars
129
130       call card_concat(controlcard,.true.)
131       call readi(controlcard,'NT',nT_h(iparm),1)
132       write (iout,*) "iparm",iparm," nt",nT_h(iparm)
133       call flush(iout)
134       if (nT_h(iparm).gt.MaxT_h) then
135         write (iout,*)  "Error: parameter out of range: NT",nT_h(iparm),
136      &    MaxT_h
137         return1
138       endif
139       replica(iparm)=index(controlcard,"REPLICA").gt.0
140       umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
141       read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
142       write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
143      &  replica(iparm)," umbrella ",umbrella(iparm),
144      &  " read_iset",read_iset(iparm)
145       call flush(iout)
146       do ib=1,nT_h(iparm)
147         call card_concat(controlcard,.true.)
148         call readi(controlcard,'NR',nR(ib,iparm),1)
149         if (umbrella(iparm)) then
150           nRR(ib,iparm)=1
151         else
152           nRR(ib,iparm)=nR(ib,iparm)
153         endif
154         if (nR(ib,iparm).gt.MaxR) then
155           write (iout,*)  "Error: parameter out of range: NR",
156      &      nR(ib,iparm),MaxR
157         return1
158         endif
159         call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
160         beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
161         call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
162      &    0.0d0)
163         do i=1,nR(ib,iparm)
164           call card_concat(controlcard,.true.)
165           call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
166      &      100.0d0)
167           call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
168      &      0.0d0)
169         enddo
170       enddo
171       do ib=1,nT_h(iparm)
172         write (iout,*) "ib",ib," beta_h",
173      &    1.0d0/(0.001987*beta_h(ib,iparm))
174         write (iout,*) "nR",nR(ib,iparm)
175         write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
176         do i=1,nR(ib,iparm)
177           write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
178      &      "q0",(q0(j,i,ib,iparm),j=1,nQ)
179         enddo
180         call flush(iout)
181       enddo
182
183       enddo
184
185       if (hamil_rep) then
186
187        do iparm=2,nParmSet
188           nT_h(iparm)=nT_h(1)
189          do ib=1,nT_h(iparm)
190            nR(ib,iparm)=nR(ib,1)
191            if (umbrella(iparm)) then
192              nRR(ib,iparm)=1
193            else
194              nRR(ib,iparm)=nR(ib,1)
195            endif
196            beta_h(ib,iparm)=beta_h(ib,1)
197            do i=1,nR(ib,iparm)
198              f(i,ib,iparm)=f(i,ib,1)
199              do j=1,nQ
200                KH(j,i,ib,iparm)=KH(j,i,ib,1) 
201                Q0(j,i,ib,iparm)=Q0(j,i,ib,1) 
202              enddo
203            enddo
204            replica(iparm)=replica(1)
205            umbrella(iparm)=umbrella(1)
206            read_iset(iparm)=read_iset(1)
207          enddo
208        enddo
209         
210       endif
211
212       return
213       end
214 c-----------------------------------------------------------------------------
215       subroutine read_protein_data(*)
216       implicit none
217       include "DIMENSIONS"
218       include "DIMENSIONS.ZSCOPT"
219       include "DIMENSIONS.FREE"
220 #ifdef MPI
221       include "mpif.h"
222       integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
223       include "COMMON.MPI"
224 #endif
225       include "COMMON.CHAIN"
226       include "COMMON.IOUNITS"
227       include "COMMON.PROT"
228       include "COMMON.PROTFILES"
229       include "COMMON.NAMES"
230       include "COMMON.FREE"
231       include "COMMON.OBCINKA"
232       character*64 nazwa
233       character*16000 controlcard
234       integer i,ii,ib,iR,iparm,ilen,iroof,nthr,npars
235       external ilen,iroof
236       if (hamil_rep) then
237         npars=1
238       else
239         npars=nparmset
240       endif
241
242       do iparm=1,npars
243
244 C Read names of files with conformation data.
245       if (replica(iparm)) then
246         nthr = 1
247       else
248         nthr = nT_h(iparm)
249       endif
250       do ib=1,nthr
251       do ii=1,nRR(ib,iparm)
252       write (iout,*) "Parameter set",iparm," temperature",ib,
253      & " window",ii
254       call flush(iout)
255       call card_concat(controlcard,.true.) 
256       write (iout,*) controlcard(:ilen(controlcard))
257       call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
258       call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
259       call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
260       call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
261       call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
262      & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
263       call reada(controlcard,"TIME_START",
264      &  time_start_collect(ii,ib,iparm),0.0d0)
265       call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
266      &  1.0d10)
267       write (iout,*) "rec_start",rec_start(ii,ib,iparm),
268      & " rec_end",rec_end(ii,ib,iparm)
269       write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
270      & " time_end",time_end_collect(ii,ib,iparm)
271       call flush(iout)
272       if (replica(iparm)) then
273         call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
274         write (iout,*) "Number of trajectories",totraj(ii,iparm)
275         call flush(iout)
276       endif
277       if (nfile_bin(ii,ib,iparm).lt.2 
278      &    .and. nfile_asc(ii,ib,iparm).eq.0
279      &    .and. nfile_cx(ii,ib,iparm).eq.0) then
280         write (iout,*) "Error - no action specified!"
281         return1
282       endif
283       if (nfile_bin(ii,ib,iparm).gt.0) then
284         call card_concat(controlcard,.false.)
285         call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
286      &   maxfile_prot,nfile_bin(ii,ib,iparm))
287 #ifdef DEBUG
288         write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
289         write(iout,*) (protfiles(i,1,ii,ib,iparm),
290      &    i=1,nfile_bin(ii,ib,iparm))
291 #endif
292       endif
293       if (nfile_asc(ii,ib,iparm).gt.0) then
294         call card_concat(controlcard,.false.)
295         call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
296      &   maxfile_prot,nfile_asc(ii,ib,iparm))
297 #ifdef DEBUG
298         write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
299         write(iout,*) (protfiles(i,2,ii,ib,iparm),
300      &    i=1,nfile_asc(ii,ib,iparm))
301 #endif
302       else if (nfile_cx(ii,ib,iparm).gt.0) then
303         call card_concat(controlcard,.false.)
304         call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
305      &   maxfile_prot,nfile_cx(ii,ib,iparm))
306 #ifdef DEBUG
307         write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
308         write(iout,*) (protfiles(i,2,ii,ib,iparm),
309      &    i=1,nfile_cx(ii,ib,iparm))
310 #endif
311       endif
312       call flush(iout)
313       enddo
314       enddo
315
316       enddo
317
318       return
319       end
320 c-------------------------------------------------------------------------------
321       subroutine opentmp(islice,iunit,bprotfile_temp)
322       implicit none
323       include "DIMENSIONS"
324       include "DIMENSIONS.ZSCOPT"
325       include "DIMENSIONS.FREE"
326 #ifdef MPI
327       include "mpif.h"
328       integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
329       include "COMMON.MPI"
330 #endif
331       include "COMMON.IOUNITS"
332       include "COMMON.PROTFILES"
333       include "COMMON.PROT"
334       include "COMMON.FREE"
335       character*64 bprotfile_temp
336       character*3 liczba,liczba2
337       character*2 liczba1
338       integer iunit,islice
339       integer ilen,iroof
340       external ilen,iroof
341       logical lerr
342
343       write (liczba1,'(bz,i2.2)') islice
344       write (liczba,'(bz,i3.3)') me
345 #ifdef MPI
346 c      write (iout,*) "separate_parset ",separate_parset,
347 c     &  " myparm",myparm
348       if (separate_parset) then
349       write (liczba2,'(bz,i3.3)') myparm
350       bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
351      &  prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
352       open (iunit,file=bprotfile_temp,status="unknown",
353      &    form="unformatted",access="direct",recl=lenrec)
354       else
355       bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
356      &  prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
357       open (iunit,file=bprotfile_temp,status="unknown",
358      &    form="unformatted",access="direct",recl=lenrec)
359       endif
360 #else
361       bprotfile_temp = scratchdir(:ilen(scratchdir))//
362      &  "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
363       open (iunit,file=bprotfile_temp,status="unknown",
364      &    form="unformatted",access="direct",recl=lenrec)
365 #endif      
366 c      write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
367 c     &  bprotfile_temp
368 c      call flush(iout)
369       return
370       end
371 c-------------------------------------------------------------------------------
372       subroutine read_database(*)
373       implicit none
374       include "DIMENSIONS"
375       include "DIMENSIONS.ZSCOPT"
376       include "DIMENSIONS.FREE"
377 #ifdef MPI
378       include "mpif.h"
379       integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
380       include "COMMON.MPI"
381 #endif
382       include "COMMON.CHAIN"
383       include "COMMON.IOUNITS"
384       include "COMMON.PROTFILES"
385       include "COMMON.NAMES"
386       include "COMMON.VAR"
387       include "COMMON.GEO"
388       include "COMMON.ENEPS"
389       include "COMMON.PROT"
390       include "COMMON.INTERACT"
391       include "COMMON.FREE"
392       include "COMMON.SBRIDGE"
393       include "COMMON.OBCINKA"
394       real*4 csingle(3,maxres2)
395       character*64 nazwa,bprotfile_temp
396       character*3 liczba
397       character*2 liczba1
398       integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
399      &  ll(maxslice),mm(maxslice),if
400       integer nrec,nlines,iscor,iunit,islice
401       double precision energ
402       integer ilen,iroof
403       external ilen,iroof
404       double precision rmsdev,energia(0:max_ene),efree,eini,temp
405       double precision prop(maxQ)
406       integer ntot_all(maxslice,0:maxprocs-1)
407       integer iparm,ib,iib,ir,nprop,nthr,npars
408       double precision etot,time
409       integer ixdrf,iret 
410       logical lerr,linit
411
412       lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
413       lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
414       lenrec=lenrec2+8
415       write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
416      &  " lenrec2",lenrec2
417
418       do i=1,nQ
419         prop(i)=0.0d0
420       enddo
421       do islice=1,nslice
422         ll(islice)=0
423         mm(islice)=0
424       enddo
425       write (iout,*) "nparmset",nparmset
426       if (hamil_rep) then
427         npars=1
428       else
429         npars=nparmset
430       endif
431       do iparm=1,npars
432
433       if (replica(iparm)) then
434         nthr = 1
435       else
436         nthr = nT_h(iparm)
437       endif
438
439       do ib=1,nthr
440       do iR=1,nRR(ib,iparm)
441
442       write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
443       do islice=1,nslice
444         jj(islice)=0
445         kk(islice)=0
446       enddo
447
448       IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
449 c Read conformations from binary DA files (one per batch) and write them to 
450 c a binary DA scratchfile.
451         write (liczba,'(bz,i3.3)') me
452         do if=1,nfile_bin(iR,ib,iparm)
453           nazwa=protfiles(if,1,iR,ib,iparm)
454      &     (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
455           open (ientin,file=nazwa,status="old",form="unformatted",
456      &     access="direct",recl=lenrec2,err=1111)
457           ii=0
458           do islice=1,nslice
459             call opentmp(islice,ientout,bprotfile_temp)
460             call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
461      &        mm(islice),iR,ib,iparm)
462             close(ientout)
463           enddo
464           close(ientin)
465         enddo
466       ENDIF ! NFILE_BIN>0
467 c
468       IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
469 c Read conformations from multiple ASCII int files and write them to a binary
470 c DA scratchfile.
471         do if=1,nfile_asc(iR,ib,iparm)
472           nazwa=protfiles(if,2,iR,ib,iparm)
473      &     (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
474           open(unit=ientin,file=nazwa,status='old',err=1111)
475           write(iout,*) "reading ",nazwa(:ilen(nazwa))
476           ii=0
477           call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
478         enddo ! if
479       ENDIF
480       IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
481 c Read conformations from cx files and write them to a binary
482 c DA scratchfile.
483         do if=1,nfile_cx(iR,ib,iparm)
484           nazwa=protfiles(if,2,iR,ib,iparm)
485      &     (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
486           write(iout,*) "reading ",nazwa(:ilen(nazwa))
487           ii=0
488           print *,"Calling cxread"
489           call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
490      &       *1111)
491           close(ientout)
492           write (iout,*) "exit cxread"
493           call flush(iout)
494         enddo
495       ENDIF
496
497       do islice=1,nslice
498         stot(islice)=stot(islice)+jj(islice)
499       enddo
500
501       enddo
502       enddo
503       write (iout,*) "IPARM",iparm
504       enddo
505
506       if (nslice.eq.1) then
507 #ifdef MPI
508         write (liczba,'(bz,i3.3)') me
509         bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
510      &    prefix(:ilen(prefix))//liczba//".xbin.tmp"
511 #else
512         bprotfile_temp = scratchdir(:ilen(scratchdir))//
513      &     "/"//prefix(:ilen(prefix))//".xbin.tmp"
514 #endif
515         write(iout,*) mm(1)," conformations read",ll(1),
516      &    " conformations written to ", 
517      &    bprotfile_temp(:ilen(bprotfile_temp))
518       else
519         do islice=1,nslice
520           write (liczba1,'(bz,i2.2)') islice
521 #ifdef MPI
522           write (liczba,'(bz,i3.3)') me
523           bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
524      &      prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
525 #else
526           bprotfile_temp = scratchdir(:ilen(scratchdir))//
527      &       "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
528 #endif
529           write(iout,*) mm(islice)," conformations read",ll(islice),
530      &    " conformations written to ", 
531      &    bprotfile_temp(:ilen(bprotfile_temp))
532         enddo
533       endif
534
535 #ifdef MPI
536 c Check if everyone has the same number of conformations
537       call MPI_Allgather(stot(1),maxslice,MPI_INTEGER,
538      &  ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
539       lerr=.false.
540       do i=0,nprocs-1
541         if (i.ne.me) then
542           do islice=1,nslice
543           if (stot(islice).ne.ntot_all(islice,i)) then
544             write (iout,*) "Number of conformations at processor",i,
545      &       " differs from that at processor",me,
546      &       stot(islice),ntot_all(islice,i)," slice",islice
547             lerr = .true.
548           endif
549           enddo
550         endif
551       enddo 
552       if (lerr) then
553         write (iout,*)
554         write (iout,*) "Numbers of conformations read by processors"
555         write (iout,*)
556         do i=0,nprocs-1
557           write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
558         enddo
559         write (iout,*) "Calculation terminated."
560         call flush(iout)
561         return1
562       endif
563       do islice=1,nslice
564         ntot(islice)=stot(islice)
565       enddo
566       return
567 #endif
568  1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
569       call flush(iout)
570       return1
571       end
572 c------------------------------------------------------------------------------
573       subroutine card_concat(card,to_upper)
574       implicit none
575       include 'DIMENSIONS.ZSCOPT'
576       include "COMMON.IOUNITS"
577       character*(*) card
578       character*80 karta,ucase
579       logical to_upper
580       integer ilen
581       external ilen
582       read (inp,'(a)') karta
583       if (to_upper) karta=ucase(karta)
584       card=' '
585       do while (karta(80:80).eq.'&')
586         card=card(:ilen(card)+1)//karta(:79)
587         read (inp,'(a)') karta
588         if (to_upper) karta=ucase(karta)
589       enddo
590       card=card(:ilen(card)+1)//karta
591       return
592       end
593 c------------------------------------------------------------------------------
594       subroutine readi(rekord,lancuch,wartosc,default)
595       implicit none
596       character*(*) rekord,lancuch
597       integer wartosc,default
598       integer ilen,iread
599       external ilen
600       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
601       if (iread.eq.0) then
602         wartosc=default
603         return
604       endif
605       iread=iread+ilen(lancuch)+1
606       read (rekord(iread:),*) wartosc
607       return
608       end
609 c----------------------------------------------------------------------------
610       subroutine reada(rekord,lancuch,wartosc,default)
611       implicit none
612       character*(*) rekord,lancuch
613       character*80 aux
614       double precision wartosc,default
615       integer ilen,iread
616       external ilen
617       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
618       if (iread.eq.0) then
619         wartosc=default
620         return
621       endif
622       iread=iread+ilen(lancuch)+1
623       read (rekord(iread:),*) wartosc
624       return
625       end
626 c----------------------------------------------------------------------------
627       subroutine multreadi(rekord,lancuch,tablica,dim,default)
628       implicit none
629       integer dim,i
630       integer tablica(dim),default
631       character*(*) rekord,lancuch
632       character*80 aux
633       integer ilen,iread
634       external ilen
635       do i=1,dim
636         tablica(i)=default
637       enddo
638       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
639       if (iread.eq.0) return
640       iread=iread+ilen(lancuch)+1
641       read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
642    10 return
643       end
644 c----------------------------------------------------------------------------
645       subroutine multreada(rekord,lancuch,tablica,dim,default)
646       implicit none
647       integer dim,i
648       double precision tablica(dim),default
649       character*(*) rekord,lancuch
650       character*80 aux
651       integer ilen,iread
652       external ilen
653       do i=1,dim
654         tablica(i)=default
655       enddo
656       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
657       if (iread.eq.0) return
658       iread=iread+ilen(lancuch)+1
659       read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
660    10 return
661       end
662 c----------------------------------------------------------------------------
663       subroutine reads(rekord,lancuch,wartosc,default)
664       implicit none
665       character*(*) rekord,lancuch,wartosc,default
666       character*80 aux
667       integer ilen,lenlan,lenrec,iread,ireade
668       external ilen
669       logical iblnk
670       external iblnk
671       lenlan=ilen(lancuch)
672       lenrec=ilen(rekord)
673       iread=index(rekord,lancuch(:lenlan)//"=")
674 c      print *,"rekord",rekord," lancuch",lancuch
675 c      print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
676       if (iread.eq.0) then
677         wartosc=default
678         return
679       endif
680       iread=iread+lenlan+1
681 c      print *,"iread",iread
682 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
683       do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
684         iread=iread+1
685 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
686       enddo
687 c      print *,"iread",iread
688       if (iread.gt.lenrec) then
689          wartosc=default
690         return
691       endif
692       ireade=iread+1
693 c      print *,"ireade",ireade
694       do while (ireade.lt.lenrec .and.
695      &   .not.iblnk(rekord(ireade:ireade)))
696         ireade=ireade+1
697       enddo
698       wartosc=rekord(iread:ireade)
699       return
700       end
701 c----------------------------------------------------------------------------
702       subroutine multreads(rekord,lancuch,tablica,dim,default)
703       implicit none
704       integer dim,i
705       character*(*) rekord,lancuch,tablica(dim),default
706       character*80 aux
707       integer ilen,lenlan,lenrec,iread,ireade
708       external ilen
709       logical iblnk
710       external iblnk
711       do i=1,dim
712         tablica(i)=default
713       enddo
714       lenlan=ilen(lancuch)
715       lenrec=ilen(rekord)
716       iread=index(rekord,lancuch(:lenlan)//"=")
717 c      print *,"rekord",rekord," lancuch",lancuch
718 c      print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
719       if (iread.eq.0) return
720       iread=iread+lenlan+1
721       do i=1,dim
722 c      print *,"iread",iread
723 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
724       do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
725         iread=iread+1
726 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
727       enddo
728 c      print *,"iread",iread
729       if (iread.gt.lenrec) return
730       ireade=iread+1
731 c      print *,"ireade",ireade
732       do while (ireade.lt.lenrec .and.
733      &   .not.iblnk(rekord(ireade:ireade)))
734         ireade=ireade+1
735       enddo
736       tablica(i)=rekord(iread:ireade)
737       iread=ireade+1
738       enddo
739       end
740 c----------------------------------------------------------------------------
741       subroutine split_string(rekord,tablica,dim,nsub)
742       implicit none
743       integer dim,nsub,i,ii,ll,kk
744       character*(*) tablica(dim)
745       character*(*) rekord
746       integer ilen
747       external ilen
748       do i=1,dim
749         tablica(i)=" "
750       enddo
751       ii=1
752       ll = ilen(rekord)
753       nsub=0
754       do i=1,dim
755 C Find the start of term name
756         kk = 0
757         do while (ii.le.ll .and. rekord(ii:ii).eq." ") 
758           ii = ii+1
759         enddo
760 C Parse the name into TABLICA(i) until blank found
761         do while (ii.le.ll .and. rekord(ii:ii).ne." ") 
762           kk = kk+1
763           tablica(i)(kk:kk)=rekord(ii:ii)
764           ii = ii+1
765         enddo
766         if (kk.gt.0) nsub=nsub+1
767         if (ii.gt.ll) return
768       enddo
769       return
770       end
771 c--------------------------------------------------------------------------------
772       integer function iroof(n,m)
773       ii = n/m
774       if (ii*m .lt. n) ii=ii+1
775       iroof = ii
776       return
777       end