source/cluster/wham/src-M/readrtns.F

   1       subroutine read_control
   2 C
   3 C Read molecular data
   4 C
   5       implicit none
   6       include 'DIMENSIONS'
   7       include 'sizesclu.dat'
   8       include 'COMMON.IOUNITS'
   9       include 'COMMON.TIME1'
  10       include 'COMMON.SBRIDGE'
  11       include 'COMMON.CONTROL'
  12       include 'COMMON.CLUSTER'
  13       include 'COMMON.CHAIN'
  14       include 'COMMON.HEADER'
  15       include 'COMMON.FFIELD'
  16       include 'COMMON.FREE'
  17       include 'COMMON.INTERACT'
  18       character*320 controlcard,ucase
  19 #ifdef MPL
  20       include 'COMMON.INFO'
  21 #endif
  22       integer i
  23
  24       read (INP,'(a80)') titel
  25       call card_concat(controlcard)
  26
  27       call readi(controlcard,'NRES',nres,0)
  28       call readi(controlcard,'RESCALE',rescale_mode,2)
  29       call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
  30       write (iout,*) "DISTCHAINMAX",distchainmax
  31       call readi(controlcard,'PDBOUT',outpdb,0)
  32       call readi(controlcard,'MOL2OUT',outmol2,0)
  33       refstr=(index(controlcard,'REFSTR').gt.0)
  34       write (iout,*) "REFSTR",refstr
  35       pdbref=(index(controlcard,'PDBREF').gt.0)
  36       iscode=index(controlcard,'ONE_LETTER')
  37       tree=(index(controlcard,'MAKE_TREE').gt.0)
  38       min_var=(index(controlcard,'MINVAR').gt.0)
  39       plot_tree=(index(controlcard,'PLOT_TREE').gt.0)
  40       punch_dist=(index(controlcard,'PUNCH_DIST').gt.0)
  41       call readi(controlcard,'NCUT',ncut,1)
  42       call readi(controlcard,'SYM',symetr,1)
  43       write (iout,*) 'sym', symetr
  44       call readi(controlcard,'NSTART',nstart,0)
  45       call readi(controlcard,'NEND',nend,0)
  46       call reada(controlcard,'ECUT',ecut,10.0d0)
  47       call reada(controlcard,'PROB',prob_limit,0.99d0)
  48       write (iout,*) "Probability limit",prob_limit
  49       lgrp=(index(controlcard,'LGRP').gt.0)
  50       caonly=(index(controlcard,'CA_ONLY').gt.0)
  51       print_dist=(index(controlcard,'PRINT_DIST').gt.0)
  52       call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0)
  53       call readi(controlcard,'IOPT',iopt,2)
  54       lside = index(controlcard,"SIDE").gt.0
  55       efree = index(controlcard,"EFREE").gt.0
  56       call readi(controlcard,'NTEMP',nT,1)
  57       write (iout,*) "nT",nT
  58       call multreada(controlcard,'TEMPER',beta_h,nT,300.0d0)
  59       write (iout,*) "nT",nT
  60       write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
  61       do i=1,nT
  62         beta_h(i)=1.0d0/(1.987D-3*beta_h(i))
  63       enddo
  64       write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
  65       lprint_cart=index(controlcard,"PRINT_CART") .gt.0
  66       lprint_int=index(controlcard,"PRINT_INT") .gt.0
  67       if (min_var) iopt=1
  68       return
  69       end
  70 c--------------------------------------------------------------------------
  71       subroutine molread
  72 C
  73 C Read molecular data.
  74 C
  75       implicit none
  76       include 'DIMENSIONS'
  77       include 'COMMON.IOUNITS'
  78       include 'COMMON.GEO'
  79       include 'COMMON.VAR'
  80       include 'COMMON.INTERACT'
  81       include 'COMMON.LOCAL'
  82       include 'COMMON.NAMES'
  83       include 'COMMON.CHAIN'
  84       include 'COMMON.FFIELD'
  85       include 'COMMON.SBRIDGE'
  86       include 'COMMON.HEADER'
  87       include 'COMMON.CONTROL'
  88       include 'COMMON.CONTACTS'
  89       include 'COMMON.TIME1'
  90 #ifdef MPL
  91       include 'COMMON.INFO'
  92 #endif
  93       character*4 sequence(maxres)
  94       character*800 weightcard
  95       integer rescode
  96       double precision x(maxvar)
  97       integer itype_pdb(maxres)
  98       logical seq_comp
  99       integer i,j
 100 C
 101 C Body
 102 C
 103 C Read weights of the subsequent energy terms.
 104       call card_concat(weightcard)
 105       call reada(weightcard,'WSC',wsc,1.0d0)
 106       call reada(weightcard,'WLONG',wsc,wsc)
 107       call reada(weightcard,'WSCP',wscp,1.0d0)
 108       call reada(weightcard,'WELEC',welec,1.0D0)
 109       call reada(weightcard,'WVDWPP',wvdwpp,welec)
 110       call reada(weightcard,'WEL_LOC',wel_loc,1.0D0)
 111       call reada(weightcard,'WCORR4',wcorr4,0.0D0)
 112       call reada(weightcard,'WCORR5',wcorr5,0.0D0)
 113       call reada(weightcard,'WCORR6',wcorr6,0.0D0)
 114       call reada(weightcard,'WTURN3',wturn3,1.0D0)
 115       call reada(weightcard,'WTURN4',wturn4,1.0D0)
 116       call reada(weightcard,'WTURN6',wturn6,1.0D0)
 117       call reada(weightcard,'WSTRAIN',wstrain,1.0D0)
 118       call reada(weightcard,'WSCCOR',wsccor,1.0D0)
 119       call reada(weightcard,'WBOND',wbond,1.0D0)
 120       call reada(weightcard,'WTOR',wtor,1.0D0)
 121       call reada(weightcard,'WTORD',wtor_d,1.0D0)
 122       call reada(weightcard,'WANG',wang,1.0D0)
 123       call reada(weightcard,'WSCLOC',wscloc,1.0D0)
 124       call reada(weightcard,'SCAL14',scal14,0.4D0)
 125       call reada(weightcard,'SCALSCP',scalscp,1.0d0)
 126       call reada(weightcard,'CUTOFF',cutoff_corr,7.0d0)
 127       call reada(weightcard,'DELT_CORR',delt_corr,0.5d0)
 128       if (index(weightcard,'SOFT').gt.0) ipot=6
 129 C 12/1/95 Added weight for the multi-body term WCORR
 130       call reada(weightcard,'WCORRH',wcorr,1.0D0)
 131       if (wcorr4.gt.0.0d0) wcorr=wcorr4
 132       weights(1)=wsc
 133       weights(2)=wscp
 134       weights(3)=welec
 135       weights(4)=wcorr
 136       weights(5)=wcorr5
 137       weights(6)=wcorr6
 138       weights(7)=wel_loc
 139       weights(8)=wturn3
 140       weights(9)=wturn4
 141       weights(10)=wturn6
 142       weights(11)=wang
 143       weights(12)=wscloc
 144       weights(13)=wtor
 145       weights(14)=wtor_d
 146       weights(15)=wstrain
 147       weights(16)=wvdwpp
 148       weights(17)=wbond
 149       weights(18)=scal14
 150       weights(19)=wsccor
 151       write (iout,10) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,
 152      &  wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wturn3,
 153      &  wturn4,wturn6,wsccor
 154    10 format (/'Energy-term weights (unscaled):'//
 155      & 'WSCC=   ',f10.6,' (SC-SC)'/
 156      & 'WSCP=   ',f10.6,' (SC-p)'/
 157      & 'WELEC=  ',f10.6,' (p-p electr)'/
 158      & 'WVDWPP= ',f10.6,' (p-p VDW)'/
 159      & 'WBOND=  ',f10.6,' (stretching)'/
 160      & 'WANG=   ',f10.6,' (bending)'/
 161      & 'WSCLOC= ',f10.6,' (SC local)'/
 162      & 'WTOR=   ',f10.6,' (torsional)'/
 163      & 'WTORD=  ',f10.6,' (double torsional)'/
 164      & 'WSTRAIN=',f10.6,' (SS bridges & dist. cnstr.)'/
 165      & 'WEL_LOC=',f10.6,' (multi-body 3-rd order)'/
 166      & 'WCORR4= ',f10.6,' (multi-body 4th order)'/
 167      & 'WCORR5= ',f10.6,' (multi-body 5th order)'/
 168      & 'WCORR6= ',f10.6,' (multi-body 6th order)'/
 169      & 'WTURN3= ',f10.6,' (turns, 3rd order)'/
 170      & 'WTURN4= ',f10.6,' (turns, 4th order)'/
 171      & 'WTURN6= ',f10.6,' (turns, 6th order)'/
 172      & 'WSCCOR= ',f10.6,' (SC-backbone torsinal correalations)')
 173
 174       if (wcorr4.gt.0.0d0) then
 175         write (iout,'(/2a/)') 'Local-electrostatic type correlation ',
 176      &   'between contact pairs of peptide groups'
 177         write (iout,'(2(a,f5.3/))')
 178      &  'Cutoff on 4-6th order correlation terms: ',cutoff_corr,
 179      &  'Range of quenching the correlation terms:',2*delt_corr
 180       else if (wcorr.gt.0.0d0) then
 181         write (iout,'(/2a/)') 'Hydrogen-bonding correlation ',
 182      &   'between contact pairs of peptide groups'
 183       endif
 184       write (iout,'(a,f8.3)')
 185      &  'Scaling factor of 1,4 SC-p interactions:',scal14
 186       write (iout,'(a,f8.3)')
 187      &  'General scaling factor of SC-p interactions:',scalscp
 188       r0_corr=cutoff_corr-delt_corr
 189       do i=1,20
 190         aad(i,1)=scalscp*aad(i,1)
 191         aad(i,2)=scalscp*aad(i,2)
 192         bad(i,1)=scalscp*bad(i,1)
 193         bad(i,2)=scalscp*bad(i,2)
 194       enddo
 195
 196       call flush(iout)
 197       print *,'indpdb=',indpdb,' pdbref=',pdbref
 198
 199 C Read sequence if not taken from the pdb file.
 200       if (iscode.gt.0) then
 201         read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
 202       else
 203         read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
 204       endif
 205 C Convert sequence to numeric code
 206       do i=1,nres
 207         itype(i)=rescode(i,sequence(i),iscode)
 208       enddo
 209       if (itype(2).eq.10.and.itype(1).eq.ntyp1) then
 210         write (iout,*)
 211      &   "Glycine is the first full residue, initial dummy deleted"
 212         do i=1,nres
 213           itype(i)=itype(i+1)
 214         enddo
 215         nres=nres-1
 216       endif
 217       if (itype(nres-1).eq.10.and.itype(nres).eq.ntyp1) then
 218         write (iout,*)
 219      &   "Glycine is the last full residue, terminal dummy deleted"
 220         nres=nres-1
 221       endif
 222       print *,nres
 223       print '(20i4)',(itype(i),i=1,nres)
 224
 225       do i=1,nres
 226 #ifdef PROCOR
 227         if (itype(i).eq.21 .or. itype(i+1).eq.21) then
 228 #else
 229         if (itype(i).eq.21) then
 230 #endif
 231           itel(i)=0
 232 #ifdef PROCOR
 233         else if (itype(i+1).ne.20) then
 234 #else
 235         else if (itype(i).ne.20) then
 236 #endif
 237           itel(i)=1
 238         else
 239           itel(i)=2
 240         endif
 241       enddo
 242       write (iout,*) "ITEL"
 243       do i=1,nres-1
 244         write (iout,*) i,itype(i),itel(i)
 245       enddo
 246
 247       print *,'Call Read_Bridge.'
 248       call read_bridge
 249       nnt=1
 250       nct=nres
 251       print *,'NNT=',NNT,' NCT=',NCT
 252       if (itype(1).eq.21) nnt=2
 253       if (itype(nres).eq.21) nct=nct-1
 254       if (nstart.lt.nnt) nstart=nnt
 255       if (nend.gt.nct .or. nend.eq.0) nend=nct
 256       write (iout,*) "nstart",nstart," nend",nend
 257       nres0=nres
 258 c      if (pdbref) then
 259 c        read(inp,'(a)') pdbfile
 260 c        write (iout,'(2a)') 'PDB data will be read from file ',pdbfile
 261 c        open(ipdbin,file=pdbfile,status='old',err=33)
 262 c        goto 34
 263 c  33    write (iout,'(a)') 'Error opening PDB file.'
 264 c        stop
 265 c  34    continue
 266 c        print *,'Begin reading pdb data'
 267 c        call readpdb
 268 c        print *,'Finished reading pdb data'
 269 c        write (iout,'(a,i3,a,i3)')'nsup=',nsup,' nstart_sup=',nstart_sup
 270 c        do i=1,nres
 271 c          itype_pdb(i)=itype(i)
 272 c        enddo
 273 c        close (ipdbin)
 274 c        write (iout,'(a,i3)') 'nsup=',nsup
 275 c        nstart_seq=nnt
 276 c        if (nsup.le.(nct-nnt+1)) then
 277 c          do i=0,nct-nnt+1-nsup
 278 c            if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),nsup)) then
 279 c              nstart_seq=nnt+i
 280 c              goto 111
 281 c            endif
 282 c          enddo
 283 c          write (iout,'(a)')
 284 c     &            'Error - sequences to be superposed do not match.'
 285 c          stop
 286 c        else
 287 c          do i=0,nsup-(nct-nnt+1)
 288 c            if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),nct-nnt+1))
 289 c     &      then
 290 c              nstart_sup=nstart_sup+i
 291 c              nsup=nct-nnt+1
 292 c              goto 111
 293 c            endif
 294 c          enddo
 295 c          write (iout,'(a)')
 296 c     &            'Error - sequences to be superposed do not match.'
 297 c        endif
 298 c  111   continue
 299 c        write (iout,*) 'nsup=',nsup,' nstart_sup=',nstart_sup,
 300 c     &                 ' nstart_seq=',nstart_seq
 301 c      endif
 302       call init_int_table
 303       call setup_var
 304       write (iout,*) "molread: REFSTR",refstr
 305       if (refstr) then
 306         if (.not.pdbref) then
 307           call read_angles(inp,*38)
 308           goto 39
 309    38     write (iout,'(a)') 'Error reading reference structure.'
 310 #ifdef MPL
 311           call mp_stopall(Error_Msg)
 312 #else
 313           stop 'Error reading reference structure'
 314 #endif
 315    39     call chainbuild
 316           nstart_sup=nnt
 317           nstart_seq=nnt
 318           nsup=nct-nnt+1
 319           do i=1,2*nres
 320             do j=1,3
 321               cref(j,i)=c(j,i)
 322             enddo
 323           enddo
 324         endif
 325         call contact(.true.,ncont_ref,icont_ref)
 326       endif
 327       return
 328       end
 329 c-----------------------------------------------------------------------------
 330       logical function seq_comp(itypea,itypeb,length)
 331       implicit none
 332       integer length,itypea(length),itypeb(length)
 333       integer i
 334       do i=1,length
 335         if (itypea(i).ne.itypeb(i)) then
 336           seq_comp=.false.
 337           return
 338         endif
 339       enddo
 340       seq_comp=.true.
 341       return
 342       end
 343 c-----------------------------------------------------------------------------
 344       subroutine read_bridge
 345 C Read information about disulfide bridges.
 346       implicit none
 347       include 'DIMENSIONS'
 348       include 'COMMON.IOUNITS'
 349       include 'COMMON.GEO'
 350       include 'COMMON.VAR'
 351       include 'COMMON.INTERACT'
 352       include 'COMMON.LOCAL'
 353       include 'COMMON.NAMES'
 354       include 'COMMON.CHAIN'
 355       include 'COMMON.FFIELD'
 356       include 'COMMON.SBRIDGE'
 357       include 'COMMON.HEADER'
 358       include 'COMMON.CONTROL'
 359       include 'COMMON.TIME1'
 360 #ifdef MPL
 361       include 'COMMON.INFO'
 362 #endif
 363       integer i,j
 364 C Read bridging residues.
 365       read (inp,*) ns,(iss(i),i=1,ns)
 366       print *,'ns=',ns
 367 C Check whether the specified bridging residues are cystines.
 368       do i=1,ns
 369         if (itype(iss(i)).ne.1) then
 370           write (iout,'(2a,i3,a)')
 371      &   'Do you REALLY think that the residue ',restyp(iss(i)),i,
 372      &   ' can form a disulfide bridge?!!!'
 373           write (*,'(2a,i3,a)')
 374      &   'Do you REALLY think that the residue ',restyp(iss(i)),i,
 375      &   ' can form a disulfide bridge?!!!'
 376 #ifdef MPL
 377          call mp_stopall(error_msg)
 378 #else
 379          stop
 380 #endif
 381         endif
 382       enddo
 383 C Read preformed bridges.
 384       if (ns.gt.0) then
 385       read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
 386       if (nss.gt.0) then
 387         nhpb=nss
 388 C Check if the residues involved in bridges are in the specified list of
 389 C bridging residues.
 390         do i=1,nss
 391           do j=1,i-1
 392             if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j)
 393      &      .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then
 394               write (iout,'(a,i3,a)') 'Disulfide pair',i,
 395      &      ' contains residues present in other pairs.'
 396               write (*,'(a,i3,a)') 'Disulfide pair',i,
 397      &      ' contains residues present in other pairs.'
 398 #ifdef MPL
 399               call mp_stopall(error_msg)
 400 #else
 401               stop
 402 #endif
 403             endif
 404           enddo
 405           do j=1,ns
 406             if (ihpb(i).eq.iss(j)) goto 10
 407           enddo
 408           write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
 409    10     continue
 410           do j=1,ns
 411             if (jhpb(i).eq.iss(j)) goto 20
 412           enddo
 413           write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
 414    20     continue
 415           dhpb(i)=dbr
 416           forcon(i)=fbr
 417         enddo
 418         do i=1,nss
 419           ihpb(i)=ihpb(i)+nres
 420           jhpb(i)=jhpb(i)+nres
 421         enddo
 422       endif
 423       endif
 424       return
 425       end
 426 c----------------------------------------------------------------------------
 427       subroutine read_angles(kanal,*)
 428       implicit none
 429       include 'DIMENSIONS'
 430       include 'COMMON.GEO'
 431       include 'COMMON.VAR'
 432       include 'COMMON.CHAIN'
 433       include 'COMMON.IOUNITS'
 434       integer i,kanal
 435       read (kanal,*,err=10,end=10) (theta(i),i=3,nres)
 436       read (kanal,*,err=10,end=10) (phi(i),i=4,nres)
 437       read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1)
 438       read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1)
 439       do i=1,nres
 440         theta(i)=deg2rad*theta(i)
 441         phi(i)=deg2rad*phi(i)
 442         alph(i)=deg2rad*alph(i)
 443         omeg(i)=deg2rad*omeg(i)
 444       enddo
 445       return
 446    10 return1
 447       end
 448 c----------------------------------------------------------------------------
 449       subroutine reada(rekord,lancuch,wartosc,default)
 450       implicit none
 451       character*(*) rekord,lancuch
 452       double precision wartosc,default
 453       integer ilen,iread
 454       external ilen
 455       iread=index(rekord,lancuch)
 456       if (iread.eq.0) then
 457         wartosc=default
 458         return
 459       endif
 460       iread=iread+ilen(lancuch)+1
 461       read (rekord(iread:),*) wartosc
 462       return
 463       end
 464 c----------------------------------------------------------------------------
 465       subroutine multreada(rekord,lancuch,tablica,dim,default)
 466       implicit none
 467       integer dim,i
 468       double precision tablica(dim),default
 469       character*(*) rekord,lancuch
 470       integer ilen,iread
 471       external ilen
 472       do i=1,dim
 473         tablica(i)=default
 474       enddo
 475       iread=index(rekord,lancuch)
 476       if (iread.eq.0) return
 477       iread=iread+ilen(lancuch)+1
 478       read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
 479    10 return
 480       end
 481 c----------------------------------------------------------------------------
 482       subroutine readi(rekord,lancuch,wartosc,default)
 483       implicit none
 484       character*(*) rekord,lancuch
 485       integer wartosc,default
 486       integer ilen,iread
 487       external ilen
 488       iread=index(rekord,lancuch)
 489       if (iread.eq.0) then
 490         wartosc=default
 491         return
 492       endif
 493       iread=iread+ilen(lancuch)+1
 494       read (rekord(iread:),*) wartosc
 495       return
 496       end
 497 c----------------------------------------------------------------------------
 498       subroutine card_concat(card)
 499       include 'DIMENSIONS'
 500       include 'COMMON.IOUNITS'
 501       character*(*) card
 502       character*80 karta,ucase
 503       external ilen
 504       read (inp,'(a)') karta
 505       karta=ucase(karta)
 506       card=' '
 507       do while (karta(80:80).eq.'&')
 508         card=card(:ilen(card)+1)//karta(:79)
 509         read (inp,'(a)') karta
 510         karta=ucase(karta)
 511       enddo
 512       card=card(:ilen(card)+1)//karta
 513       return
 514       end
 515 c----------------------------------------------------------------------------
 516       subroutine openunits
 517       implicit none
 518       include 'DIMENSIONS'
 519 #ifdef MPI
 520       include "mpif.h"
 521       character*3 liczba
 522       include "COMMON.MPI"
 523 #endif
 524       include 'COMMON.IOUNITS'
 525       include 'COMMON.CONTROL'
 526       integer lenpre,lenpot,ilen
 527       external ilen
 528       character*16 cformat,cprint
 529       character*16 ucase
 530       integer lenint,lenout
 531       call getenv('INPUT',prefix)
 532       call getenv('OUTPUT',prefout)
 533       call getenv('INTIN',prefintin)
 534       call getenv('COORD',cformat)
 535       call getenv('PRINTCOOR',cprint)
 536       call getenv('SCRATCHDIR',scratchdir)
 537       from_bx=.true.
 538       from_cx=.false.
 539       if (index(ucase(cformat),'CX').gt.0) then
 540         from_cx=.true.
 541         from_bx=.false.
 542       endif
 543       from_cart=.true.
 544       lenpre=ilen(prefix)
 545       lenout=ilen(prefout)
 546       lenint=ilen(prefintin)
 547 C Get the names and open the input files
 548       open (inp,file=prefix(:ilen(prefix))//'.inp',status='old')
 549 #ifdef MPI
 550       write (liczba,'(bz,i3.3)') me
 551       outname=prefout(:lenout)//'_clust.out_'//liczba
 552 #else
 553       outname=prefout(:lenout)//'_clust.out'
 554 #endif
 555       if (from_bx) then
 556         intinname=prefintin(:lenint)//'.bx'
 557       else if (from_cx) then
 558         intinname=prefintin(:lenint)//'.cx'
 559       else
 560         intinname=prefintin(:lenint)//'.int'
 561       endif
 562       rmsname=prefintin(:lenint)//'.rms'
 563       open (jplot,file=prefout(:ilen(prefout))//'.tex',
 564      &       status='unknown')
 565       open (jrms,file=rmsname,status='unknown')
 566       open(iout,file=outname,status='unknown')
 567 C Get parameter filenames and open the parameter files.
 568       call getenv('BONDPAR',bondname)
 569       open (ibond,file=bondname,status='old')
 570       call getenv('THETPAR',thetname)
 571       open (ithep,file=thetname,status='old')
 572       call getenv('ROTPAR',rotname)
 573       open (irotam,file=rotname,status='old')
 574       call getenv('TORPAR',torname)
 575       open (itorp,file=torname,status='old')
 576       call getenv('TORDPAR',tordname)
 577       open (itordp,file=tordname,status='old')
 578       call getenv('FOURIER',fouriername)
 579       open (ifourier,file=fouriername,status='old')
 580       call getenv('ELEPAR',elename)
 581       open (ielep,file=elename,status='old')
 582       call getenv('SIDEPAR',sidename)
 583       open (isidep,file=sidename,status='old')
 584       call getenv('SIDEP',sidepname)
 585       open (isidep1,file=sidepname,status="old")
 586       call getenv('SCCORPAR',sccorname)
 587       open (isccor,file=sccorname,status="old")
 588 #ifndef OLDSCP
 589 C
 590 C 8/9/01 In the newest version SCp interaction constants are read from a file
 591 C Use -DOLDSCP to use hard-coded constants instead.
 592 C
 593       call getenv('SCPPAR',scpname)
 594       open (iscpp,file=scpname,status='old')
 595 #endif
 596       return
 597       end