How to extract the co-ordinates of an ECG signal from a scanned image

Question

I checked the earlier question on "How do I recover the signal from an ECG image"

Currently I am trying to write a program to extract the x and y co-ordinates of an ECG signal embedded within the background grid . I was able to use the image from the earlier post and used some inbuilt functions to crop the image and convert the image to grayscale and then to binary.

This is how I did it:

% Cropping the image and saving the co-ordinates 
img=imcrop(imread('ecg_test.jpg'));
% select the boundries in the pop up window , hit on crop image to save the
% cropped image x,y co-ordinates
img1=rgb2gray(img);  % Converting the cropped image to grayscale
Bw=im2bw(img1,0.8);  % converting the grayscale image to binary image 
imshow(Bw) ![cropped and binary image of the actual ecg_test.jpg][4]

Now I want to clean the image and then extract the co-ordinates of this image. How do I implement the MATLAB part of:

img = ColorConvert[
   ImagePad[
    Import["https://i.stack.imgur.com/500Kg.jpg"], {{0, -20}, {0, 0}}],
    "Grayscale"];

Image[
 Transpose[Function[With[{m = Min[#]},
     Map[Function[{v}, If[v == m, 1, 0]], #]]] /@ 
   Transpose[ImageData[img, "Real"]]
  ]
 ]

Answer 1

I believe this is an optimisation problem. E.g. what is the path from the left hand side of the image to the right hand of the image that is 'most black'.

As such 'rectangular shortest path' method would be ideal. You need to resample the image something like 10x to ensure that the method gets all the peaks properly, but in the end you will get something like this:

Update:

Here is a paper on circular shortest paths, but it is basically the same thing:

http://vision-cdc.csiro.au/changs/doc/sun-pallottino03pr2.pdf

Update 2:

MATLAB script:

IM = imread('500Kg.jpg');
IM = rgb2gray(IM);
IM = double(IM);

% Over sample 8 times in horizontal direction
% Amount to oversample is determined by maximum slope, e.g. if max slope is
% 16, oversample by 16.
IM2 = interp2(IM,linspace(1,size(IM,2),size(IM,2)*8),linspace(1,size(IM,1),size(IM,1))');

% Calculate shortest path using 255 - image intensity as energy function
[spath,senergy,slastindex,im] = shortestPath(255-IM2',255,0,1,100);

% Down sample path to account for original oversampling
p_x = 1:size(IM2,2);
p_x = p_x ./ 8;
p_y = spath(:,2);

% Plot on top of original image
imagesc(IM); colormap gray; hold on; plot(p_x,p_y,'r'); hold off;

Shortest path function:

function [ p, e, l, im ] = shortestPath( m, opt, a, ax, bx  )
%SHORTESTPATH Summary of this function goes here
%   Detailed explanation goes here

%Get size of matrix
[sy, sx] = size(m);
im = zeros(sy,sx);

%Convert to value different from optimal
value = abs(opt - m);

%Preallocate index matrix and energy matrix
lastindex = zeros(sy, sx);
energy = zeros(sy, sx);

%Initialise first rows
energy(1,:) = value(1,:);
lastindex(1,:) = 1:sx;

%Initialise temp matrix
t = zeros(1,3);

%Loop through remaining rows
for row = 2:sy
    for col = ax:bx
        %Get the last energies and current values
        if col == ax 
            t = [inf energy(row-1,col:(col+1))];
            cv = [inf value(row,col:(col+1))];
        elseif col == bx
            t = [energy(row-1,(col-1):col) inf];
            cv = [value(row,(col-1):col) inf];
        else
            t = energy(row-1,(col-1):(col+1));
            cv = value(row,(col-1):(col+1));
        end

        %Add energy for moving
        t = t + [a 0 a];

        %Add energy from difference from optimum
        t = t + cv;

        %Find minimum
        [v,i] = min(t);

        %Save new values
        energy(row,col) = v;
        lastindex(row,col) = col + i - 2;

    end
end

[v,li] = min(energy(sy,ax:bx));
li = ax + li - 1;
p(sy,:) = [sy li];
im(sy,li) = 1;
for row = (sy-1):-1:1
    i = lastindex(row,li);
    p(row,:) = [row,i];
    im(row,i) = 1;
    li = i;
end

e = energy;
l = lastindex;

end