%psiav 
% average PSI for empirical water calculation
% reads files named <fname>.<XXXX>.PSI where fname is given below and
% XXXX is (Tinker) frame number, from 0001 to Nmax
% It calculates 'running' averages of N2av frames and saves them
% as 'psi<N2av><XXXX>.PSI' where XXXX is the number of the saved
% average.
% 'step' defines how many frames are left out(skipped) between the individual
% running averages (step=1 will do all the running averages, step=10
% from 1 to N2av, from 11 to N2av+10, from 21 to N2av+20 etc.)
%
clear all;
%
Nmax=10000; % number of time slices/files
N2av=59; % number(s) to average
step=1; % step for averaging
fn='p22full.'; % filename
%
nstartstart=1; % start for each
%
nstart=nstartstart;
nend=nstart+N2av-1;
nav=nstart; % counter for ea. average block
%
while nend<=Nmax
   disp('avg. # out of');
   disp([nav,floor((Nmax-N2av+1)/step)+1]);
   %
   for iav=nstart:nend
      %disp('@ file out of');
      %disp([iav nend]);
      %
      % put together file index
      if iav<1000
         if iav<10
          ifn=strcat('000',int2str(iav));
         elseif iav<100
          ifn=strcat('00',int2str(iav));
         else 
          ifn=strcat('0',int2str(iav));
         end %if
      else
          ifn=int2str(iav);
      end %if
      % load PSI file
      psi=load(strcat(fn,ifn,'.PSI')); 
      % average
      if iav==nstartstart
         psiav=psi./N2av;
      else
      psiav=psiav+psi./N2av;
      end %if

   end %for

   %
   % write the avg out
   cd('av59');
   %
   fid=fopen(strcat('psi',int2str(N2av),'_',int2str(nav),'.PSI'),'w');
   fprintf(fid,'%12.6f %12.6f %12.6f %12.6f %12.6f %12.6f \n',psiav');
   fclose(fid);
   %
   cd('../');
   %
   % move to the next
   nstart=nstart+step;
   nend=nend+step;
   nav=nav+1;
   % zero out average fields
   psiav=zeros(size(psi));

end %while

fclose all;
%

