algoLib/sourceCode/WD_watershed.cpp
2025-11-09 22:54:23 +08:00

506 lines
15 KiB
C++
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

#if 0
#include <iostream>
#include <vector>
#include <queue>
#include <cmath>
#include <cstring>
#include <fstream>
using namespace std;
// 图像尺寸
const int WIDTH = 256;
const int HEIGHT = 256;
// 8邻域方向
const int dx[] = { -1, -1, -1, 0, 0, 1, 1, 1 };
const int dy[] = { -1, 0, 1, -1, 1, -1, 0, 1 };
// 像素结构
struct Pixel {
int x, y; // 坐标
int value; // 像素值(灰度)
Pixel(int x_, int y_, int v_) : x(x_), y(y_), value(v_) {}
// 优先队列需要的比较函数(小顶堆)
bool operator<(const Pixel& other) const {
return value > other.value; // 注意:优先队列默认是大顶堆,这里反转比较
}
};
// 读取PGM格式图像
bool readPGM(const string& filename, int image[HEIGHT][WIDTH]) {
ifstream file(filename, ios::binary);
if (!file) return false;
string magic;
int w, h, max_val;
file >> magic >> w >> h >> max_val;
if (magic != "P5" || w != WIDTH || h != HEIGHT) return false;
// 跳过换行符
file.ignore(1);
// 读取像素值
unsigned char pixel;
for (int i = 0; i < HEIGHT; ++i) {
for (int j = 0; j < WIDTH; ++j) {
file.read((char*)&pixel, 1);
image[i][j] = (int)pixel;
}
}
return true;
}
// 保存分割结果为PGM图像
void writePGM(const string& filename, int labels[HEIGHT][WIDTH], int num_labels) {
ofstream file(filename, ios::binary);
file << "P5\n" << WIDTH << " " << HEIGHT << "\n255\n";
// 将标签映射到0-255范围
int step = 255 / (num_labels + 1);
for (int i = 0; i < HEIGHT; ++i) {
for (int j = 0; j < WIDTH; ++j) {
unsigned char val = (labels[i][j] + 1) * step;
if (labels[i][j] == -1) val = 0; // 分水岭边界
file.write((char*)&val, 1);
}
}
}
// 分水岭算法实现
int watershed(int image[HEIGHT][WIDTH], int labels[HEIGHT][WIDTH]) {
// 初始化标签:-1表示未标记0表示边界
memset(labels, -1, sizeof(int) * HEIGHT * WIDTH);
// 优先队列(按像素值从小到大处理)
priority_queue<Pixel> pq;
// 标记所有局部最小值作为初始盆地
int current_label = 0;
bool is_min;
// 第一步:找到所有局部最小值并标记
for (int i = 0; i < HEIGHT; ++i) {
for (int j = 0; j < WIDTH; ++j) {
is_min = true;
for (int d = 0; d < 8; ++d) {
int ni = i + dx[d];
int nj = j + dy[d];
if (ni >= 0 && ni < HEIGHT && nj >= 0 && nj < WIDTH) {
if (image[ni][nj] < image[i][j]) {
is_min = false;
break;
}
}
}
if (is_min) {
labels[i][j] = current_label++;
pq.push(Pixel(i, j, image[i][j]));
}
}
}
// 第二步:模拟淹没过程
while (!pq.empty()) {
Pixel p = pq.top();
pq.pop();
int x = p.x;
int y = p.y;
for (int d = 0; d < 8; ++d) {
int nx = x + dx[d];
int ny = y + dy[d];
if (nx >= 0 && nx < HEIGHT && ny >= 0 && ny < WIDTH) {
// 如果邻域像素未标记
if (labels[nx][ny] == -1) {
// 标记为当前区域
labels[nx][ny] = labels[x][y];
pq.push(Pixel(nx, ny, image[nx][ny]));
}
// 如果邻域像素已标记且属于不同区域,则标记为边界
else if (labels[nx][ny] != labels[x][y]) {
labels[x][y] = -1; // 当前像素设为边界
}
}
}
}
return current_label;
}
#if 0
算法说明
这个实现遵循分水岭算法的基本原理:
地形表示:将图像的灰度值视为地形高度,灰度值越低表示地势越低
初始盆地:首先识别图像中的所有局部最小值,这些是初始的 "积水盆地"
淹没过程:使用优先队列(最小堆)模拟从低到高的淹没过程,按像素值从小到大处理
区域生长:每个像素被分配给其相邻的最低盆地
边界确定:当两个不同盆地的水相遇时,形成分水岭边界
使用方法
准备一张 256x256 的 PGM 格式灰度图像命名为input.pgm
编译代码g++ watershed.cpp - o watershed
运行程序:. / watershed
结果将保存为output.pgm其中不同区域用不同灰度表示边界为黑色
注意事项
这个实现是基础版本,没有加入预处理步骤(如去噪),可能对复杂图像分割效果不佳
实际应用中通常需要先对图像进行平滑处理,减少过度分割
可以通过调整邻域处理方式4 邻域或 8 邻域)来改变分割效果
对于彩色图像,需要先转换为灰度图再处理
如果需要处理更大尺寸的图像可以修改WIDTH和HEIGHT常量并相应调整内存分配
int main() {
int image[HEIGHT][WIDTH];
int labels[HEIGHT][WIDTH];
// 读取输入图像
if (!readPGM("input.pgm", image)) {
cerr << "无法读取图像文件!" << endl;
return -1;
}
// 执行分水岭算法
int num_labels = watershed(image, labels);
cout << "分割完成,共得到 " << num_labels << " 个区域" << endl;
// 保存结果
writePGM("output.pgm", labels, num_labels);
cout << "结果已保存为 output.pgm" << endl;
return 0;
}
#endif
#else
#include <iostream>
#include <vector>
#include <queue>
#include <algorithm>
#include <fstream>
#include <cmath>
#include "SG_baseDataType.h"
#include "SG_baseAlgo_Export.h"
using namespace std;
#if 0
// 读取PPM格式图像简化版
bool readPPM(const string& filename, Image& img) {
ifstream file(filename, ios::binary);
if (!file) return false;
string magic;
int maxVal;
file >> magic >> img.width >> img.height >> maxVal;
if (magic != "P6") return false; // 仅支持二进制PPM
file.ignore(); // 忽略换行符
vector<unsigned char> data(img.width * img.height * 3);
file.read((char*)data.data(), data.size());
// 转为灰度图Y = 0.299*R + 0.587*G + 0.114*B
img.gray.resize(img.height, vector<int>(img.width));
img.markers.resize(img.height, vector<int>(img.width, 0)); // 初始化标记图为0
for (int i = 0; i < img.height; ++i) {
for (int j = 0; j < img.width; ++j) {
int idx = (i * img.width + j) * 3;
int r = data[idx];
int g = data[idx + 1];
int b = data[idx + 2];
img.gray[i][j] = (int)(0.299 * r + 0.587 * g + 0.114 * b);
}
}
return true;
}
// 保存分割结果为PPM格式
void saveResult(const string& filename, const Image& img) {
ofstream file(filename, ios::binary);
file << "P6\n" << img.width << " " << img.height << "\n255\n";
for (int i = 0; i < img.height; ++i) {
for (int j = 0; j < img.width; ++j) {
if (img.markers[i][j] == -1) { // 分水岭边界(红色)
file.put(255); file.put(0); file.put(0);
}
else { // 区域(根据标记值生成不同颜色)
int color = (img.markers[i][j] * 50) % 256;
file.put(color);
file.put((color + 85) % 256);
file.put((color + 170) % 256);
}
}
}
}
#endif
// 8邻域坐标偏移
const int dx[] = { -1, -1, -1, 0, 0, 1, 1, 1 };
const int dy[] = { -1, 0, 1, -1, 1, -1, 0, 1 };
const int dx4[] = { -1, 0, 0, 1};
const int dy4[] = { 0, -1, 1, 0};
// 计算局部极小值作为初始种子点
void findMinima(SWD_waterShedImage& img, int& markerCount) {
markerCount = 1; //背景对应的ID
// 遍历每个像素判断是否为局部极小值8邻域内最小
for (int i = 1; i < img.height - 1; ++i) {
for (int j = 1; j < img.width - 1; ++j) {
if (1 == img.markers[i][j]) //背景
continue;
bool isMin = true;
#if 0
for (int d = 0; d < 8; ++d) {
int y = i + dy[d];
int x = j + dx[d];
if (img.gray[y][x] < img.gray[i][j]) {
isMin = false;
break;
}
}
#else
for (int dy = -7; dy < 7; dy++)
{
for (int dx = -7; dx < 7; dx++)
{
int y = i + dy;
int x = j + dx;
if ((x >= 0) && (x < img.width) && (y >= 0) && (y < img.height))
{
if (img.gray[y][x] < img.gray[i][j]) {
isMin = false;
break;
}
}
}
}
#endif
if (isMin) {
// 避免重复标记同一区域的极小值
bool hasMarker = false;
for (int d = 0; d < 8; ++d) {
int y = i + dy[d];
int x = j + dx[d];
if (img.markers[y][x] > 1) {
hasMarker = true;
img.markers[y][x] = img.markers[i][j];
break;
}
}
if (!hasMarker) {
img.markers[i][j] = ++markerCount; // 新种子点
}
}
}
}
}
// 分水岭算法主函数
void watershed(SWD_waterShedImage& img)
{
int markerCount;
findMinima(img, markerCount); // 自动标记种子点
// 按灰度值排序所有像素(模拟水位上升)
vector<pair<int, pair<int, int>>> pixels; // (灰度值, (x,y))
for (int i = 0; i < img.height; ++i) {
for (int j = 0; j < img.width; ++j) {
pixels.emplace_back(img.gray[i][j], make_pair(j, i));
}
}
sort(pixels.begin(), pixels.end());
// 遍历排序后的像素,执行注水过程
for (const auto& p : pixels) {
int x = p.second.first;
int y = p.second.second;
if (img.markers[x][y] > 0) continue; // 已标记的种子点
// 收集邻域已标记的区域
vector<int> neighborMarkers;
for (int d = 0; d < 8; ++d) {
int nx = x + dx[d];
int ny = y + dy[d];
if (nx >= 0 && nx < img.height && ny >= 0 && ny < img.width) {
int m = img.markers[nx][ny];
if (m > 0) {
neighborMarkers.push_back(m);
}
}
}
if (neighborMarkers.empty()) {
continue; // 无邻域标记,暂不处理
}
else if (neighborMarkers.size() == 1) {
img.markers[x][y] = neighborMarkers[0]; // 属于同一区域
}
else {
// 多区域交汇,标记为分水岭
img.markers[x][y] = -1;
}
}
}
//模拟注水,从注水口开始,将同水位注满
// inlet: 注水口
// level: 水位
void waterInjection(SWD_waterShedImage& img, SVzNL2DPoint inlet, int level)
{
int py = inlet.y;
int px = inlet.x;
int marker = img.markers[py][px];
std::vector<SVzNL2DPoint> levelPts;
levelPts.push_back(inlet);
int readPtr = 0;
while (readPtr < levelPts.size())
{
SVzNL2DPoint seed = levelPts[readPtr];
readPtr++;
//4邻接
for (int d = 0; d < 4; ++d)
{
int nx = seed.x + dx4[d];
int ny = seed.y + dy4[d];
if (ny >= 0 && ny < img.height && nx >= 0 && nx < img.width)
{
if ((img.gray[ny][nx] == level) && (img.markers[ny][nx] == 0))
{
img.markers[ny][nx] = marker;
SVzNL2DPoint nxt_seed = { nx, ny };
levelPts.push_back(nxt_seed);
}
}
}
}
}
bool _checkMarkerExist(std::vector<int>& MarkerLst, int marker)
{
bool exist = false;
for (int i = 0, i_max = (int)MarkerLst.size(); i < i_max; i++)
{
if (MarkerLst[i] == marker)
{
exist = true;
break;
}
}
return exist;
}
// 根据输入的种子点进行分水岭算法
// watershedSeeds:种子点
// maxLevel:最大水位值
void wd_seedWatershed(SWD_waterShedImage& img, std::vector<SSG_2DValueI>& watershedSeeds, int maxLevel, int startMakerID)
{
int markerCount = startMakerID;
int seedSize = (int)watershedSeeds.size();
for (int i = 0; i < seedSize; i++)
{
int px = watershedSeeds[i].x;
int py = watershedSeeds[i].y;
watershedSeeds[i].value = markerCount;
img.markers[py][px] = markerCount; // 新种子点
int greyValue = img.gray[py][px]; //水位
SVzNL2DPoint inlet = { px, py }; //注水口
waterInjection(img, inlet, greyValue);//注水,将同一水位连接的部分注满
markerCount++;
}
// 按灰度值排序所有像素(模拟水位上升)
std::vector<std::vector<SVzNL2DPoint>> levelPtList;
levelPtList.resize(maxLevel + 1);
for (int y = 0; y < img.height; ++y)
{
for (int x = 0; x < img.width; ++x)
{
int level = img.gray[y][x];
if (level <= maxLevel)
{
SVzNL2DPoint a_pt = { x, y };
levelPtList[level].push_back(a_pt);
}
}
}
while (1)
{
bool allDone = true;
// 遍历排序后的像素,执行注水过程
for (int i = 0; i <= maxLevel; i++)
{
if (levelPtList[i].size() == 0)
continue;
//对同一水位level)执行迭代注水
bool doInjection = true;
while (doInjection)
{
int injectionNum = 0;
std::vector<SVzNL2DPoint> resiPts;
int ptSize = (int)levelPtList[i].size();
for (int pi = 0; pi < ptSize; pi++)
{
int x = levelPtList[i][pi].x;
int y = levelPtList[i][pi].y;
if (img.markers[y][x] > 0) continue; // 已标记的种子点
// 收集邻域已标记的区域
vector<int> neighborMarkers;
for (int d = 0; d < 8; ++d)
{
int nx = x + dx[d];
int ny = y + dy[d];
if (ny >= 0 && ny < img.height && nx >= 0 && nx < img.width)
{
int m = img.markers[ny][nx];
if (m > 0)
{
bool exist = _checkMarkerExist(neighborMarkers, m);
if (false == exist)
neighborMarkers.push_back(m);
}
}
}
if (neighborMarkers.empty())
{
resiPts.push_back(levelPtList[i][pi]);// 无邻域标记,等待迭代
}
else if (neighborMarkers.size() == 1) {
img.markers[y][x] = neighborMarkers[0]; // 属于同一区域
injectionNum++;
}
else {
// 多区域交汇,标记为分水岭
img.markers[y][x] = -1;
injectionNum++;
}
}
levelPtList[i].clear();
levelPtList[i].insert(levelPtList[i].end(), resiPts.begin(), resiPts.end());
if ((injectionNum == 0) || (levelPtList[i].size() == 0))
{
doInjection = false;
if (levelPtList[i].size() > 0)
allDone = false;
}
}
}
if (true == allDone)
break;
}
}
#endif