如何使用MATLAB计算图像中曲线下的面积？

Question

alt text http://internationalpropertiesregistry.com/Server/showFile.php?file=%2FUpload%2Fstatistics.gifc49ca28823a561a41d09ef9adbb5e0c5.gif

x轴的单位是小时(h),总共有24小时.

y轴的单位是数百万(m).

如何以图像为单位计算图像中红色曲线下面积m*h？

重要更新

只有图像可用(而不是数据),我想以编程方式计算区域.

Answer 1

这是一个有趣的解决方案:).顺便说一句,它使用bwfill(类似于imfill)需要一些用户交互.

码

%# Constants
gray_value_curve = 2;
gray_value_box = 3;
area_box_in_units = 10;

%# Read the image
I = imread('C:\p23\graph.gif');
%# Find the area of a unit block
figure(1);
imshow(I,[]);
[BS sq_elem] = bwfill;   
imshow(BS,[]);
%# Get the dimensions to make the estimate more accurate
X = zeros(size(BS));
X(sq_elem) = 1;
s = regionprops(X,'Area','BoundingBox');
block_area = s.Area + 2*(s.BoundingBox(3)-1) + 2*(s.BoundingBox(4)-1) + 4;

%#Find the area under the curve
I( ~(I == gray_value_curve | I == gray_value_box) ) = 0;
figure(2);
imshow(I,[]);
[BA area_curve_elem] = bwfill;
imshow(BA,[]);
%# Area under the curve
curve_area = numel(area_curve_elem);

%# Display the area in the required units
area = area_box_in_units*curve_area/block_area;
disp(area);

产量

113.5259

图1 图2

Answer 2

创建全自动解决方案的困难在于，它需要您将有关要处理的输入图像的某些假设硬编码到解决方案中。如果这些假设不适用于您可能遇到的所有潜在图像，则全自动解决方案将不会给出值得信赖的结果，并且尝试扩展全自动解决方案以处理所有可能的输入可能会导致其膨胀难以理解且复杂的代码。

当对输入图像特征的可变性有疑问时，像Jacob这样带有一些用户交互的解决方案通常是最好的。如果您可以确定输入图像的特征遵循一组严格的规则，那么可以考虑自动化解决方案。

作为示例，下面是我编写的一些自动代码，用于估算图表中红色曲线下的面积。由于我使用上图作为指导，因此必须满足许多条件才能使其发挥作用：

• 绘制线的红色像素必须在图像中唯一地描述为包含等于 0 的绿色和蓝色分量以及等于 1 的红色分量。
• 网格线的绿色像素必须在图像中唯一地描述为包含小于 1 的红色和蓝色分量以及等于 1 的绿色分量。
• 轴线的蓝色像素必须在图像中唯一地描述为包含等于 0 的红色和绿色分量以及等于 1 的蓝色分量。
• 网格线和轴线必须始终在水平或垂直方向上精确对齐。
• 网格线的长度必须远远超过图像宽度和高度的一半。
• x 轴必须是图像中最长的水平蓝线。
• 网格线必须始终为 1 像素粗。

根据输入图像的上述条件，可以使用以下代码来近似红色曲线下的面积，而无需用户输入：

[img,map] = imread('original_chart.gif');  %# Read the indexed image
[r,c] = size(img);                         %# Get the image size

redIndex = find((map(:,1) == 1) & ...    %# Find the red index value
                (map(:,2) == 0) & ...
                (map(:,3) == 0))-1;
greenIndex = find((map(:,1) < 1) & ...   %# Find the green index value
                  (map(:,2) == 1) & ...
                  (map(:,3) < 1))-1;
blueIndex = find((map(:,1) == 0) & ...   %# Find the blue index value
                 (map(:,2) == 0) & ...
                 (map(:,3) == 1))-1;

redLine = (img == redIndex);      %# A binary image to locate the red line
greenLine = (img == greenIndex);  %# A binary image to locate the grid lines
blueLine = (img == blueIndex);    %# A binary image to locate the axes lines

w = mean(diff(find(sum(greenLine,1) > r/2)));  %# Compute unit square width
h = mean(diff(find(sum(greenLine,2) > c/2)));  %# Compute unit square height
squareArea = w*h;                              %# Compute unit square area

[maxValue,maxIndex] = max(redLine);          %# Find top edge of red line
x = find(maxValue > 0);                      %# Find x coordinates of red line
y = maxIndex(maxValue > 0);                  %# Find y coordinates of red line
[maxValue,maxIndex] = max(sum(blueLine,2));  %# Find row index of x axis
y = maxIndex-y;                              %# Zero the y coordinate
totalArea = trapz(x,y)/squareArea;           %# Compute the area under the curve

给出以下结果：

squareArea = 460.6 square pixels
totalArea = 169.35 m*h

解释：

我将详细说明计算中涉及的步骤w：

1. greenLine使用函数SUM对二值图像沿每列求和，给出一个1-by-c向量，其中每个元素是图像每列中有多少网格线像素的计数。
2. 该向量中大于r/2（图像行数一半）的元素表示图像中包含垂直网格线的列。使用函数FIND可以找到这些列的索引。
3. 使用函数DIFF可以找到这些列索引之间的成对差异。这给出了一个包含网格线之间的空间宽度（以像素为单位）的向量。
4. 最后，函数MEAN用于计算图像中所有网格线之间的空间的平均宽度。

计算 时h，唯一的区别是沿每一行求和r/2并替换为c/2（图像中列数的一半）。