Duplicate Points In Overlapping Section of Two Overlapping Spheres

I have two overlapping spheres, and I am using the Monte Carlo method to calculate the volume of the combined spheres. Everything is working except for one part. I need to find a way to NOT double count points placed in the overlapping region so that I can correctly calculate the ratio of points inside the spheres to the points outside in the cube. Can someone help me with how not to double count points in the overlapping region? Here is my code so far:

clc 
clear all

S=input('Number of Simulations: ');
N=input('Number of Points: ');
fileID = fopen('Volume.txt','w');

for k = 1:S

[X, Y, Z]=sphere;
a = [1 1 1 1.2; 2.293,1,1,1.8];

%Draw sphere #1
figure('color','white');
s1=surf(X*a(1,4)+a(1,1), Y*a(1,4)+a(1,2), Z*a(1,4)+a(1,3),'FaceColor', [1 1 1],'edgecolor','none','FaceAlpha',0.6);
light
axis equal
set(gca,'Color','k')
hold on
%Draw sphere #2
s2=surf(X*a(2,4)+a(2,1), Y*a(2,4)+a(2,2), Z*a(2,4)+a(2,3),'FaceColor', [0 1 0],'edgecolor','none','FaceAlpha',0.5);
daspect([1 1 1])
view(30,10)
xlabel('x')
ylabel('y')
zlabel('z')

x1range1 = -0.2;
x1range2 = 4.2;
y1range1 = -0.9;
y1range2 = 2.9;
z1range1 = -0.9;
z1range2 = 2.9;

 hit1 = 0; %Points inside spheres
 hit2 = 0;
 
 for i = 1:N
     
     %Define Box Range and Randomize Points
     x = (x1range2-x1range1).*rand(N,1)+x1range1;
     y = (y1range2-y1range1).*rand(N,1)+y1range1;
     z = (z1range2-z1range1).*rand(N,1)+z1range1;
     
     x1=x-1;
     y1=y-1;
     z1=z-1;
     r1=sqrt(x1.^2+y1.^2+z1.^2);
     
     x2=x-2.293;
     y2=y-1;
     z2=z-1;
     r2=sqrt(x2.^2+y2.^2+z2.^2);
  
     
ii = r1<=1.2; 
jj = r2<=1.8;

hit1 = sum(ii);
hit2 = sum(jj);

plot3(x(ii), y(ii), z(ii), 'w+');
%plot3(x(~ii), y(~ii), z(~ii), 'r+');
plot3(x(jj), y(jj), z(jj), 'g+');
drawnow

     
 
 end
 
 SpherePercentVolume = (hit1+hit2)/N;
 Volume = ((x1range2-x1range1)*(y1range2-y1range1)*(z1range2-z1range1))*SpherePercentVolume
 fprintf(fileID,'%6.3f\r\n',Volume);
end

fclose(fileID);

load Volume.txt
average = mean(Volume);
stdev = std(Volume);
fileID = fopen('AvgVolume.txt','w');
fprintf(fileID,'%6s','Average Volume: ');
fprintf(fileID,'%6.3f\r\n',average);
fprintf(fileID,'%6s','Standard Deviation: ');
fprintf(fileID,'%6.3f',stdev);

[David Goodmanson]

Hi Kelly,

I don't know what the for loop from 1 to N does, but aside from that there are various ways to accomplish this with some boolian operations.

ii & jj indexes the points in both spheres, which you can use to correct for double counting.

ii & ~jj indexes the points in sphere 1 and not in sphere 2, etc.

[Kelly McGuire]

Hey David, thanks for the response. I don't quite have the experience yet with programming to use boolean operations in this way. Would you be able to show me in the posted code? The for loop from 1 to N just defines the cube ranges, the centers of the spheres, and then with ii and jj defines which points fall inside the spheres and then plots them as green and white. The hit1 and hit2 summations tell me how many points are in the white sphere and how many are in the green sphere so I can take the ratio of those to the total number of points. N is just the number of simulations I want to run. If N is 10, then it will create 10 plots of 10 different simulations.