Extract individual curve from a multiple curve plot
7 views (last 30 days)
Show older comments
I have plottedthe attached data using
figure; scatter(data(:,1),data(:,2));
It gives me plot as shown below

As can be seen from the figure that we have 11 curves in this plot. How can I trace each curve and get (X,Y) value and so that I have X and Y values for all 11 curves?
2 Comments
Torsten
on 10 Jan 2025
Can't you modify the way the data are saved ? I don't see a simple pattern in your data matrix that would allow to separate the 10 and 1/10 curves.
Answers (3)
DGM
on 10 Jan 2025
This is probably fragile and overcomplicated, but eh.
% the data
load data.mat
x = data(:,1);
y = data(:,2);
%scatter(x,y,'.');
% try to figure out where everything is supposed to go
nlines = 12;
samplewid = diff([0; find(diff(y) < 0)]); % how many traces are being sampled at each xpos
breakpoints = [0; find(diff(samplewid)~=0); numel(samplewid)]; % the points where samplewid changes
blocklen = diff(breakpoints); % the run length of each block of constant sample length
blockwid = samplewid(breakpoints(1:end-1)+1); % the characteristic number of traces in each block
sampleseq = [{3:12}; {2:12}; {[2:7 9:12]}; {2:12}; {2:11}; {[2 1 3:11]}]; % the sequence of traces in each block
% disassemble everything and put it back together
sorteddata = repmat(struct('x',{[]},'y',{[]}),[nlines 1]); % allocate
baseidx = 0;
for kb = 1:numel(blocklen)
% get index range
bw = blockwid(kb);
bl = blocklen(kb);
blockrange = baseidx + (1:bl*bw);
% extract and reshape [traces xsamples]
thisblockx = reshape(x(blockrange),bw,bl);
thisblocky = reshape(y(blockrange),bw,bl);
% distribute points to the corresponding traces
thisseq = sampleseq{kb};
for kt = 1:blockwid(kb)
sorteddata(thisseq(kt)).x = [sorteddata(thisseq(kt)).x; thisblockx(kt,:).'];
sorteddata(thisseq(kt)).y = [sorteddata(thisseq(kt)).y; thisblocky(kt,:).'];
end
% increment
baseidx = blockrange(end);
end
% plot it?
clf
for k = 1:nlines
plot(sorteddata(k).x,sorteddata(k).y); hold on
end
At first it looked like those lines crossed, but the closer I look, it really doesn't. I have no idea of what created this data, so maybe that local offset is an artifact of how they were captured.
6 Comments
DGM
on 11 Jan 2025
I don't know of any simple/robust way to automatically sort the lines and find their implied intersections. If it's possible to get the sorted (but non-intersecting lines as in the first example, then maybe you can find places where adjacent lines are maximally-close and have high curvature, but that's probably easily defeated. The problem in that assumption is that we're relying on being able to automatically sort them to begin with.
I obtained sampleseq just by visual observation. The breakpoints used for swapping were found using a datatip on the sorted data in a plot. The last example with the interpolated data makes the process easier since the xdata is the same between both curves.

% in the second example, line 8 and 9 are offset by 1 (19 vs 20)
% in the interpolated example, both 8 and 9 use the same index (20)
find(abs(sorteddata(8).x - 561)<1) % find the index for this datatip
ans =
20
Like I said, I feel like what I wrote is fragile and probably not the ideal way to do it, but I don't know what is. I was assuming someone would have had a better idea by now.
DGM
on 11 Jan 2025
Edited: DGM
on 11 Jan 2025
I was trying to put together a crude line-following routine, but I haven't had much luck getting it to produce usable results except on some of the simpler, more isolated lines. It only handles one of the two crossings, and dealing with gaps is problematic.
Image Analyst
on 11 Jan 2025
I think the best solution is not to try to fix a mess of data but to try to prevent the mess in the first place. Was the .mat file created in MATLAB or by some other program? I'd see if there was some way to write out the .mat file in a more orderly way by writing oukt separate variables to it. Or even writing out multiple mat files with one curve in each.
0 Comments
DGM
on 12 Jan 2025
Edited: DGM
on 12 Jan 2025
This is the garbage I put together. It could stand to be cleaned up, but it still has problems. It seems to sort data and fill small gaps okay, but the crossed traces (8 & 9) have too much curvature for it to fix. In the end, this requires less manual work, but swapping traces still needs to be done manually.
% the data
load data.mat
% for demonstration, make some new gaps
data(100,:) = [];
data(181,:) = [];
% parameters
seglen = 2; % max length of trailing segment over which to average slope
angletol = 15; % only look within +-angletol of current trajectory
maxbridgesteps = 2; % number of mean x-steps to bridge gaps
bridgeangletol = 15; % same as before, but for the bridging routine
% normalize the data columns temporarily
datan = data; % work on a copy
normval = [imrange(datan(:,1)); imrange(datan(:,2))];
datan(:,1) = rescalelin(datan(:,1),normval(1,:),[0 1]);
datan(:,2) = rescalelin(datan(:,2),normval(2,:),[0 1]);
% plot
scatter(datan(:,1),datan(:,2),'.'); hold on
% initial data sorting
% this will still have issues with gaps or crossings
allx = unique(datan(:,1));
kt = 1;
while ~isempty(datan)
xci = find(allx == datan(1,1)); % find the position of the smallest remaining x-value
thistrace = datan(1,:); % pick a new leftmost point to start a new trace
datan(1,:) = []; % remove it from the set
firstpoint = true; % initialize
followingtrace = true; % initialize
while followingtrace
nextx = allx(min(xci+1,numel(allx))); % get the next x-value
Pn = datan(any(datan(:,1) == nextx,2),:); % these are candidates in the next stripes
if firstpoint
% the next point is the closest point
D = sqrt(sum((thistrace(end,:) - Pn).^2,2)); % find distance to all candidate points
[~,idx] = min(D); % minimize
thistrace = [thistrace; Pn(idx,:)]; %#ok<*AGROW>
firstpoint = false;
% plot
hp = plot(thistrace(:,1),thistrace(:,2));
else
% the next point is the closest point
% in the direction of the current trajectory
%d = thistrace(end,:) - thistrace(end-1,:);
endseg = thistrace(max(end-seglen,1):end,:);
d = mean(diff(endseg,1,1),1);
th = atan2d(d(2),d(1)); % angle between last two points
dn = Pn - thistrace(end,:);
thn = atan2d(dn(:,2),dn(:,1)); % angle between last point and all candidates
% get rid of candidates which are off-trajectory
goodangle = abs(th - thn) <= angletol;
Pn = Pn(goodangle,:);
% find the closest candidate of those left
if isempty(Pn)
% terminate this trace
followingtrace = false;
continue;
else
D = sqrt(sum((thistrace(end,:) - Pn).^2,2)); % find distance to all candidate points
[~,idx] = min(D); % minimize
thistrace = [thistrace; Pn(idx,:)];
end
% update the plot
hp.XData = thistrace(:,1);
hp.YData = thistrace(:,2);
end
% prune this point from the data
datan(all(datan == thistrace(end,:),2),:) = [];
% increment
xci = find(allx == thistrace(end,1));
end
% store the completed trace
sorteddata{kt} = thistrace;
kt = kt + 1;
end
% try to stitch together any obvious gaps
traceindex = 1;
while traceindex <= numel(sorteddata)
% setup
ntraces = numel(sorteddata);
startpoints = zeros(ntraces,2);
startpointidx = 1:ntraces;
% get current end segment and all candidate start points
thisendseg = sorteddata{traceindex}(end-1:end,:);
for k = 1:ntraces
startpoints(k,:) = sorteddata{k}(1,:);
end
% look for traces which begin nearby
% in the direction of the end trajectory
d = diff(thisendseg,1,1);
th = atan2d(d(2),d(1)); % angle between last two points
dn = startpoints - thisendseg(end,:);
thn = atan2d(dn(:,2),dn(:,1)); % angle between last point and all start points
% get rid of candidates which are off-trajectory
goodangle = abs(th - thn) <= bridgeangletol;
startpoints = startpoints(goodangle,:);
startpointidx = startpointidx(goodangle);
% find the closest of the points which are within allowable x-distance
xD = abs(thisendseg(end,1) - startpoints(:,1));
gooddistance = xD <= (maxbridgesteps + 1)*mean(diff(allx));
xD = xD(gooddistance);
[~,matchidx] = min(xD); % minimize
matchidx = startpointidx(matchidx);
if isempty(matchidx)
% start building a new composite trace
traceindex = traceindex + 1;
else
% append the selected trace and remove it from the set
sorteddata{traceindex} = [sorteddata{traceindex}; sorteddata{matchidx}];
sorteddata(matchidx) = [];
end
end
ntraces = numel(sorteddata);
% denormalize the sorted data
for k = 1:ntraces
sorteddata{k}(:,1) = rescalelin(sorteddata{k}(:,1),[0 1],normval(1,:));
sorteddata{k}(:,2) = rescalelin(sorteddata{k}(:,2),[0 1],normval(2,:));
end
% plot final sorted & bridged data
figure(2)
clf
%scatter(data(:,1),data(:,2),'.')
hold on
for k = 1:ntraces
plot(sorteddata{k}(:,1),sorteddata{k}(:,2))
end
There are a couple utility functions attached.
If you want, you can manually fix the failed crossing like so:
%% manually fix crossing
radtol = 0.1; % normalized x-distance local to click
% click on the figure where you think an intersection should be
[x y] = ginput(1);
% find the two traces which are closest to the click
mindist2click = zeros(ntraces,1);
for k = 1:ntraces
D = sqrt(sum((sorteddata{k} - [x y]).^2,2));
mindist2click(k) = min(D);
end
[~,tidx] = mink(mindist2click,2);
% sample the traces local to the click
% this avoids issues caused by traces which may be close on their tails
xtol = radtol*diff(normval(1,:));
localsamples = cell(2,1);
for k = 1:2
thissample = sorteddata{tidx(k)};
mk = abs(thissample(:,1) - x) <= xtol;
localsamples{k} = thissample(mk,:);
end
% assume that we can discard unique x-values which may be caused by an offset
% the pointlists should have the same length now
commonx = intersect(localsamples{1}(:,1),localsamples{2}(:,1));
mk = ismember(localsamples{1}(:,1),commonx);
localsamples{1} = localsamples{1}(mk,:);
mk = ismember(localsamples{2}(:,1),commonx);
localsamples{2} = localsamples{2}(mk,:);
% minimize the y-distance
yD = abs(localsamples{1}(:,2) - localsamples{2}(:,2));
[~,idx] = min(yD);
% find the corresponding points in the original-length data
% the data is not interpolated, so they might not be the same index
for k = 1:2
bp(k) = find(sorteddata{tidx(k)} == commonx(idx),1,'first');
end
% swap the curves at the nearest point
temp1 = sorteddata{tidx(1)}(bp(1):end,:);
temp2 = sorteddata{tidx(2)}(bp(2):end,:);
sorteddata{tidx(1)} = [sorteddata{tidx(1)}(1:bp(1)-1,:); temp2];
sorteddata{tidx(2)} = [sorteddata{tidx(2)}(1:bp(2)-1,:); temp1];
% replot final sorted & bridged data
figure(2)
clf
%scatter(data(:,1),data(:,2),'.')
hold on
for k = 1:ntraces
plot(sorteddata{k}(:,1),sorteddata{k}(:,2))
end
That's all I've got. These bits should probably all be rolled up into functions of some sort so that they're more practical to use.
2 Comments
DGM
on 13 Jan 2025
That's what this code does. The problem is that the intersections are not distinct. The local slope and curvature favors not making the crossing. Trying to make the averaging window larger will start creating connections between traces long before the desired crossing is made. Is there a way to improve it? Probably.
See Also
Categories
Find more on Time Series Events in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!