1 subroutine ana_rdf_init() 2c 3c $Id$ 4c 5 implicit none 6c 7#include "ana_common.fh" 8#include "global.fh" 9#include "mafdecls.fh" 10#include "msgids.fh" 11#include "rtdb.fh" 12c 13 if(.not.ma_push_get(mt_int,nrdf*nsel*mwa,'irdf',l_rdf,i_rdf)) 14 + call md_abort('Could not allocate irdf',0) 15 print*,'rdf allocated in rdfhdr ',nrdf*nsel*mwa 16c 17 return 18 end 19 subroutine ana_rdfhdr(irdf) 20c 21 implicit none 22c 23#include "ana_common.fh" 24#include "global.fh" 25#include "mafdecls.fh" 26#include "msgids.fh" 27#include "rtdb.fh" 28c 29 integer irdf(nsel,mwa,nrdf) 30c 31 character*255 fname 32 integer i,j,k,lq,l,m 33c 34 fname=filtrj 35 lq=index(filtrj,'?') 36 if(lq.gt.0) then 37 fname=filtrj(1:lq-1)//cnum//filtrj(lq+1:index(filtrj,' ')-1) 38 endif 39 lq=index(fname,'.trj') 40 fname(lq:lq+3)='.rdf' 41c 42 open(unit=lfnrdf,file=fname(1:index(fname,' ')-1), 43 + status='unknown') 44c 45 write(*,3333) fname(1:index(fname,' ')-1) 46 3333 format(/,' Opening rdf file ',a) 47 rewind(lfnrdf) 48c 49 numrdf=0 50c 51 do 1 i=1,nrdf 52 do 2 j=1,mwa 53 do 3 k=1,nsel 54 irdf(k,j,i)=0 55 3 continue 56 2 continue 57 1 continue 58c 59 print*,'rdf set to zero in rdfhdr' 60c 61 return 62 end 63 subroutine ana_rdf(isel,xs,xw,irdf) 64c 65 implicit none 66c 67#include "ana_common.fh" 68#include "global.fh" 69#include "mafdecls.fh" 70#include "msgids.fh" 71#include "rtdb.fh" 72c 73 integer isel(msa),irdf(nsel,mwa,nrdf) 74 real*8 xs(msa,3),xw(mwm,mwa,3) 75 integer i,j,k,l,m 76 real*8 d 77c 78 i=0 79 do 1 j=1,nsa 80 if(isel(j).gt.0) then 81 i=i+1 82 do 2 k=1,nwm 83 do 3 l=1,nwa 84 d=sqrt((xs(j,1)-xw(k,l,1))**2+(xs(j,2)-xw(k,l,2))**2+ 85 + (xs(j,3)-xw(k,l,3))**2) 86 m=int(dble(nrdf*rrdf)/d) 87 if(m.le.nrdf) irdf(i,l,m)=irdf(i,l,m)+1 88 3 continue 89 2 continue 90 endif 91 1 continue 92c 93 numrdf=numrdf+1 94c 95 return 96 end 97 subroutine ana_rdfwrt(isel,irdf) 98c 99 implicit none 100c 101#include "ana_common.fh" 102#include "global.fh" 103#include "mafdecls.fh" 104#include "msgids.fh" 105#include "rtdb.fh" 106c 107 integer isel(msa) 108 integer irdf(nsel,mwa,nrdf) 109c 110 integer i,j,k,l,m 111 real*8 d 112c 113 write(*,'(a)') 'ANA_RDFWRT' 114c 115 do 4 i=1,nsel 116 do 5 m=1,nrdf 117 write(*,'(2i5,10i10)') i,m,(irdf(i,l,m),l=1,nwa) 118 5 continue 119 4 continue 120c 121c i=0 122c do 1 j=1,nsa 123c if(isel(j).gt.0) then 124c i=i+1 125c do 2 k=1,nrdf 126c d=dble(k)*rrdf/dble(nrdf) 127c write(*,1000) d,(irdf(i,l,k),l=1,nwa) 128c 1000 format(f12.6,10i5) 129c 2 continue 130c endif 131c 1 continue 132c 133 close(unit=lfnrdf,status='keep') 134 write(*,'(a)') ' Closing rdf file ' 135c 136 return 137 end 138