%Import the files
%assume you are in ~/protomol/examples/alanDipeptideSol
clear;
DELIMITER = ' ';
HEADERLINES = 0;
ad = importdata('ESM/alanine.EW.force.xyz', DELIMITER, HEADERLINES);  % explicit solvent, Hassan05 charges
c  = importdata('ISM_BORN/alanine.EW.forceBorn.xyz', DELIMITER, HEADERLINES); % dielectric D = 53.1,S=0.25
c1  = importdata('ISM_NOBORN/alanine.EW.forceSCPISM.xyz', DELIMITER, HEADERLINES); % dielectric D=53,S=0.20
m  = ad(:,4); 
ad = ad(:,1:3);
c = c(:,1:3);
c1 = c1(:,1:3);
l = length(m);

% denominator
d=zeros(l,1);
for i=1:l
d(i)=norm(ad(i,:)) * m(i);
end
n=sum(d)/l;

% numerator
diff = c - ad;
diffn = zeros(l,1);
for i=1:l
diffn(i)=norm(diff(i,:)) * m(i);
end

% max error for coulomb dielec
maxfc = max(diffn) / n;

%avg error for coulomb dielec
avgfc = mean(diffn) / n;

maxfc

avgfc

% numerator
diff = c1 - ad;
diffn = zeros(l,1);
for i=1:l
diffn(i)=norm(diff(i,:)) * m(i);
end

% max error for coulomb dielec
maxfc1 = max(diffn) / n;

%avg error for coulomb dielec
avgfc1 = mean(diffn) / n;

maxfc1

avgfc1


