QA0001
ZIP with notes and MATLAB script:
Q: A start point script is needed to simulate metal sheet damage location as described in B.Zima M.Rucka article
notes:
A:
the previous ellipse related question was about drawing ellipses given 3 random points, all 3 belonging to the group of points actually drawing ellipse curve, and we saw that given 3 random ellipse points, there are many different ellipses meeting an ellipse equation and containing those 3 points, because major and minor axes are not actually defined.
However, for this question, given 2 ellipse centre points and 1 ellipse curve point, there's only 1 ellipse
meeting the now uniquely defined ellipse equation, that has a well defined single pair of major and minor axes.
% ASSUMPTIONS
% homogeneous metal
% single scratch
% sensors do not overlap
% no sensor right on scratch, this is why manual input of sensor locations
% here d=dt*cg in Zima Rucka app note
% DEFINE PLATE
clear all;clc;close all
hf1=figure(1);grid on
ax1=hf1.CurrentAxes
ax1.DataAspectRatio=[1 1 1]
hold all
x_span=100;y_span=50 % half size of [0 0] centred workspace
x_res=.1 % spatial resolution,on both x ad y
x_range=[-x_span:x_res:x_span] % test area
y_range=[-y_span:x_res:y_span]
plot([x_range(1) x_range(end) x_range(end) x_range(1) x_range(1)],...
[y_range(1) y_range(1) y_range(end) y_range(end) y_range(1)],...
'LineWidth',2,...
'Color','b')
ax1.XTick=[x_range(1):10*x_res:x_range(end)]
ax1.YTick=[y_range(1):10*x_res:y_range(end)]
% SIMULATE DAMAGE
Nd=randi([3 25],1,1) % let the damage be a single straight scratch of random length range [5 25]
d_ang_res=2*pi/360 % let be the scratch attitude a random angle with resolution 1 degree
d_ang_range=[0:d_ang_res:2*pi-d_ang_res];
d_ang_target=d_ang_range(randi(numel([d_ang_range]),1,1))
x_d=x_range(randi(numel(x_range),1,1)) % plate damage location
y_d=y_range(randi(numel(y_range),1,1))
plot([x_d-Nd*x_res x_d+Nd*x_res],[y_d-Nd*x_res y_d+Nd*x_res],...
'LineWidth',1,...
'Color','k')
S=.1*floor(10*[x_d-Nd*x_res:x_res:x_d+Nd*x_res;y_d-Nd*x_res:x_res:y_d+Nd*x_res]) % scratch points
NS=5 % amount sensors, there are NS sensors
% x_s=x_range(randi(numel(x_range),1,NS)) % sensors locations
% y_s=y_range(randi(numel(y_range),1,NS))
% INPUT SENSORS
disp(['Place sensors with mouse: ' ' sensors.'])
p_s=ginput(NS)
plot(p_s(:,1),p_s(:,2),...
'LineStyle','none',...
'Marker','o',...
'MarkerSize',4,...
'MarkerEdgeColor',[.5 .5 .5],...
'MarkerFaceColor',[1,0.3921,0])
% round p_s
% pending
% just for plotting purpose, add safety frame on x axis max(x_s(i),x_d) and y axis max(y_s(i),y_d)
% pending
Detect crack locations on metal sheet
The 'field' test consists of a cross formation array of sensors right in front of the target crack. Needless to elaborate that the location of the damage is unknown, it may be quite smaller than the one used in the shown test, and it may convenient to first scatter the sensors as much as possible, for then once pin-pointed the where-about of the damage, then use an tighter array like the close formation shown in the application note, to increase resolution of the damage/target.
The following picture is the key figure of the application note that highlights the successive overlapping of ellipses that eventually reveals the location of the target.
% MEASURE DISTANCES
L=combnk([1:1:NS],2) % L contains all possible combinations of paired sensors
x_c=[];y_c=[];phi=[]
for k=1:1:size(L,1)
pair_1=L(k,1); pair_2=L(k,2)
x_1=p_s(pair_1,1); y_1=p_s(pair_1,2)
x_2=p_s(pair_2,1); y_2=p_s(pair_2,2)
x_c=.5*(x_1+x_2); y_c=.5*(y_1+y_2)
% ell(k).p1=[x_1 y_1] % update ell centres, no needed
% ell(k).p2=[x_2 y_2]
ell(k).centroid=[x_c y_c] % update ell centroid (mid point between ell centres)
c_dist=.5*((x_1-x_2)^2+(y_2-y_2)^2)^.5 % update ell c distance (half distance between ell centres)
% same as pdist([x_1 x_2;y_1 y_2]')
ell(k).c_dist=c_dist
plot([x_1 x_2],[y_1 y_2]) % check
phi0=atan((y_2-y_1)/(x_2-x_1)) % update ell attitude = angle phi
if phi0<0 phi0=2*pi+phi0; end
% phi0/(2*pi)*180 % check
ell(k).phi=phi0
d=[]
for q=1:1:size(S,2)
d=[d ((S(1,q)-x_1)^2+(S(2,q)-y_1)^2)^.5+((S(1,q)-x_2)^2+(S(2,q)-y_2)^2)^.5]
% d is like in bistatic radar: wave travels from p1 [x_1 y_1] to target (scratch) and then to p2 [x_2 y_2]
% this d is in Zima Rucka app note = dt*cg
% d=[d pdist([S(1,q) x_1;S(2,q) y_1]')+pdist([S(1,q) x_2;S(2,q) y_2]')] % same as above line
plot([x_1 S(1,q) x_2],[y_1 S(2,q) y_2],'r-') % check
end
ell(k).d=d
a=d/2 % ell major set of axes for given ell generated for pair of sensors
ell(k).a=a
b=(a.^2.-c_dist.^2).^.5 % ell minor axes
ell(k).b=b
end
plot([x_d-Nd*x_res x_d+Nd*x_res],[y_d-Nd*x_res y_d+Nd*x_res],...
'LineWidth',1,...
'Color','k')
And the closing line, that in MATLAB, there's no need for any other for loop:
% APPROXIMATING TARGET WITH ELLIPSE POINTS
% d_ang_res2=d_ang_res/10 % you need very accurate angle resolution, or you lose
% resolution when target not orbiting, this is, when approaching or departing
% When sensors too far, aliases show up
d_ang_res2=180/(2*pi)*d_ang_res
ang2=[0:d_ang_res2:360-d_ang_res2]
x=[];y=[];
for k=1:1:size(L,1)
x=[x ell(k).centroid(1)+(ell(k).a)'.*cosd(ang2).*cosd(180/(2*pi)*ell(k).phi)-(ell(k).b)'.*sind(ang2).*sind(180/(2*pi)*ell(k).phi)];
y=[y ell(k).centroid(2)+(ell(k).a)'.*cosd(ang2).*sind(180/(2*pi)*ell(k).phi)+(ell(k).b)'.*sind(ang2).*cosd(180/(2*pi)*ell(k).phi)];
end
figure(2);histogram(x);title('histogram(x)');grid on
figure(3);histogram(y);;title('histogram(y)');grid on
as the research article shows with the nice blue graphs of overlapped ellipses, the location of teh damage is eventually revealed. yet depending upon how far are the sensors, the formation used to array them, and the size of the damage, it may be that false targets show up, that may be, or may be not easily discarded.
Increasing the amount of sensors while keeping them in tight formation improves resolution.
For instance repeating with
NS=10
..
and using the following :
nbins1=numel(x_range)
nbins2=numel(y_range)
figure(2);histogram(x,...
'NumBins',nbins1,...
'Normalization','count',...
'FaceColor',[0 0 .5],...
'EdgeColor',[0 0 .5]);
title('histogram(x)');grid on
figure(3);histogram(y,...
'NumBins',nbins2,...
'Normalization','count',...
'FaceColor',[0 0 .5],...
'EdgeColor',[0 0 .5]);
title('histogram(y)');grid on