C++实现D8算法，提取流向和河道点

C++实现D8算法，提取流向和河道点

#include <fstream>
#include <iostream>
#include <string>
#include<vector>
#include<sstream>
#include<ctime>
#include"math.h"
using namespace std;

#define DEM_SIZE 100		//size要小于等于dem数据的长和宽，可以采用vector容器存储数据，复用性更强



void main()
{
	const char *fileName = "F:/yym/zunhua/DEM/dem.txt";//dem数据准备

	//打开流对象
	ifstream ifs;
	//打开方式并判断是否打开成功
	ifs.open(fileName,ios::in);
	if (!ifs.is_open())
	{
		cout << "文件打开失败" << endl;
		return;
	}	
	//读取数据，并存储到数组当中
	string buf;
	int src[DEM_SIZE][DEM_SIZE] = {0};//数据存储到二维数组中
	int count = 0;
	while (getline(ifs,buf))
	{
		//cout << buf << endl;
		
		stringstream ss;
		ss << buf;
		int tmp;
		int num = 0;
		
		
		//从txt文件中（以逗号或者空格分隔）提取dem数值
		while (ss >> tmp)
		{
			//data[num] = tmp;
			src[count][num] = tmp;//dem数据写入二维数组
			num++;
			if (num == DEM_SIZE)//指定每行读取的数据个数（列）
				break;
			if (ss.peek() == ',' || ss.peek() == ' ')
				ss.ignore();
		}
		
		count++;

		//测试
		/*for (int i = 0; i < DEM_SIZE; i++)
		{
			cout << data[i] << " " ;
			
		}
		cout << endl;*/

	}

	//打印dem提取结果

	for (int i = 0; i < 100; i++)
	{
		for (int j = 0; j < 100; j++)
		{
			cout << src[i][j] << " ";

		}
		cout << endl;
	}
	
	//关闭文件
	ifs.close();
	//D8算法

	system("pause");
	system("cls");
	
	//对dem数据进行流向和河道提取
	int Vector[DEM_SIZE][DEM_SIZE] = { 0 }, Result[DEM_SIZE][DEM_SIZE] = { 0 };//初始化流向和河道数组变量
	int row = DEM_SIZE, col = DEM_SIZE;
	double Sqrt2;
	double S = 0, N = 0, E = 0, SE = 0, NE = 0, NW = 0, W = 0, SW = 0;
	Sqrt2 = sqrt(2);//sqrt只支持float和double类型数据，计算平方根
	time_t begin, end;
	begin = clock();
	//d8算法实体，三目运算符
	for (int i = 0; i < row; i++)
	{
		for (int j = 0; j < col; j++)
		{
			S = (i != (row - 1)) ? (src[i][j] - src[i + 1][j]) : -1;
			SE = (i != (row - 1) && j != (col - 1)) ? (src[i][j] - src[i + 1][j + 1]) / sqrt(2) : -1;
			N = (i != 0) ? (src[i][j] - src[i - 1][j]) : -1;
			E = (j != (col - 1)) ? (src[i][j] - src[i][j + 1]) : -1;
			NE = (i != 0 && j != (col - 1)) ? (src[i][j] - src[i - 1][j + 1]) / sqrt(2) : -1;
			NW = (i != 0 && j != 0) ? (src[i][j] - src[i - 1][j - 1]) / sqrt(2) : -1;
			W = (j != 0) ? (src[i][j] - src[i][j - 1]) : -1;
			SW = (i != (row - 1) && j != 0) ? (src[i][j] - src[i + 1][j - 1]) / sqrt(2) : -1;
			//坡降计算
			double M = 0;
			M = (M > S) ? M : S;
			M = (M > SE) ? M : SE;
			M = (M > N) ? M : N;
			M = (M > E) ? M : E;
			M = (M > NE) ? M : NE;
			M = (M > NW) ? M : NW;
			M = (M > W) ? M : W;
			M = (M > SW) ? M : SW;
			//取最大的坡降
			if (M > 0)
			{
				if (M == S)
				{
					Vector[i][j] = 4;
				}
				else if (M == SE)
				{
					Vector[i][j] = 2;
				}
				else if (M == N)
				{
					Vector[i][j] = 64;
				}
				else if (M == E)
				{
					Vector[i][j] = 1;
				}
				else if (M == NE)
				{
					Vector[i][j] = 128;
				}
				else if (M == NW)
				{
					Vector[i][j] = 32;
				}
				else if (M == W)
				{
					Vector[i][j] = 16;
				}
				else if (M == SW)
				{
					Vector[i][j] = 8;
				}
				else
				{
					Vector[i][j] = 0;
				}
			}
		}
	}
	int i1 = 0, j1 = 0;
	//测试流向数据
	/*for(i=0;i<row;i++)
	{
	for(int j=0;j<col;j++)
	cout<<Vector[i][j]<<" ";
	cout<<endl;
	}*/
	for (int i = 0; i < row; i++)
	{
		for (int j = 0; j < col; j++)
		{
			i1 = i;
			j1 = j;
			while (1)// 方向为0 时退出
			{
				int x = Vector[i1][j1];// 取出方向值
				if (Vector[i1][j1] == 0 || i1 >= row || j1 >= col)
					break;
				switch (x)
				{
					// temp=0;
					Vector[i1][j1] |= 0x100000;// 设置该点已被遍历
				case 1:
				{
					Result[i1][j1 + 1]++;
					j1 = j1 + 1;
					break;
				}
				case 2:
				{
					// if(!(V ector[i1][j1]>>16))
					// temp=1;
					Result[i1 + 1][j1 + 1]++;
					i1 = i1 + 1;
					j1 = j1 + 1;
					break;
				}
				case 4:
				{
					//if(!(V ector[i1][j1]>>16))
					// temp=1;
					Result[i1 + 1][j1]++;
					i1 = i1 + 1;
					break;
				}
				case 8:
				{
					// if(!(V ector[i1][j1]>>16))
					// temp=1;
					Result[i1 + 1][j1 - 1]++;
					i1 = i1 + 1;
					j1 = j1 - 1; break;
				}
				case 16:
				{
					// if(!(V ector[i1][j1]>>16))
					// temp=1;
					Result[i1][j1 - 1]++;
					j1 = j1 - 1; break;
				}
				case 32:
				{
					//if(!(V ector[i1][j1]>>16))
					// temp=1;
					Result[i1 - 1][j1 - 1]++;
					i1 = i1 - 1;
					j1 = j1 - 1; break;
				}
				case 64:
				{
					//if(!(V ector[i1][j1]>>16))
					// temp=1;
					Result[i1 - 1][j1]++;
					i1 = i1 - 1; break;
				}
				case 128:
				{
					//if(!(V ector[i1][j1]>>16))
					// temp=1;
					Result[i1 - 1][j1 + 1]++;
					i1 = i1 - 1;
					j1 = j1 + 1; break;
				}
				default:
					break;
				}
				if (Vector[i1][j1] == 0 || i1 >= row || j1 >= col)
					break;
			}
		}
	}
	for (int i = 0; i < row; i++)
	{
		for (int j = 0; j < col; j++)
			cout << Result[i][j] << " ";
		cout << endl;
	}

	//流向结果和汇流结果保存

	//创建流对象
	ofstream ofs;
	//打开方式
	ofs.open("F:/yym/zunhua/DEM/direction.txt",ios::out);

	//写数据
	for (int i = 0; i < 100; i++)
	{
		for (int j = 0; j < 100; j++)
		{
			ofs<<Vector[i][j]<<" " ;
		}
		ofs << endl;
	}
	
	//关闭文件
	ofs.close();
	//创建流对象
	ofstream ofs02;
	//打开方式
	ofs02.open("F:/yym/zunhua/DEM/riverpoint.txt", ios::out);

	//写数据
	for (int i = 0; i < 100; i++)
	{
		for (int j = 0; j < 100; j++)
		{
			ofs02 << Result[i][j] << " ";
		}
		ofs02 << endl;
	}

	//关闭文件
	ofs02.close();
	
	end = clock();
	cout << "runtime: " << double(end - begin) / CLOCKS_PER_SEC << endl;
	getchar();

	//return 0;

	system("pause");
	
}