commit 1093d8febc1d6b679083d2f75d280fd6af7c6d99 Author: 陈增辉 <3045316072@qq.com> Date: Fri Nov 15 09:38:46 2024 +0800 代码初始化 diff --git a/BaseConstVariable.h b/BaseConstVariable.h new file mode 100644 index 0000000..a4d7878 --- /dev/null +++ b/BaseConstVariable.h @@ -0,0 +1,129 @@ +#pragma once + +#ifndef BASECONSTVARIABLE_H +#define BASECONSTVARIABLE_H +#define EIGEN_USE_MKL_ALL +//#define EIGEN_NO_DEBUG + + +#define EIGEN_USE_BLAS +#define EIGEN_USE_LAPACK +#define EIGEN_VECTORIZE_SSE2 + +//#define DEBUGSHOWDIALOG + +//#include +#include +#include +#include + +#define PI_180 180/3.141592653589793238462643383279 +#define T180_PI 3.141592653589793238462643383279/180 +#define LIGHTSPEED 299792458 +#define PRECISIONTOLERANCE 1e-9 +#define Radians2Degrees(Radians) Radians*PI_180 +#define Degrees2Radians(Degrees) Degrees*T180_PI + +#define MATPLOTDRAWIMAGE + + +#define earthRoute 0.000072292115 +#define Memory1GB 1073741824 +#define Memory1MB 1048576 +#define Memory1KB 1024 + +const std::complex imagI(0, 1); + +const double PI = 3.141592653589793238462643383279; +const double epsilon = 0.000000000000001; +const double pi = 3.14159265358979323846; +const double d2r = pi / 180; +const double r2d = 180 / pi; + +const double a = 6378137.0; //椭球长半轴 +const double ae = 6378137.0; //椭球长半轴 +const double ee = 0.0818191910428;// 第一偏心率 +const double f_inverse = 298.257223563; //扁率倒数 +const double b = a - a / f_inverse; +const double eSquare = (a * a - b * b) / (a * a); +const double e = sqrt(eSquare); +const double earth_Re = 6378136.49; +const double earth_Rp = (1 - 1 / f_inverse) * earth_Re; +const double earth_We = 0.000072292115; + + +/*********************************************** 基础枚举 ********************************************************************/ + +enum POLARTYPEENUM { // 极化类型 + POLARHH, + POLARHV, + POLARVH, + POLARVV +}; + + +/*********************************************** 基础结构体区域 ********************************************************************/ + +/// +/// 地理坐标点 +/// +struct Landpoint // 点 SAR影像的像素坐标； +{ + double lon; // 经度x lon pixel_col + double lat; // 纬度y lat pixel_row + double ati; // 高程z ati pixel_time +}; + +struct Point2 { + double x = 0; + double y = 0; +}; + + +struct Point3 { + double x = 0; + double y = 0; + double z = 0; +}; + +/*********************************************** FEKO仿真参数 ********************************************************************/ + + +struct SatelliteAntPos { + double time; + double Px; + double Py; + double Pz; + double Vx; + double Vy; + double Vz; + double AntDirectX; + double AntDirectY; + double AntDirectZ; + double AVx; + double AVy; + double AVz; + double lon; + double lat; + double ati; +}; + + +/*********************************************** 指针回收区域 ********************************************************************/ + +inline void delArrPtr(void* p) +{ + delete[] p; + p = nullptr; +} + +inline void delPointer(void* p) +{ + delete p; + p = nullptr; +} + + + + +#endif \ No newline at end of file diff --git a/BaseTool.cpp b/BaseTool.cpp new file mode 100644 index 0000000..602747a --- /dev/null +++ b/BaseTool.cpp @@ -0,0 +1,542 @@ +#include "stdafx.h" +#define EIGEN_USE_BLAS +#define EIGEN_VECTORIZE_SSE4_2 +//#include + +#include +#include +#include +#include +////#include +#include +#include +#include +#include < io.h > +#include < stdio.h > +#include < stdlib.h > +#include +#include +#include +//#include +#include //#include "ogrsf_frmts.h" + +#include +#include +#include "GeoOperator.h" + +#include +#include +#include "baseTool.h" + +#include // 多项式拟合 +#include /* 提供了 gammaq 函数 */ +#include /* 提供了向量结构*/ +#include +#include + + + + + +QString getCurrentTimeString() { + struct tm ConversionTime; + std::time_t t = std::time(NULL); + char mbstr[100]; + _localtime64_s(&ConversionTime, &t); + std::strftime(mbstr, sizeof(mbstr), "%Y-%m-%d %H:%M:%S", &ConversionTime); + QString strTime = mbstr; + return strTime; +} + +QString getCurrentShortTimeString() { + struct tm ConversionTime; + std::time_t t = std::time(NULL); + char mbstr[100]; + _localtime64_s(&ConversionTime, &t); + std::strftime(mbstr, sizeof(mbstr), "%Y-%m-%d %H:%M:%S", &ConversionTime); + QString strTime = mbstr; + return strTime; +} + +std::vector splitString(const QString& str, char delimiter) +{ + QStringList tokens = str.split(delimiter); + return convertQStringListToStdVector(tokens); +} + + +std::complex Cubic_Convolution_interpolation(double u, double v, Eigen::MatrixX> img) +{ + if (img.rows() != 4 || img.cols() != 4) { + throw std::exception("the size of img's block is not right"); + } + // 斤拷锟斤拷模锟斤拷 + Eigen::MatrixX> wrc(1, 4);// 使锟斤拷 complex 斤拷锟斤拷要原斤拷为锟剿伙拷取值 + Eigen::MatrixX> wcr(4, 1);// + for (int i = 0; i < 4; i++) { + wrc(0, i) = Cubic_kernel_weight(u + 1 - i); // u+1,u,u-1,u-2 + wcr(i, 0) = Cubic_kernel_weight(v + 1 - i); + } + + Eigen::MatrixX> interValue = wrc * img * wcr; + return interValue(0, 0); +} + +std::complex Cubic_kernel_weight(double s) +{ + s = abs(s); + if (s <= 1) { + return std::complex(1.5 * s * s * s - 2.5 * s * s + 1, 0); + } + else if (s <= 2) { + return std::complex(-0.5 * s * s * s + 2.5 * s * s - 4 * s + 2, 0); + } + else { + return std::complex(0, 0); + } +} + +/// +/// p11 p12 -- x +/// p0(u,v) +/// p21 p22 +/// | +/// y +/// p11(0,0) +/// p21(0,1) +/// P12(1,0) +/// p22(1,1) +/// +/// x,y,z +/// x,y,z +/// x,y,z +/// x,y,z +/// x,y,z +/// +double Bilinear_interpolation(Landpoint p0, Landpoint p11, Landpoint p21, Landpoint p12, Landpoint p22) +{ + + return p11.ati * (1 - p0.lon) * (1 - p0.lat) + + p12.ati * p0.lon * (1 - p0.lat) + + p21.ati * (1 - p0.lon) * p0.lat + + p22.ati * p0.lon * p0.lat; +} + + + +bool onSegment(Point3 Pi, Point3 Pj, Point3 Q) +{ + if ((Q.x - Pi.x) * (Pj.y - Pi.y) == (Pj.x - Pi.x) * (Q.y - Pi.y) //叉乘 + //保证Q点坐标在pi,pj之间 + && std::min(Pi.x, Pj.x) <= Q.x && Q.x <= std::max(Pi.x, Pj.x) + && std::min(Pi.y, Pj.y) <= Q.y && Q.y <= std::max(Pi.y, Pj.y)) + return true; + else + return false; +} + +Point3 invBilinear(Point3 p, Point3 a, Point3 b, Point3 c, Point3 d) +{ + Point3 res; + + Point3 e = b - a; + Point3 f = d - a; + Point3 g = a - b + c - d; + Point3 h = p - a; + + double k2 = cross2d(g, f); + double k1 = cross2d(e, f) + cross2d(h, g); + double k0 = cross2d(h, e); + double u, v; + // if edges are parallel, this is a linear equation + if (abs(k2) < 0.001) + { + v = -k0 / k1; + u = (h.x - f.x * v) / (e.x + g.x * v); + p.z = a.z + (b.z - a.z) * u + (d.z - a.z) * v + (a.z - b.z + c.z - d.z) * u * v; + return p; + } + // otherwise, it's a quadratic + else + { + float w = k1 * k1 - 4.0 * k0 * k2; + if (w < 0.0) { + // 可能在边界上 + if (onSegment(a, b, p)) { + Point3 tt = b - a; + Point3 ttpa = p - a; + double scater = ttpa / tt; + if (scater < 0 || scater>1) { return { -9999,-9999,-9999 }; } + p.z = a.z + scater * tt.z; + return p; + } + else if (onSegment(b, c, p)) { + Point3 tt = c - b; + Point3 ttpa = p - b; + double scater = ttpa / tt; + if (scater < 0 || scater>1) { return { -9999,-9999,-9999 }; } + p.z = b.z + scater * tt.z; + return p; + } + else if (onSegment(c, d, p)) { + Point3 tt = d - c; + Point3 ttpa = p - c; + double scater = ttpa / tt; + if (scater < 0 || scater>1) { return { -9999,-9999,-9999 }; } + p.z = c.z + scater * tt.z; + return p; + } + else if (onSegment(d, a, p)) { + Point3 tt = a - d; + Point3 ttpa = p - d; + double scater = ttpa / tt; + if (scater < 0 || scater>1) { return { -9999,-9999,-9999 }; } + p.z = d.z + scater * tt.z; + return p; + } + + return { -9999,-9999,-9999 }; + } + else { + w = sqrt(w); + + float ik2 = 0.5 / k2; + float v = (-k1 - w) * ik2; + float u = (h.x - f.x * v) / (e.x + g.x * v); + + if (u < 0.0 || u>1.0 || v < 0.0 || v>1.0) + { + v = (-k1 + w) * ik2; + u = (h.x - f.x * v) / (e.x + g.x * v); + } + p.z = a.z + (b.z - a.z) * u + (d.z - a.z) * v + (a.z - b.z + c.z - d.z) * u * v; + return p; + } + } + p.z = a.z + (b.z - a.z) * u + (d.z - a.z) * v + (a.z - b.z + c.z - d.z) * u * v; + return p; +} + +double sind(double degree) +{ + return sin(degree * d2r); +} + +double cosd(double d) +{ + return cos(d * d2r); +} + + +std::string Convert(float Num) +{ + std::ostringstream oss; + oss << Num; + std::string str(oss.str()); + return str; +} + +QString JoinPath(const QString& path, const QString& filename) +{ + QDir dir(path); + + // Ensure the path ends with the appropriate separator + if (!QDir::isAbsolutePath(path)) + dir.makeAbsolute(); + + return dir.filePath(filename); +} +std::vector convertQStringListToStdVector(const QStringList& qStringList) +{ + std::vector stdVector; + + for (const QString& str : qStringList) { + stdVector.push_back(str); + } + + return stdVector; +}; + + +ErrorCode polyfit(const double* x, const double* y, int xyLength, int poly_n, std::vector& out_factor, double& out_chisq) { + /* + * x：自变量，视差 + * y：因变量，距离 + * xyLength: x、y长度 + * poly_n：拟合的阶次 + * out_factor：拟合的系数结果，从0阶到poly_n阶的系数 + * out_chisq：拟合曲线与数据点的优值函数最小值 ,χ2 检验 + */ + + gsl_matrix* XX = gsl_matrix_alloc(xyLength, poly_n + 1); + gsl_vector* c = gsl_vector_alloc(poly_n + 1); + gsl_matrix* cov = gsl_matrix_alloc(poly_n + 1, poly_n + 1); + gsl_vector* vY = gsl_vector_alloc(xyLength); + + for (size_t i = 0; i < xyLength; i++) + { + gsl_matrix_set(XX, i, 0, 1.0); + gsl_vector_set(vY, i, y[i]); + for (int j = 1; j <= poly_n; j++) + { + gsl_matrix_set(XX, i, j, pow(x[i], j)); + } + } + gsl_multifit_linear_workspace* workspace = gsl_multifit_linear_alloc(xyLength, poly_n + 1); + int r = gsl_multifit_linear(XX, vY, c, cov, &out_chisq, workspace); + gsl_multifit_linear_free(workspace); + out_factor.resize(c->size, 0); + for (size_t i = 0; i < c->size; ++i) + { + out_factor[i] = gsl_vector_get(c, i); + } + + gsl_vector_free(vY); + gsl_matrix_free(XX); + gsl_matrix_free(cov); + gsl_vector_free(c); + + return GSLState2ErrorCode(r); +} + +Point3 crossProduct(const Point3& a, const Point3& b) +{ + return Point3{ + a.y * b.z - a.z * b.y, // C_x + a.z * b.x - a.x * b.z, // C_y + a.x * b.y - a.y * b.x // C_z + }; +} + + + + +// 计算绕任意轴旋转的旋转矩阵 +Eigen::Matrix3d rotationMatrix(const Eigen::Vector3d& axis, double angle) { + // 确保旋转轴是单位向量 + Eigen::Vector3d u = axis.normalized(); + + double cos_theta = cos(angle); + double sin_theta = sin(angle); + + // 计算旋转矩阵 + Eigen::Matrix3d R; + R << + cos_theta + u.x() * u.x() * (1 - cos_theta), u.x()* u.y()* (1 - cos_theta) - u.z() * sin_theta, u.x()* u.z()* (1 - cos_theta) + u.y() * sin_theta, + u.y()* u.x()* (1 - cos_theta) + u.z() * sin_theta, cos_theta + u.y() * u.y() * (1 - cos_theta), u.y()* u.z()* (1 - cos_theta) - u.x() * sin_theta, + u.z()* u.x()* (1 - cos_theta) - u.y() * sin_theta, u.z()* u.y()* (1 - cos_theta) + u.x() * sin_theta, cos_theta + u.z() * u.z() * (1 - cos_theta); + + return R; +} + + + + + +long double convertToMilliseconds(const std::string& dateTimeStr) { + // 定义日期时间字符串 + std::string dateTimeString = dateTimeStr; + // 使用 Boost.Date_Time 解析日期时间字符串 + boost::posix_time::ptime dateTime = boost::posix_time::time_from_string(dateTimeString); + + // 将 ptime 转换为自 epoch (1970年1月1日) 以来的毫秒数 + boost::posix_time::ptime epoch(boost::gregorian::date(1970, 1, 1)); + boost::posix_time::time_duration duration = dateTime - epoch; + long double milliseconds = duration.total_milliseconds() / 1000.0; + return milliseconds; + +} + +long FindValueInStdVector(std::vector& list, double& findv) +{ + if (list.size() == 0) { + return -1; + } + else if (list.size() == 1) { + if (abs(list[0] - findv) < 1e-9) { + return 0; + } + else { + return -1; + } + } + else {} + + if (abs(list[0] - findv) < 1e-9) { + return 0; + } + else if (abs(list[list.size() - 1] - findv) < 1e-9) { + return list.size() - 1; + } + else if (list[0] > findv) { + return -1; + } + else if (list[list.size() - 1] < findv) { + return -1; + } + else {} + + + + long low = 0; + long high = list.size() - 1; + while (low <= high) { + long mid = (low + high) / 2; + if (abs(list[mid] - findv) < 1e-9) { + return mid; + } + else if (list[mid] < findv) { + low = mid + 1; + } + else if (list[mid] > findv) { + high = mid - 1; + } + else {} + } + return -1; +} + +long InsertValueInStdVector(std::vector& list, double insertValue, bool repeatValueInsert) +{ + if (list.size() == 0) { + list.insert(list.begin(), insertValue); + return 0; + } + else if (list.size() == 1) { + if (std::abs(list[0] - insertValue) < PRECISIONTOLERANCE) { + if (repeatValueInsert) { + list.push_back(insertValue); + return 1; + } + else { + return -1; + } + } + else if (insertValue > list[0]) { + list.push_back(insertValue); + return 1; + } + else if (insertValue < list[0]) { + list.push_back(insertValue); + return 0; + } + } + else { + long low = 0; + long high = list.size() - 1; + long mid = -1; + // 处理边界 + if (list[high] <= insertValue) + { + if (std::abs(list[high] - insertValue) < PRECISIONTOLERANCE) { + if (repeatValueInsert) { + list.push_back(insertValue); + } + else {} + } + else { + list.push_back(insertValue); + } + return list.size() - 1; + } + + if (list[low] >= insertValue) { + + if (std::abs(list[low] - insertValue) < PRECISIONTOLERANCE) { + if (repeatValueInsert) { + list.insert(list.begin(), insertValue); + } + else {} + } + else { + list.insert(list.begin(), insertValue); + } + return 0; + } + // low 1) { + mid = (low + high) / 2; + if (std::abs(list[mid] - insertValue) < PRECISIONTOLERANCE) { + if (repeatValueInsert) { + list.insert(list.begin() + mid, insertValue); + } + return mid; + } + else if (insertValue < list[mid]) { + high = mid; + } + else if (list[mid] < insertValue) { + low = mid; + } + + } + if (list[low] <= insertValue && insertValue <= list[high]) + { + + if (std::abs(list[high] - insertValue) < PRECISIONTOLERANCE) { + if (repeatValueInsert) { + list.insert(list.begin() + high, insertValue); + } + else {} + } + else { + list.insert(list.begin() + high, insertValue); + } + return low; + } + else { + return -1; + } + } + return -1; +} + +long FindValueInStdVectorLast(std::vector& list, double& insertValue) +{ + if (list.size() == 0 || list.size() == 1) { + return -1; + } + else { + long low = 0; + long high = list.size() - 1; + long mid = -1; + // 处理边界 + if (list[high] <= insertValue) + { + return -1; + } + else if (std::abs(list[high] - insertValue) < PRECISIONTOLERANCE) + { + return high; + } + else if (list[low] > insertValue) { + return -1; + } + else if (std::abs(list[low] - insertValue) < PRECISIONTOLERANCE) { + return low; + } + else {} + // low 1) { + mid = (low + high) / 2; + if (insertValue < list[mid]) { + high = mid; + } + else if (list[mid] < insertValue) { + low = mid; + } + else {// insertValue==list[mid] , list[mid]===insertValue + return mid;// + } + } + if (list[low] <= insertValue && insertValue <= list[high]) + { + return low; + } + else { + return -1; + } + } + return -1; +} + diff --git a/BaseTool.h b/BaseTool.h new file mode 100644 index 0000000..86153c5 --- /dev/null +++ b/BaseTool.h @@ -0,0 +1,112 @@ +#pragma once +#ifndef BASETOOL_H +#define BASETOOL_H +#include "BaseConstVariable.h" + +/// +/// 基本类、基本函数 +/// + + +// //#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include "GeoOperator.h" +#include +#include +#include +#include + + +#include "LogInfoCls.h" + +///////////////////////////////////// 运行时间打印 +///////////////////////////////////////////////////////////// + + +QString getCurrentTimeString(); +QString getCurrentShortTimeString(); + +std::vector splitString(const QString& str, char delimiter); +std::vector convertQStringListToStdVector(const QStringList& qStringList); +/////////////////////////////// 基本图像类 结束 +///////////////////////////////////////////////////////////// + +std::string Convert(float Num); +QString JoinPath(const QString& path, const QString& filename); + +////////////////////////////// 坐标部分基本方法 ////////////////////////////////////////// + +////////////////////////////// 坐标部分基本方法 ////////////////////////////////////////// + +////////////////////////////// 插值 //////////////////////////////////////////// + +std::complex Cubic_Convolution_interpolation(double u, double v, + Eigen::MatrixX> img); + +std::complex Cubic_kernel_weight(double s); + +double Bilinear_interpolation(Landpoint p0, Landpoint p11, Landpoint p21, Landpoint p12, + Landpoint p22); + +bool onSegment(Point3 Pi, Point3 Pj, Point3 Q); + +Point3 invBilinear(Point3 p, Point3 a, Point3 b, Point3 c, Point3 d); + +// +// WGS84 到J2000 坐标系的变换 +// 参考网址：https://blog.csdn.net/hit5067/article/details/116894616 +// 资料网址：http://celestrak.org/spacedata/ +// 参数文件： +// a. Earth Orientation Parameter 文件： http://celestrak.org/spacedata/EOP-Last5Years.csv +// b. Space Weather Data 文件： http://celestrak.org/spacedata/SW-Last5Years.csv +// 备注：上述文件是自2017年-五年内 +/** +在wgs84 坐标系转到J2000 坐标系 主要 涉及到坐标的相互转换。一般给定的WGS坐标为 给定时刻的 t ,BLH +转换步骤： +step 1: WGS 84 转换到协议地球坐标系 +step 2: 协议地球坐标系 转换为瞬时地球坐标系 +step 3: 瞬时地球坐标系 转换为 瞬时真天球坐标系 +step 4: 瞬时真天球坐标系 转到瞬时平天球 坐标系 +step 5: 瞬时平天球坐标系转换为协议天球坐标系（J2000） +**/ + +double sind(double degree); + +double cosd(double d); + +// 插值 +ErrorCode polyfit(const double* x, const double* y, int xyLength, int poly_n, std::vector& out_factor, double& out_chisq); + +// 叉乘 +Point3 crossProduct(const Point3& a, const Point3& b); + +Eigen::Matrix3d rotationMatrix(const Eigen::Vector3d& axis, double angle); + +long double convertToMilliseconds(const std::string& dateTimeStr); + + +/// +/// list 应该是按照从小到大的顺序排好 +/// +/// +/// +/// +long FindValueInStdVector(std::vector& list,double& findv); + +long InsertValueInStdVector(std::vector& list, double insertValue, bool repeatValueInsert = false); + +long FindValueInStdVectorLast(std::vector& list, double& findv); + + +#endif \ No newline at end of file diff --git a/EchoDataFormat.cpp b/EchoDataFormat.cpp new file mode 100644 index 0000000..2e24a7b --- /dev/null +++ b/EchoDataFormat.cpp @@ -0,0 +1,613 @@ +#include "stdafx.h" +#include "EchoDataFormat.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include "ImageOperatorBase.h" + + +std::shared_ptr CreatePluseAntPosArr(long pluseCount) +{ + return std::shared_ptr(new PluseAntPos[pluseCount],delArrPtr); +} + +long getPluseDataSize(PluseData& pluseData) +{ + long datasize = sizeof(long) + sizeof(double) + sizeof(double) * 6 + sizeof(long);// PRFId,time,px-vz,pluseCount + datasize = datasize + pluseData.plusePoints * sizeof(double) * 2; // ݿ + return datasize; +} + +ErrorCode getPluseDataFromBuffer(char* buffer, PluseData & data) +{ + long seekid = 0; + std::memcpy(&data.id, buffer + seekid, sizeof(data.id)); seekid += sizeof(data.id); + std::memcpy(&data.time, buffer + seekid, sizeof(data.time)); seekid += sizeof(data.time); + std::memcpy(&data.Px, buffer + seekid, sizeof(data.Px)); seekid += sizeof(data.Px); + std::memcpy(&data.Py, buffer + seekid, sizeof(data.Py)); seekid += sizeof(data.Py); + std::memcpy(&data.Pz, buffer + seekid, sizeof(data.Pz)); seekid += sizeof(data.Pz); + std::memcpy(&data.Vx, buffer + seekid, sizeof(data.Vx)); seekid += sizeof(data.Vx); + std::memcpy(&data.Vy, buffer + seekid, sizeof(data.Vy)); seekid += sizeof(data.Vy); + std::memcpy(&data.Vz, buffer + seekid, sizeof(data.Vz)); seekid += sizeof(data.Vz); + std::memcpy(&data.plusePoints, buffer + seekid, sizeof(data.plusePoints)); seekid += sizeof(data.plusePoints); + + Eigen::MatrixXd echoData_real = Eigen::MatrixXd::Zero(1, data.plusePoints); + Eigen::MatrixXd echoData_imag = Eigen::MatrixXd::Zero(1, data.plusePoints); + std::memcpy(echoData_real.data(), buffer + seekid, sizeof(data.plusePoints)); seekid += data.plusePoints * sizeof(double); + std::memcpy(echoData_imag.data(), buffer + seekid, sizeof(data.plusePoints)); + + std::shared_ptr echoData = std::make_shared(1, data.plusePoints); + echoData->real() = echoData_real.array(); + echoData->imag() = echoData_imag.array(); + return ErrorCode::SUCCESS; +} + +std::shared_ptr CreatePluseDataArr(long pluseCount) +{ + return std::shared_ptr(new PluseData[pluseCount],delArrPtr); +} + +std::shared_ptr> CreateEchoData(long plusePoints) +{ + return std::shared_ptr>(new std::complex[plusePoints],delArrPtr); +} + +EchoL0Dataset::EchoL0Dataset() +{ + this->folder=""; + this->filename=""; + this->xmlname=""; + this->GPSPointFilename=""; + this->echoDataFilename=""; + this->xmlFilePath=""; + this->GPSPointFilePath=""; + this->echoDataFilePath=""; + this->simulationTaskName = ""; + + this->PluseCount = 0; + this->PlusePoints = 0; + this->NearRange = 0; + this->FarRange = 0; + this->centerFreq = 0; + this->Fs = 0; + + this->folder.clear(); + this->filename.clear(); + this->xmlname.clear(); + this->GPSPointFilename.clear(); + this->echoDataFilename.clear(); + this->xmlFilePath.clear(); + this->GPSPointFilePath.clear(); + this->echoDataFilePath.clear(); + this->simulationTaskName.clear(); + +} + +EchoL0Dataset::~EchoL0Dataset() +{ +} + +ErrorCode EchoL0Dataset::OpenOrNew(QString folder, QString filename, long PluseCount, long PlusePoints) +{ + qDebug() << "--------------Echo Data OpenOrNew ---------------------------------------"; + qDebug() << "Folder: " << folder; + qDebug() << "Filename: " << filename; + QDir dir(folder); + if (dir.exists() == false) + { + dir.mkpath(folder); + } + else {} + + if (dir.exists() == false) { + return ErrorCode::FOLDER_NOT_EXIST; + } + else {} + QString filePath = dir.filePath(filename); // ļ· + + this->folder = folder; + this->filename = filename; + this->simulationTaskName = filename; + + this->xmlname=filename+".xml"; + this->GPSPointFilename=filename+".gpspos.data"; + this->echoDataFilename = filename + ".L0echo.data"; + + this->xmlFilePath = dir.filePath(this->xmlname); + this->GPSPointFilePath = dir.filePath(this->GPSPointFilename); + this->echoDataFilePath = dir.filePath(this->echoDataFilename); + + this->PluseCount = PluseCount; + this->PlusePoints = PlusePoints; + + + // + if (QFile(this->xmlFilePath).exists()) + { + this->loadFromXml(); + } + else + { + this->saveToXml(); + } + + if (QFile(this->GPSPointFilePath).exists() == false) { + // ļ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + std::shared_ptr poDstDS(poDriver->Create(this->GPSPointFilePath.toUtf8().constData(), 16, PluseCount, 1, GDT_Float64, NULL)); + GDALFlushCache((GDALDatasetH)poDstDS.get()); + poDstDS.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + + } + + if (QFile(this->echoDataFilePath).exists() == false) { + + // ļ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + std::shared_ptr poDstDS(poDriver->Create(this->echoDataFilePath.toUtf8().constData(), PlusePoints, PluseCount, 1, GDT_CFloat64, NULL)); + GDALFlushCache((GDALDatasetH)poDstDS.get()); + poDstDS.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + + } + + + + return ErrorCode::SUCCESS; +} + +ErrorCode EchoL0Dataset::Open(QString xmlfilepath) +{ + QFileInfo fileInfo(xmlfilepath); + QString fileName = fileInfo.fileName(); // ȡļ + QString fileSuffix = fileInfo.suffix(); // ȡ׺ + QString fileNameWithoutSuffix = fileInfo.baseName(); // ȡ׺ļ + QString directoryPath = fileInfo.path(); // ȡļ· + if (fileSuffix.toLower() == "xml" || fileSuffix.toLower() == ".xml") { + return this->Open(directoryPath,fileNameWithoutSuffix); + } + else { + return ErrorCode::ECHO_L0DATA_XMLNAMEERROR; + } + + return ErrorCode::SUCCESS; +} + +ErrorCode EchoL0Dataset::Open(QString folder, QString filename) +{ + QDir dir(folder); + if (dir.exists() == false) + { + return ErrorCode::FOLDER_NOT_EXIST; + } + else {} + + if (dir.exists() == false) { + return ErrorCode::FOLDER_NOT_EXIST; + } + else {} + QString filePath = dir.filePath(filename); // ļ· + + this->folder = folder; + this->filename = filename; + this->simulationTaskName = filename; + + this->xmlname = filename + ".xml"; + this->GPSPointFilename = filename + ".gpspos.data"; + this->echoDataFilename = filename + ".L0echo.data"; + + this->xmlFilePath = dir.filePath(this->xmlname); + this->GPSPointFilePath = dir.filePath(this->GPSPointFilename); + this->echoDataFilePath = dir.filePath(this->echoDataFilename); + + this->PluseCount = PluseCount; + this->PlusePoints = PlusePoints; + + // + if (QFile(this->xmlFilePath).exists()) + { + this->loadFromXml(); + } + else + { + return ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR; + } +} + +QString EchoL0Dataset::getxmlName() +{ + return xmlname; +} + +QString EchoL0Dataset::getGPSPointFilename() +{ + return GPSPointFilename; +} + +QString EchoL0Dataset::getEchoDataFilename() +{ + return GPSPointFilePath; +} + +// Getter Setter ʵ +long EchoL0Dataset::getPluseCount() { return this->PluseCount; } +void EchoL0Dataset::setPluseCount(long pulseCount) { this->PluseCount = pulseCount; } + +long EchoL0Dataset::getPlusePoints() { return this->PlusePoints; } +void EchoL0Dataset::setPlusePoints(long pulsePoints) { this->PlusePoints = pulsePoints; } + + +double EchoL0Dataset::getNearRange() { return this->NearRange; } +void EchoL0Dataset::setNearRange(double nearRange) { this->NearRange = nearRange; } + +double EchoL0Dataset::getFarRange() { return this->FarRange; } +void EchoL0Dataset::setFarRange(double farRange) { this->FarRange = farRange; } + +double EchoL0Dataset::getCenterFreq() { return this->centerFreq; } +void EchoL0Dataset::setCenterFreq(double freq) { this->centerFreq = freq; } + +double EchoL0Dataset::getFs() { return this->Fs; } +void EchoL0Dataset::setFs(double samplingFreq) { this->Fs = samplingFreq; } + +QString EchoL0Dataset::getSimulationTaskName() { return this->simulationTaskName; } +void EchoL0Dataset::setSimulationTaskName(const QString& taskName) { this->simulationTaskName = taskName; } + +double EchoL0Dataset::getCenterAngle() +{ + return this->CenterAngle; +} + +void EchoL0Dataset::setCenterAngle(double angle) +{ + this->CenterAngle = angle; +} + +QString EchoL0Dataset::getLookSide() +{ + return this->LookSide; +} + +void EchoL0Dataset::setLookSide(QString lookside) +{ + this->LookSide = lookside; +} + +SatelliteAntPos EchoL0Dataset::getSatelliteAntPos(long prf_id) +{ + std::shared_ptr antpos = this->getAntPos(); + SatelliteAntPos prfpos{}; + prfpos.time = antpos.get()[prf_id * 16 + 0]; + prfpos.Px = antpos.get()[prf_id * 16 + 1]; + prfpos.Py = antpos.get()[prf_id * 16 + 2]; + prfpos.Pz = antpos.get()[prf_id * 16 + 3]; + prfpos.Vx = antpos.get()[prf_id * 16 + 4]; + prfpos.Vy = antpos.get()[prf_id * 16 + 5]; + prfpos.Vz = antpos.get()[prf_id * 16 + 6]; + prfpos.AntDirectX = antpos.get()[prf_id * 16 + 7]; + prfpos.AntDirectY = antpos.get()[prf_id * 16 + 8]; + prfpos.AntDirectZ = antpos.get()[prf_id * 16 + 9]; + prfpos.AVx = antpos.get()[prf_id * 16 + 10]; + prfpos.AVy = antpos.get()[prf_id * 16 + 11]; + prfpos.AVz =antpos.get()[prf_id * 16 + 12]; + prfpos.lon =antpos.get()[prf_id * 16 + 13]; + prfpos.lat =antpos.get()[prf_id * 16 + 14]; + prfpos.ati =antpos.get()[prf_id * 16 + 15]; + return prfpos; +} + +// ӡϢʵ +void EchoL0Dataset::printInfo() { + std::cout << "Simulation Task Name: " << this->simulationTaskName.toStdString() << std::endl; + std::cout << "Pulse Count: " << this->PluseCount << std::endl; + std::cout << "Pulse Points: " << this->PlusePoints << std::endl; + std::cout << "Near Range: " << this->NearRange << std::endl; + std::cout << "Far Range: " << this->FarRange << std::endl; + std::cout << "Center Frequency: " << this->centerFreq << std::endl; + std::cout << "Sampling Frequency: " << this->Fs << std::endl; +} + +// xmlļд +void EchoL0Dataset::saveToXml() { + + QString filePath = this->xmlFilePath; + + QFile file(filePath); + if (!file.open(QIODevice::WriteOnly)) { + qWarning() << "Cannot open file for writing:" << filePath; + return; + } + + QXmlStreamWriter xmlWriter(&file); + xmlWriter.setAutoFormatting(true); + xmlWriter.writeStartDocument(); + xmlWriter.writeStartElement("SimulationConfig"); + + xmlWriter.writeTextElement("PluseCount", QString::number(this->PluseCount)); + xmlWriter.writeTextElement("PlusePoints", QString::number(this->PlusePoints)); + xmlWriter.writeTextElement("NearRange", QString::number(this->NearRange)); + xmlWriter.writeTextElement("FarRange", QString::number(this->FarRange)); + xmlWriter.writeTextElement("CenterFreq", QString::number(this->centerFreq)); + xmlWriter.writeTextElement("Fs", QString::number(this->Fs)); + xmlWriter.writeTextElement("CenterAngle", QString::number(this->CenterAngle)); + xmlWriter.writeTextElement("LookSide", this->LookSide); + xmlWriter.writeTextElement("SimulationTaskName", this->simulationTaskName); + xmlWriter.writeTextElement("Filename", this->filename); + xmlWriter.writeTextElement("Xmlname", this->xmlname); + xmlWriter.writeTextElement("GPSPointFilename", this->GPSPointFilename); + xmlWriter.writeTextElement("EchoDataFilename", this->echoDataFilename); + + xmlWriter.writeEndElement(); // SimulationConfig + xmlWriter.writeEndDocument(); + + file.close(); +} + +ErrorCode EchoL0Dataset::loadFromXml() { + QString filePath = this->xmlFilePath; + QFile file(filePath); + if (!file.open(QIODevice::ReadOnly)) { + qWarning() << "Cannot open file for reading:" << filePath; + return ErrorCode::FILEOPENFAIL; + } + + + bool PluseCountflag = false; + bool PlusePointsflag = false; + bool NearRangeflag = false; + bool FarRangeflag = false; + bool CenterFreqflag = false; + bool Fsflag = false; + + + QXmlStreamReader xmlReader(&file); + while (!xmlReader.atEnd() && !xmlReader.hasError()) { + xmlReader.readNext(); + + if (xmlReader.isStartElement()) { + QString elementName = xmlReader.name().toString(); + if (elementName == "PluseCount") { + this->PluseCount = xmlReader.readElementText().toLong(); + PluseCountflag = true; + } + else if (elementName == "PlusePoints") { + this->PlusePoints = xmlReader.readElementText().toLong(); + PlusePointsflag = true; + } + else if (elementName == "NearRange") { + this->NearRange = xmlReader.readElementText().toDouble(); + NearRangeflag = true; + } + else if (elementName == "FarRange") { + this->FarRange = xmlReader.readElementText().toDouble(); + FarRangeflag = true; + } + else if (elementName == "CenterFreq") { + this->centerFreq = xmlReader.readElementText().toDouble(); + CenterFreqflag = true; + } + else if (elementName == "Fs") { + this->Fs = xmlReader.readElementText().toDouble(); + Fsflag = true; + } + else if (elementName == "SimulationTaskName") { + this->simulationTaskName = xmlReader.readElementText(); + } + else if (elementName == "Filename") { + this->filename = xmlReader.readElementText(); + } + else if (elementName == "Xmlname") { + this->xmlname = xmlReader.readElementText(); + } + else if (elementName == "GPSPointFilename") { + this->GPSPointFilename = xmlReader.readElementText(); + } + else if (elementName == "EchoDataFilename") { + this->echoDataFilename = xmlReader.readElementText(); + } + else if (elementName == "CenterAngle") { + this->CenterAngle = xmlReader.readElementText().toDouble(); + } + else if (elementName == "LookSide") { + this->LookSide = xmlReader.readElementText(); + } + } + } + + if (xmlReader.hasError()) { + qWarning() << "XML error:" << xmlReader.errorString(); + file.close(); + return ErrorCode::ECHO_L0DATA_XMLFILENOTOPEN; + } + + if (!(PlusePointsflag && PlusePointsflag && FarRangeflag && NearRangeflag && CenterFreqflag && Fsflag)) { + file.close(); + return ErrorCode::ECHO_L0DATA_XMLFILENOTOPEN; + } + + file.close(); + + return ErrorCode::SUCCESS; +} + +std::shared_ptr EchoL0Dataset::getAntPos() +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + // ȡļ + std::shared_ptr rasterDataset = OpenDataset(this->GPSPointFilePath, GDALAccess::GA_ReadOnly); + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* demBand = rasterDataset->GetRasterBand(1); + + long width = rasterDataset->GetRasterXSize(); + long height = rasterDataset->GetRasterYSize(); + long band_num = rasterDataset->GetRasterCount(); + + std::shared_ptr temp = nullptr; + + if (gdal_datatype == GDT_Float64) { + temp=std::shared_ptr(new double[this->PluseCount * 16],delArrPtr); + demBand->RasterIO(GF_Read, 0, 0, 16, this->PluseCount, temp.get(), 16, this->PluseCount, gdal_datatype, 0, 0); + } + else { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_GPSFILEFORMATERROR)) ; + } + rasterDataset.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return temp; +} + +std::shared_ptr> EchoL0Dataset::getEchoArr(long startPRF, long PRFLen) +{ + if (!(startPRF < this->PluseCount)) { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_PRFIDXOUTRANGE))<PluseCount; + return nullptr; + } + else {} + + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + std::shared_ptr rasterDataset = OpenDataset(this->echoDataFilePath, GDALAccess::GA_ReadOnly); + + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* poBand = rasterDataset->GetRasterBand(1); + if (NULL==poBand||nullptr == poBand) { + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + return nullptr; + } + else {} + long width = rasterDataset->GetRasterXSize(); + long height = rasterDataset->GetRasterYSize(); + long band_num = rasterDataset->GetRasterCount(); + std::shared_ptr> temp = nullptr; + if (height != this->PluseCount || width != this->PlusePoints) { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + } + else { + if (gdal_datatype == GDT_CFloat64) { + temp= std::shared_ptr>(new std::complex[PRFLen * width ], delArrPtr); + poBand->RasterIO(GF_Read, 0, startPRF, width, PRFLen, temp.get(), width, PRFLen, GDT_CFloat64, 0, 0); + GDALFlushCache((GDALDatasetH)rasterDataset.get()); + } + else { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + } + } + + rasterDataset.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return temp; +} + +std::shared_ptr> EchoL0Dataset::getEchoArr() +{ + return this->getEchoArr(0,this->PluseCount); +} + +ErrorCode EchoL0Dataset::saveAntPos(std::shared_ptr ptr) +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + // ȡļ + std::shared_ptr rasterDataset = OpenDataset(this->GPSPointFilePath, GDALAccess::GA_Update); + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* demBand = rasterDataset->GetRasterBand(1); + + long width = rasterDataset->GetRasterXSize(); + long height = rasterDataset->GetRasterYSize(); + long band_num = rasterDataset->GetRasterCount(); + + if (gdal_datatype == GDT_Float64) { + demBand->RasterIO(GF_Write, 0, 0, 16, this->PluseCount, ptr.get(), 16, this->PluseCount, gdal_datatype, 0, 0); + } + else { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_GPSFILEFORMATERROR)); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return ErrorCode::ECHO_L0DATA_GPSFILEFORMATERROR; + } + rasterDataset.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return ErrorCode::SUCCESS; +} + +ErrorCode EchoL0Dataset::saveEchoArr(std::shared_ptr> echoPtr, long startPRF, long PRFLen) +{ + if (!(startPRF < this->PluseCount)) { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_PRFIDXOUTRANGE)); + return ErrorCode::ECHO_L0DATA_PRFIDXOUTRANGE; + } + else {} + + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + std::shared_ptr rasterDataset = OpenDataset(this->echoDataFilePath, GDALAccess::GA_Update); + + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* poBand = rasterDataset->GetRasterBand(1); + if (NULL == poBand || nullptr == poBand) { + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + return ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR; + } + else {} + long width = rasterDataset->GetRasterXSize(); + long height = rasterDataset->GetRasterYSize(); + long band_num = rasterDataset->GetRasterCount(); + + if (height != this->PluseCount || width != this->PlusePoints) { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR; + } + else { + if (gdal_datatype == GDT_CFloat64) { + poBand->RasterIO(GF_Write, 0, startPRF, width, PRFLen, echoPtr.get(), width, PRFLen, GDT_CFloat64, 0, 0); + GDALFlushCache((GDALDatasetH)rasterDataset.get()); + } + else { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + } + } + rasterDataset.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return ErrorCode::SUCCESS; +} + \ No newline at end of file diff --git a/EchoDataFormat.h b/EchoDataFormat.h new file mode 100644 index 0000000..887baa4 --- /dev/null +++ b/EchoDataFormat.h @@ -0,0 +1,190 @@ +#pragma once +/*****************************************************************//** + * \file EchoDataFormat.h + * \brief 洢еĸƷݸʽҪزݸʽݿ + * + * \author 30453 + * \date October 2024 + *********************************************************************/ + + +#include "BaseConstVariable.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +#include "BaseTool.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include // for CPLMalloc() + +//======================================================================== +// زʽ +// file type: +// PRFCOUNT +// PRFPOINT +// nearRange +// farRange +// PRF1,time,Px,Py,Pz,Vx,Vy,Vz +// PRF2,time,Px,Py,Pz,Vx,Vy,Vz +// ------------------ ļ -------------------------------------- +// PRF1,time,Px,Py,Pz,Vx,Vy,Vz,PRFPOINT,RealData,imagData +// PRF2,time,Px,Py,Pz,Vx,Vy,Vz,PRFPOINT,RealData,imagData +// -------------------------------------------------------------- +// עBp˳ֻǹעλãĬϰĽ˳֯ +//======================================================================== + +// õĻزļʽ + +/// +/// ز-- +/// +struct PluseData { + long id; // PRF id + double time;// ʱ + double Px; // + double Py; + double Pz; + double Vx; // ٶ + double Vy; + double Vz; + long plusePoints; // + std::shared_ptr> echoData; // ز +}; + +long getPluseDataSize(PluseData& pluseData); +ErrorCode getPluseDataFromBuffer(char* buffer, PluseData& data); +std::shared_ptr CreatePluseDataArr(long pluseCount); +std::shared_ptr> CreateEchoData(long plusePoints); + + + +/// +/// ̬ +/// +struct PluseAntPos { + long id; // PRF id + double time;// ʱ + double Px; // + double Py; + double Pz; + double Vx; // ٶ + double Vy; + double Vz; +}; +std::shared_ptr CreatePluseAntPosArr(long pluseCount); + + +// L0 +class EchoL0Dataset { + +public: + EchoL0Dataset(); + ~EchoL0Dataset(); + + +public: + ErrorCode OpenOrNew(QString folder, QString filename,long PluseCount,long PlusePoints); + ErrorCode Open(QString xmlfilepath); + ErrorCode Open(QString folder, QString filename); + QString getxmlName(); + QString getGPSPointFilename(); + QString getEchoDataFilename(); + +private: // Ʒ + QString folder; + QString filename; + QString xmlname; + QString GPSPointFilename; + QString echoDataFilename; + // + QString xmlFilePath; + QString GPSPointFilePath; + QString echoDataFilePath; + + + +public: // ļ + // Getter Setter + long getPluseCount() ; + void setPluseCount(long pulseCount); + + long getPlusePoints() ; + void setPlusePoints(long pulsePoints); + + double getNearRange() ; + void setNearRange(double nearRange); + + double getFarRange() ; + void setFarRange(double farRange); + + double getCenterFreq() ; + void setCenterFreq(double freq); + + double getFs() ; + void setFs(double samplingFreq); + + QString getSimulationTaskName() ; + void setSimulationTaskName(const QString& taskName); + + double getCenterAngle(); + void setCenterAngle(double angle); + + QString getLookSide(); + void setLookSide(QString lookside); + + SatelliteAntPos getSatelliteAntPos(long plusePRFID); + // ӡϢijԱ + void printInfo() ; + +private: // ز + long PluseCount; + long PlusePoints; + double NearRange; + double FarRange; + double centerFreq; // Ƶ + double Fs; // ЧƵ + QString simulationTaskName; + + double CenterAngle; + QString LookSide; + +public: // д XML ĺ + void saveToXml(); + ErrorCode loadFromXml(); + +public: + // ȡļ + std::shared_ptr getAntPos(); + std::shared_ptr> getEchoArr(long startPRF, long PRFLen); + std::shared_ptr> getEchoArr(); + //ļ + ErrorCode saveAntPos(std::shared_ptr ptr); // עΣգдǰǷȷ + ErrorCode saveEchoArr(std::shared_ptr> echoPtr, long startPRF, long PRFLen); + + + + + +}; + + + + diff --git a/FileOperator.cpp b/FileOperator.cpp new file mode 100644 index 0000000..5ca06f9 --- /dev/null +++ b/FileOperator.cpp @@ -0,0 +1,195 @@ +#include "stdafx.h" +#include "FileOperator.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +std::vector getFilelist(const QString& folderpath, const QString& filenameExtension, int (*logfun)(QString logtext, int value)) +{ + QString filenameExtensionStr = filenameExtension; + filenameExtensionStr = filenameExtensionStr.remove(0, 1); + std::vector filenames(0); + if (isExists(folderpath) && isDirectory(folderpath)) { + QDir directory(folderpath); + if (directory.exists() && directory.isReadable()) { + QFileInfoList fileList = directory.entryInfoList(QDir::Files | QDir::NoDotAndDotDot); + for (const QFileInfo& fileInfo : fileList) { +// qDebug() << fileInfo.filePath() + "\nExtension: (" + filenameExtensionStr + ", " + fileInfo.suffix() + ")"; + if (filenameExtensionStr == u8"*" || filenameExtensionStr == fileInfo.suffix()) { + filenames.push_back(fileInfo.filePath()); + } + if (logfun) { + logfun(fileInfo.filePath() + "\nExtension: (" + filenameExtensionStr + ", " + fileInfo.suffix() + ")", 4); + } + } + } + else { + qWarning() << "Folder does not exist or is not readable: " << folderpath; + } + return filenames; + } + else { + return std::vector(0); + } + +} + +QString getParantFolderNameFromPath(const QString& path) +{ + QDir directory(path); + directory.cdUp(); + QString parentPath = directory.absolutePath(); + return directory.dirName(); +} + +QString getParantFromPath(const QString& path) +{ + //qDebug() << path; + QDir directory(path); + directory.cdUp(); + QString parentPath = directory.absolutePath(); + //qDebug() << parentPath; + return parentPath; +} + +QString getFileNameFromPath(const QString& path) +{ + QFileInfo fileInfo(path); + return fileInfo.fileName(); +} + +bool isDirectory(const QString& path) +{ + QFileInfo fileinfo(path); + return fileinfo.isDir(); +} + +bool isExists(const QString& path) +{ + QFileInfo fileinfo(path); + return fileinfo.exists(); +} + +bool isFile(const QString& path) +{ + QFileInfo fileinfo(path); + return fileinfo.isFile(); +} + +int write_binfile(char* filepath, char* data, size_t data_len) +{ + FILE* pd = fopen(filepath, "w"); + if (NULL == pd) { + return 2; + } + //数据块首地址: "&a"，元素大小: "sizeof(unsigned __int8)"， 元素个数: "10"， 文件指针："pd" + fwrite(data, sizeof(char), data_len, pd); + fclose(pd); + return -1; +} + + char* read_textfile(char* text_path, int* length) +{ + char* data = NULL; + FILE* fp1 = fopen(text_path, "r"); + if (fp1 == NULL) { + return NULL; + } + else {} + // 读取文件 + fseek(fp1, 0, SEEK_END); + int data_length = ftell(fp1); + data = (char*)malloc((data_length + 1) * sizeof(char)); + rewind(fp1); + if (data_length == fread(data, sizeof(char), data_length, fp1)) { + data[data_length] = '\0'; // 文件尾 + } + else { + free(data); + fclose(fp1); + return NULL; + } + fclose(fp1); + *length = data_length + 1; + return data; +} + +bool exists_test(const QString& name) +{ + return isExists(name); +} + +size_t fsize(FILE* fp) +{ + size_t n; + fpos_t fpos; // 当前位置 + fgetpos(fp, &fpos); // 获取当前位置 + fseek(fp, 0, SEEK_END); + n = ftell(fp); + fsetpos(fp, &fpos); // 恢复之前的位置 + return n; +} + + +void removeFile(const QString& filePath) +{ + QFile file(filePath); + + if (file.exists()) { + if (file.remove()) { + qDebug() << "File removed successfully: " << filePath; + } + else { + qWarning() << "Failed to remove file: " << filePath; + } + } + else { + qDebug() << "File does not exist: " << filePath; + } +} + +unsigned long convertToULong(const QString& input) { + bool ok; // Used to check if the conversion was successful + unsigned long result = input.toULong(&ok); + + if (!ok) { + qWarning() << "Conversion to unsigned long failed for input: " << input; + } + + return result; +} + +void copyFile(const QString& sourcePath, const QString& destinationPath) { + QFile sourceFile(sourcePath); + QFile destinationFile(destinationPath); + qDebug() << QString("copy file ready !! from ") + sourcePath+" to "+destinationPath ; + if (sourceFile.exists()) { + if (sourceFile.copy(destinationPath)) { + qDebug() << QString("copy file sucessfully !! from ") + sourcePath+" to "+destinationPath ; + // 复制成功 + //QMessageBox::information(nullptr, u8"成功", u8"文件复制成功"); + } + else { + // 复制失败 + if(sourceFile.exists()){ + QMessageBox::critical(nullptr, QObject::tr("error"), QObject::tr("The file already exists !!")); + } + else{ + QMessageBox::critical(nullptr, QObject::tr("error"), QObject::tr("file copy error")); + } + } + } + else { + // 源文件不存在 + QMessageBox::warning(nullptr, QObject::tr("warning"), QObject::tr("Source file not found")); + } +} \ No newline at end of file diff --git a/FileOperator.h b/FileOperator.h new file mode 100644 index 0000000..b6b936c --- /dev/null +++ b/FileOperator.h @@ -0,0 +1,49 @@ +#pragma once + +#ifndef FILEOPERATOR_H +#define FILEOPERATOR_H +#include "BaseConstVariable.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +bool isDirectory(const QString& path); +bool isExists(const QString& path); +bool isFile(const QString& path); +void removeFile(const QString& filePath); +unsigned long convertToULong(const QString& input); +/// +/// 获取文件(绝对路径） +/// +/// +/// +/// +std::vector getFilelist(const QString& folderpath, const QString& FilenameExtension = ".*",int (*logfun)(QString logtext,int value)=nullptr); + +QString getParantFolderNameFromPath(const QString& path); + +QString getFileNameFromPath(const QString& path); + +int write_binfile(char* filepath, char* data, size_t data_len); + + char* read_textfile(char* text_path, int* length); + +bool exists_test(const QString& name); + +size_t fsize(FILE* fp); + +QString getParantFromPath(const QString& path); +void copyFile(const QString& sourcePath, const QString& destinationPath); +// QT FileOperator +#endif \ No newline at end of file diff --git a/GeoOperator.cpp b/GeoOperator.cpp new file mode 100644 index 0000000..6ae1182 --- /dev/null +++ b/GeoOperator.cpp @@ -0,0 +1,410 @@ +#include "stdafx.h" +#include "GeoOperator.h" +#include +#include +#include +#include +////#include +#include +#include +#include +#include +#include +#include +#include +#include +//#include +#include //#include "ogrsf_frmts.h" +#include +#include +#include "GeoOperator.h" + + + +Landpoint operator +(const Landpoint& p1, const Landpoint& p2) +{ + return Landpoint{ p1.lon + p2.lon,p1.lat + p2.lat,p1.ati + p2.ati }; +} + +Landpoint operator -(const Landpoint& p1, const Landpoint& p2) +{ + return Landpoint{ p1.lon - p2.lon,p1.lat - p2.lat,p1.ati - p2.ati }; +} + +bool operator ==(const Landpoint& p1, const Landpoint& p2) +{ + return p1.lat == p2.lat && p1.lon == p2.lon && p1.ati == p2.ati; +} + + + +Landpoint operator *(const Landpoint& p, double scale) +{ + return Landpoint{ + p.lon * scale, + p.lat * scale, + p.ati * scale + }; +} + + +Landpoint LLA2XYZ(const Landpoint& LLA) { + double L = LLA.lon * d2r; + double B = LLA.lat * d2r; + double H = LLA.ati; + + double sinB = sin(B); + double cosB = cos(B); + + //double N = a / sqrt(1 - e * e * sin(B) * sin(B)); + double N = a / sqrt(1 - eSquare * sinB * sinB); + Landpoint result = { 0,0,0 }; + result.lon = (N + H) * cosB * cos(L); + result.lat = (N + H) * cosB * sin(L); + //result.z = (N * (1 - e * e) + H) * sin(B); + result.ati = (N * (1 - 1 / f_inverse) * (1 - 1 / f_inverse) + H) * sinB; + return result; +} + +void LLA2XYZ(const Landpoint& LLA, Point3& XYZ) +{ + double L = LLA.lon * d2r; + double B = LLA.lat * d2r; + double H = LLA.ati; + + double sinB = sin(B); + double cosB = cos(B); + + //double N = a / sqrt(1 - e * e * sin(B) * sin(B)); + double N = a / sqrt(1 - eSquare * sinB * sinB); + Landpoint result = { 0,0,0 }; + XYZ.x = (N + H) * cosB * cos(L); + XYZ.y = (N + H) * cosB * sin(L); + //result.z = (N * (1 - e * e) + H) * sin(B); + XYZ.z = (N * (1 - 1 / f_inverse) * (1 - 1 / f_inverse) + H) * sinB; + +} + + +Eigen::MatrixXd LLA2XYZ(Eigen::MatrixXd landpoint) +{ + landpoint.col(0) = landpoint.col(0).array() * d2r; // lon + landpoint.col(1) = landpoint.col(1).array() * d2r; // lat + + Eigen::MatrixXd sinB = (landpoint.col(1).array().sin());//lat + Eigen::MatrixXd cosB = (landpoint.col(1).array().cos());//lat + + Eigen::MatrixXd N = a / ((1 - sinB.array().pow(2) * eSquare).array().sqrt()); + Eigen::MatrixXd result(landpoint.rows(), 3); + + result.col(0) = (N.array() + landpoint.col(2).array()) * cosB.array() * Eigen::cos(landpoint.col(0).array()).array(); //x + result.col(1) = (N.array() + landpoint.col(2).array()) * cosB.array() * Eigen::sin(landpoint.col(0).array()).array(); //y + + result.col(2) = (N.array() * (1 - 1 / f_inverse) * (1 - 1 / f_inverse) + landpoint.col(2).array()) * sinB.array(); //z + + return result; +} + + +Landpoint XYZ2LLA(const Landpoint& XYZ) { + double tmpX = XYZ.lon;// + double temY = XYZ.lat;// + double temZ = XYZ.ati; + + double curB = 0; + double N = 0; + double sqrtTempXY = sqrt(tmpX * tmpX + temY * temY); + double calB = atan2(temZ, sqrtTempXY); + int counter = 0; + double sinCurB = 0; + while (abs(curB - calB) * r2d > epsilon && counter < 25) + { + curB = calB; + sinCurB = sin(curB); + N = a / sqrt(1 - eSquare * sinCurB * sinCurB); + calB = atan2(temZ + N * eSquare * sinCurB, sqrtTempXY); + counter++; + } + + Landpoint result = { 0,0,0 }; + result.lon = atan2(temY, tmpX) * r2d; + result.lat = curB * r2d; + result.ati = temZ / sinCurB - N * (1 - eSquare); + return result; +} + + + + +double getAngle(const Landpoint& a, const Landpoint& b) +{ + double c = dot(a, b) / (getlength(a) * getlength(b)); + if (a.lon * b.lat - a.lat * b.lon >= 0) { + return acos(c > 1 ? 1 : c < -1 ? -1 : c) * r2d; + } + else { + return 360 - acos(c > 1 ? 1 : c < -1 ? -1 : c) * r2d; + } +} + +double dot(const Landpoint& p1, const Landpoint& p2) +{ + return p1.lat * p2.lat + p1.lon * p2.lon + p1.ati * p2.ati; +} + +double getlength(const Landpoint& p1) { + + return sqrt(dot(p1, p1)); +} + +Landpoint crossProduct(const Landpoint& a, const Landpoint& b) { + return Landpoint{ + a.lat * b.ati - a.ati * b.lat,//x + a.ati * b.lon - a.lon * b.ati,//y + a.lon * b.lat - a.lat * b.lon//z + }; +} + +float cross2d(Point3 a, Point3 b) +{ + return a.x * b.y - a.y * b.x; +} + +Point3 operator -(Point3 a, Point3 b) +{ + return Point3{ a.x - b.x, a.y - b.y, a.z - b.z }; +} + +Point3 operator +(Point3 a, Point3 b) +{ + return Point3{ a.x + b.x, a.y + b.y, a.z + b.z }; +} + +double operator /(Point3 a, Point3 b) +{ + return sqrt(pow(a.x, 2) + pow(a.y, 2)) / sqrt(pow(b.x, 2) + pow(b.y, 2)); +} + + +Landpoint getSlopeVector(const Landpoint& p0, const Landpoint& p1, const Landpoint& p2, const Landpoint& p3, const Landpoint& p4,bool inLBH) { + Landpoint n0, n1, n2, n3, n4; + if (inLBH) { + n0 = LLA2XYZ(p0); + n1 = LLA2XYZ(p1); + n2 = LLA2XYZ(p2); + n3 = LLA2XYZ(p3); + n4 = LLA2XYZ(p4); + } + else { + n0 = p0; + n1 = p1; + n2 = p2; + n3 = p3; + n4 = p4; + } + + Landpoint n01 = n1 - n0, n02 = n2 - n0, n03 = n3 - n0, n04 = n4 - n0; + // ��n01Ϊ�������� + Landpoint np01 = p1 - p0, np02 = p2 - p0, np03 = p3 - p0, np04 = p4 - p0; + double a2 = getAngle(Landpoint{ np01.lon,np01.lat,0 }, Landpoint{ np02.lon,np02.lat,0 });// 01->02 + double a3 = getAngle(Landpoint{ np01.lon,np01.lat,0 }, Landpoint{ np03.lon,np03.lat,0 });// 01->03 + double a4 = getAngle(Landpoint{ np01.lon,np01.lat,0 }, Landpoint{ np04.lon,np04.lat,0 });// 01->04 + //qDebug() << a2 << "\t" << a3 << "\t" << a4 << endl; + a2 = 360 - a2; + a3 = 360 - a3; + a4 = 360 - a4; + Landpoint N, W, S, E; + N = n01; + if (a2 >= a3 && a2 >= a4) { + W = n02; + if (a3 >= a4) { + S = n03; + E = n04; + } + else { + S = n04; + E = n03; + } + } + else if (a3 >= a2 && a3 >= a4) { + W = n03; + if (a2 >= a4) { + S = n02; + E = n04; + } + else { + S = n04; + E = n02; + } + } + else if (a4 >= a2 && a4 >= a3) + { + W = n04; + if (a2 >= a3) { + S = n02; + E = n03; + } + else { + S = n03; + E = n02; + } + } + return (crossProduct(N, W) + crossProduct(W, S) + crossProduct(S, E) + crossProduct(E, N)) * 0.25; +} + + + +double distance(const Vector3D& p1, const Vector3D& p2) +{ + double dx = p1.x - p2.x; + double dy = p1.y - p2.y; + double dz = p1.z - p2.z; + return std::sqrt(dx * dx + dy * dy + dz * dz); +} + +double pointToLineDistance(const Vector3D& point, const Vector3D& linePoint, const Vector3D& lineDirection) +{ + Vector3D pointToLine = { point.x - linePoint.x, point.y - linePoint.y, point.z - linePoint.z }; + + // 计算点到直线的投影点的位? + double t = (pointToLine.x * lineDirection.x + pointToLine.y * lineDirection.y + pointToLine.z * lineDirection.z) / + (lineDirection.x * lineDirection.x + lineDirection.y * lineDirection.y + lineDirection.z * lineDirection.z); + + // 计算投影? + Vector3D projection = { linePoint.x + t * lineDirection.x, linePoint.y + t * lineDirection.y, linePoint.z + t * lineDirection.z }; + + // 计算点到直线的距? + return distance(point, projection); +} + +Vector3D operator+(const Vector3D& p1, const Vector3D& p2) +{ + return Vector3D{ p1.x + p2.x,p1.y + p2.y,p1.z + p2.z }; +} + +Vector3D operator-(const Vector3D& p1, const Vector3D& p2) +{ + return Vector3D{ p1.x - p2.x,p1.y - p2.y,p1.z - p2.z }; +} + +bool operator==(const Vector3D& p1, const Vector3D& p2) +{ + return p1.x == p2.x && p1.y == p2.y && p1.z == p2.z; +} + +Vector3D operator*(const Vector3D& p, double scale) +{ + return Vector3D{ + p.x * scale, + p.y * scale, + p.z * scale + }; +} + +Vector3D operator*(double scale, const Vector3D& p) +{ + return Vector3D{ + p.x * scale, + p.y * scale, + p.z * scale + }; +} + +double getAngle(const Vector3D& a, const Vector3D& b) +{ + double c = dot(a, b) / (getlength(a) * getlength(b)); + if (a.x * b.y - a.y * b.x >= 0) { + return acos(c > 1 ? 1 : c < -1 ? -1 : c) * r2d; + } + else { + return 360 - acos(c > 1 ? 1 : c < -1 ? -1 : c) * r2d; + } +} + +double getCosAngle(const Vector3D& a, const Vector3D& b) +{ + double c = dot(a, b) / (getlength(a) * getlength(b)); + return acos(c > 1 ? 1 : c < -1 ? -1 : c) * r2d; +} + +double dot(const Vector3D& p1, const Vector3D& p2) +{ + return p1.y * p2.y + p1.x * p2.x + p1.z * p2.z; +} + +double getlength(const Vector3D& p1) +{ + return std::sqrt(std::pow(p1.x, 2) + std::pow(p1.y, 2) + std::pow(p1.z, 2)); +} + +Vector3D crossProduct(const Vector3D& a, const Vector3D& b) +{ + return Vector3D{ + a.y * b.z - a.z * b.y,//x + a.z * b.x - a.x * b.z,//y + a.x * b.y - a.y * b.x//z + }; +} + + +/// +/// n1 +/// n4 n0 n2 +/// n3 +/// +Vector3D getSlopeVector(const Vector3D& n0, const Vector3D& n1, const Vector3D& n2, const Vector3D& n3, const Vector3D& n4) +{ + Vector3D n01 = n1 - n0, n02 = n2 - n0, n03 = n3 - n0, n04 = n4 - n0; + return (crossProduct(n01, n04) + crossProduct(n04, n03) + crossProduct(n03, n02) + crossProduct(n02, n01)) * 0.25; +} + + +SphericalCoordinates cartesianToSpherical(const CartesianCoordinates& cartesian) +{ + SphericalCoordinates spherical; + + spherical.r = std::sqrt(cartesian.x * cartesian.x + cartesian.y * cartesian.y + cartesian.z * cartesian.z); + spherical.theta = std::acos(cartesian.z / spherical.r); + spherical.phi = std::atan2(cartesian.y, cartesian.x); + + return spherical; +} + +CartesianCoordinates sphericalToCartesian(const SphericalCoordinates& spherical) +{ + CartesianCoordinates cartesian; + + cartesian.x = spherical.r * std::sin(spherical.theta) * std::cos(spherical.phi); + cartesian.y = spherical.r * std::sin(spherical.theta) * std::sin(spherical.phi); + cartesian.z = spherical.r * std::cos(spherical.theta); + + return cartesian; +} + + +double getlocalIncAngle(Landpoint& satepoint, Landpoint& landpoint, Landpoint& slopeVector) { + Landpoint Rsc = satepoint - landpoint; // AB=B-A + //double R = getlength(Rsc); + //std::cout << R << endl; + double angle = getAngle(Rsc, slopeVector); + if (angle >= 180) { + return 360 - angle; + } + else { + return angle; + } +} + + +double getlocalIncAngle(Vector3D& satepoint, Vector3D& landpoint, Vector3D& slopeVector){ + Vector3D Rsc = satepoint - landpoint; // AB=B-A + //double R = getlength(Rsc); + //std::cout << R << endl; + double angle = getAngle(Rsc, slopeVector); + if (angle >= 180) { + return 360 - angle; + } + else { + return angle; + } +} diff --git a/GeoOperator.h b/GeoOperator.h new file mode 100644 index 0000000..fd1a684 --- /dev/null +++ b/GeoOperator.h @@ -0,0 +1,126 @@ +#pragma once + + +#ifndef _GEOOPERATOR_H +#define _GEOOPERATOR_H + +#include "BaseConstVariable.h" +#include +#include +#include +#include +#include +#include +#include + + +/// +/// 将经纬度转换为地固参心坐标系 +/// +/// 经纬度点--degree +/// 投影坐标系点 +Landpoint LLA2XYZ(const Landpoint& LLA); +void LLA2XYZ(const Landpoint& LLA,Point3& XYZ); +Eigen::MatrixXd LLA2XYZ(Eigen::MatrixXd landpoint); + +/// +/// 将地固参心坐标系转换为经纬度 +/// +/// 固参心坐标系 +/// 经纬度--degree +Landpoint XYZ2LLA(const Landpoint& XYZ); + + +Landpoint operator +(const Landpoint& p1, const Landpoint& p2); + +Landpoint operator -(const Landpoint& p1, const Landpoint& p2); + +bool operator ==(const Landpoint& p1, const Landpoint& p2); + +Landpoint operator *(const Landpoint& p, double scale); + +double getAngle(const Landpoint& a, const Landpoint& b); + +double dot(const Landpoint& p1, const Landpoint& p2); + +double getlength(const Landpoint& p1); + +Landpoint crossProduct(const Landpoint& a, const Landpoint& b); + + +Landpoint getSlopeVector(const Landpoint& p0, const Landpoint& p1, const Landpoint& p2, const Landpoint& p3, const Landpoint& p4, bool inLBH=true); + +double getlocalIncAngle(Landpoint& satepoint, Landpoint& landpoint, Landpoint& slopeVector); + +float cross2d(Point3 a, Point3 b); + +Point3 operator -(Point3 a, Point3 b); + +Point3 operator +(Point3 a, Point3 b); + +double operator /(Point3 a, Point3 b); + + + +// 矢量计算 +struct Vector3D { + double x, y, z; +}; + +// 计算两点之间的距离 +double distance(const Vector3D& p1, const Vector3D& p2); +// 计算点到直线的最短距离 +double pointToLineDistance(const Vector3D& point, const Vector3D& linePoint, const Vector3D& lineDirection); + +Vector3D operator +(const Vector3D& p1, const Vector3D& p2); + +Vector3D operator -(const Vector3D& p1, const Vector3D& p2); + +bool operator ==(const Vector3D& p1, const Vector3D& p2); + +Vector3D operator *(const Vector3D& p, double scale); + +Vector3D operator *(double scale,const Vector3D& p ); + +double getAngle(const Vector3D& a, const Vector3D& b); + +double getCosAngle(const Vector3D& a, const Vector3D& b); + +double dot(const Vector3D& p1, const Vector3D& p2); + +double getlength(const Vector3D& p1); + +Vector3D crossProduct(const Vector3D& a, const Vector3D& b); + + + +/// +/// n1 +/// n4 n0 n2 +/// n3 +/// +/// +/// +/// +/// +/// +/// +Vector3D getSlopeVector(const Vector3D& n0, const Vector3D& n1, const Vector3D& n2, const Vector3D& n3, const Vector3D& n4); + + +struct CartesianCoordinates { + double x, y, z; +}; + +struct SphericalCoordinates { + double r, theta, phi; +}; + +SphericalCoordinates cartesianToSpherical(const CartesianCoordinates& cartesian); + +CartesianCoordinates sphericalToCartesian(const SphericalCoordinates& spherical); + +double getlocalIncAngle(Vector3D& satepoint, Vector3D& landpoint, Vector3D& slopeVector); + + +#endif \ No newline at end of file diff --git a/ImageOperatorBase.cpp b/ImageOperatorBase.cpp new file mode 100644 index 0000000..80abc78 --- /dev/null +++ b/ImageOperatorBase.cpp @@ -0,0 +1,1512 @@ +#include "stdafx.h" +#include "ImageOperatorBase.h" +#include "BaseTool.h" +#include "GeoOperator.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include "FileOperator.h" +#include +#include +#include + + + +/** + * 输入数据是ENVI格式数据 + */ + +std::shared_ptr OpenDataset(const QString& in_path, GDALAccess rwmode) +{ + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDataset* dataset_ptr = (GDALDataset*)(GDALOpen(in_path.toUtf8().constData(), rwmode)); + std::shared_ptr rasterDataset(dataset_ptr, CloseDataset); + return rasterDataset; +} + +void CloseDataset(GDALDataset* ptr) +{ + GDALClose(ptr); + ptr = NULL; +} + +int TIFF2ENVI(QString in_tiff_path, QString out_envi_path) +{ + std::shared_ptr ds = OpenDataset(in_tiff_path); + const char* args[] = { "-of", "ENVI", NULL }; + GDALTranslateOptions* psOptions = GDALTranslateOptionsNew((char**)args, NULL); + GDALClose(GDALTranslate(out_envi_path.toUtf8().constData(), ds.get(), psOptions, NULL)); + GDALTranslateOptionsFree(psOptions); + return 0; +} + +int ENVI2TIFF(QString in_envi_path, QString out_tiff_path) +{ + std::shared_ptr ds = OpenDataset(in_envi_path); + const char* args[] = { "-of", "Gtiff", NULL }; + GDALTranslateOptions* psOptions = GDALTranslateOptionsNew((char**)args, NULL); + GDALClose(GDALTranslate(out_tiff_path.toUtf8().constData(), ds.get(), psOptions, NULL)); + GDALTranslateOptionsFree(psOptions); + return 0; +} + +int CreateDataset(QString new_file_path, int height, int width, int band_num, double* gt, + QString projection, GDALDataType gdal_dtype, bool need_gt) +{ + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + std::shared_ptr poDstDS(poDriver->Create(new_file_path.toUtf8().constData(), width, + height, band_num, gdal_dtype, NULL)); + if(need_gt) { + poDstDS->SetProjection(projection.toUtf8().constData()); + poDstDS->SetGeoTransform(gt); + } else { + } + GDALFlushCache((GDALDatasetH)poDstDS.get()); + return 0; +} + +int saveDataset(QString new_file_path, int start_line, int start_cols, int band_ids, int datacols, + int datarows, void* databuffer) +{ + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + std::shared_ptr poDstDS = OpenDataset(new_file_path, GA_Update); + GDALDataType gdal_datatype = poDstDS->GetRasterBand(1)->GetRasterDataType(); + poDstDS->GetRasterBand(band_ids)->RasterIO(GF_Write, start_cols, start_line, datacols, datarows, + databuffer, datacols, datarows, gdal_datatype, 0, 0); + GDALFlushCache(poDstDS.get()); + return 0; +} + +int block_num_pre_memory(int block_width, int height, GDALDataType gdal_datatype, double memey_size) +{ + // 计算大小 + int size_meta = 0; + + if(gdal_datatype == GDT_Byte) { + size_meta = 1; + } else if(gdal_datatype == GDT_UInt16) { + size_meta = 2; // 只有双通道才能构建 复数矩阵 + } else if(gdal_datatype == GDT_UInt16) { + size_meta = 2; + } else if(gdal_datatype == GDT_Int16) { + size_meta = 2; + } else if(gdal_datatype == GDT_UInt32) { + size_meta = 4; + } else if(gdal_datatype == GDT_Int32) { + size_meta = 4; + } + // else if (gdal_datatype == GDT_UInt64) { + // size_meta = 8; + // } + // else if (gdal_datatype == GDT_Int64) { + // size_meta = 8; + // } + else if(gdal_datatype == GDT_Float32) { + size_meta = 4; + } else if(gdal_datatype == GDT_Float64) { + size_meta = 4; + } else if(gdal_datatype == GDT_CInt16) { + size_meta = 2; + } else if(gdal_datatype == GDT_CInt32) { + size_meta = 2; + } else if(gdal_datatype == GDT_CFloat32) { + size_meta = 4; + } else if(gdal_datatype == GDT_CFloat64) { + size_meta = 8; + } else { + } + int block_num = int(memey_size / (size_meta * block_width)); + block_num = block_num > height ? height : block_num; // 行数 + block_num = block_num < 1 ? 1 : block_num; + return block_num; +} + +Eigen::Matrix +ReadComplexMatrixData(int start_line, int width, int line_num, + std::shared_ptr rasterDataset, GDALDataType gdal_datatype) +{ + int band_num = rasterDataset->GetRasterCount(); + if(gdal_datatype == 0) { + return Eigen::Matrix(0, 0); + } else if(gdal_datatype < 8) { + if(band_num != 2) { + return Eigen::Matrix(0, 0); + } + } else if(gdal_datatype < 12) { + if(band_num != 1) { + return Eigen::Matrix(0, 0); + } + + } else { + } + bool _flag = false; + Eigen::Matrix data_mat( + line_num * width, 2); // 必须强制行优先 + if(gdal_datatype == GDT_Byte) { + Eigen::MatrixX real_mat(line_num * width, 1); + Eigen::MatrixX imag_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, gdal_datatype, + 0, 0); // real + rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + imag_mat.data(), width, line_num, gdal_datatype, + 0, 0); // imag + data_mat.col(0) = (real_mat.array().cast()).array(); + data_mat.col(1) = (imag_mat.array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_UInt16) { + Eigen::MatrixX real_mat(line_num * width, 1); + Eigen::MatrixX imag_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, gdal_datatype, + 0, 0); // real + rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + imag_mat.data(), width, line_num, gdal_datatype, + 0, 0); // imag + data_mat.col(0) = (real_mat.array().cast()).array(); + data_mat.col(1) = (imag_mat.array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_Int16) { + Eigen::MatrixX real_mat(line_num * width, 1); + Eigen::MatrixX imag_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, gdal_datatype, + 0, 0); // real + rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + imag_mat.data(), width, line_num, gdal_datatype, + 0, 0); // imag + data_mat.col(0) = (real_mat.array().cast()).array(); + data_mat.col(1) = (imag_mat.array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_UInt32) { + Eigen::MatrixX real_mat(line_num * width, 1); + Eigen::MatrixX imag_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, gdal_datatype, + 0, 0); // real + rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + imag_mat.data(), width, line_num, gdal_datatype, + 0, 0); // imag + data_mat.col(0) = (real_mat.array().cast()).array(); + data_mat.col(1) = (imag_mat.array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_Int32) { + Eigen::MatrixX real_mat(line_num * width, 1); + Eigen::MatrixX imag_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, gdal_datatype, + 0, 0); // real + rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + imag_mat.data(), width, line_num, gdal_datatype, + 0, 0); // imag + data_mat.col(0) = (real_mat.array().cast()).array(); + data_mat.col(1) = (imag_mat.array().cast()).array(); + _flag = true; + } + // else if (gdal_datatype == GDT_UInt64) { + // Eigen::MatrixX real_mat(line_num * width, 1); + // Eigen::MatrixX imag_mat(line_num * width, 1); + // rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + //real_mat.data(), width, line_num, gdal_datatype, 0, 0); // real + // rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + //imag_mat.data(), width, line_num, gdal_datatype, 0, 0); // imag data_mat.col(0) = + //(real_mat.array().cast()).array(); data_mat.col(1) = + //(imag_mat.array().cast()).array(); _flag = true; + // } + // else if (gdal_datatype == GDT_Int64) { + // Eigen::MatrixX real_mat(line_num * width, 1); + // Eigen::MatrixX imag_mat(line_num * width, 1); + // rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + //real_mat.data(), width, line_num, gdal_datatype, 0, 0); // real + // rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + //imag_mat.data(), width, line_num, gdal_datatype, 0, 0); // imag data_mat.col(0) = + //(real_mat.array().cast()).array(); data_mat.col(1) = + //(imag_mat.array().cast()).array(); _flag = true; + // } + else if(gdal_datatype == GDT_Float32) { + Eigen::MatrixX real_mat(line_num * width, 1); + Eigen::MatrixX imag_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, gdal_datatype, + 0, 0); // real + rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + imag_mat.data(), width, line_num, gdal_datatype, + 0, 0); // imag + data_mat.col(0) = (real_mat.array().cast()).array(); + data_mat.col(1) = (imag_mat.array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_Float64) { + Eigen::MatrixX real_mat(line_num * width, 1); + Eigen::MatrixX imag_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, gdal_datatype, + 0, 0); // real + rasterDataset->GetRasterBand(2)->RasterIO(GF_Read, 0, start_line, width, line_num, + imag_mat.data(), width, line_num, gdal_datatype, + 0, 0); // imag + data_mat.col(0) = (real_mat.array().cast()).array(); + data_mat.col(1) = (imag_mat.array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_CInt16) { + Eigen::MatrixX> complex_short_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + complex_short_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = (complex_short_mat.real().array().cast()).array(); + data_mat.col(1) = (complex_short_mat.imag().array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_CInt32) { + Eigen::MatrixX> complex_short_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + complex_short_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = (complex_short_mat.real().array().cast()).array(); + data_mat.col(1) = (complex_short_mat.imag().array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_CFloat32) { + Eigen::MatrixX> complex_short_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + complex_short_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = (complex_short_mat.real().array().cast()).array(); + data_mat.col(1) = (complex_short_mat.imag().array().cast()).array(); + _flag = true; + } else if(gdal_datatype == GDT_CFloat64) { + Eigen::MatrixX> complex_short_mat(line_num * width, 1); + rasterDataset->GetRasterBand(1)->RasterIO(GF_Read, 0, start_line, width, line_num, + complex_short_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = (complex_short_mat.real().array().cast()).array(); + data_mat.col(1) = (complex_short_mat.imag().array().cast()).array(); + _flag = true; + } else { + } + // 保存数据 + + if(_flag) { + return data_mat; + } else { + return Eigen::Matrix( + 0, 0); // 必须强制行优先; + } +} + +Eigen::Matrix +ReadMatrixDoubleData(int start_line, int width, int line_num, + std::shared_ptr rasterDataset, GDALDataType gdal_datatype, + int band_idx) +{ + // 构建矩阵块，使用eigen 进行矩阵计算，加速计算 + bool _flag = false; + Eigen::Matrix data_mat( + line_num * width, 1); // 必须强制行优先 + if(gdal_datatype == GDT_Byte) { + Eigen::MatrixX real_mat(line_num * width, 1); + rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = ((real_mat.array().cast()).array().pow(2)).log10() * 10.0; + _flag = true; + } else if(gdal_datatype == GDT_UInt16) { + Eigen::MatrixX real_mat(line_num * width, 1); + rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = ((real_mat.array().cast()).array().pow(2)).log10() * 10.0; + _flag = true; + } else if(gdal_datatype == GDT_Int16) { + Eigen::MatrixX real_mat(line_num * width, 1); + rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = ((real_mat.array().cast()).array().pow(2)).log10() * 10.0; + _flag = true; + } else if(gdal_datatype == GDT_UInt32) { + Eigen::MatrixX real_mat(line_num * width, 1); + rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = ((real_mat.array().cast()).array().pow(2)).log10() * 10.0; + _flag = true; + } else if(gdal_datatype == GDT_Int32) { + Eigen::MatrixX real_mat(line_num * width, 1); + rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = ((real_mat.array().cast()).array().pow(2)).log10() * 10.0; + _flag = true; + } + // else if (gdal_datatype == GDT_UInt64) { + // Eigen::MatrixX real_mat(line_num * width, 1); + // rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + //real_mat.data(), width, line_num, gdal_datatype, 0, 0); // real data_mat.col(0) = + //((real_mat.array().cast()).array().pow(2)).log10() * 10.0; _flag = true; + // } + // else if (gdal_datatype == GDT_Int64) { + // Eigen::MatrixX real_mat(line_num * width, 1); + // rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + //real_mat.data(), width, line_num, gdal_datatype, 0, 0); // real data_mat.col(0) = + //((real_mat.array().cast()).array().pow(2)).log10() * 10.0; _flag = true; + // } + else if(gdal_datatype == GDT_Float32) { + Eigen::MatrixX real_mat(line_num * width, 1); + rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = ((real_mat.array().cast()).array().pow(2)).log10() * 10.0; + _flag = true; + } else if(gdal_datatype == GDT_Float64) { + Eigen::MatrixX real_mat(line_num * width, 1); + rasterDataset->GetRasterBand(band_idx)->RasterIO(GF_Read, 0, start_line, width, line_num, + real_mat.data(), width, line_num, + gdal_datatype, 0, 0); // real + data_mat.col(0) = ((real_mat.array().cast()).array().pow(2)).log10() * 10.0; + _flag = true; + } else { + } + + return data_mat; +} + +Eigen::MatrixXd getGeoTranslationArray(QString in_path) +{ + return Eigen::MatrixXd(); +} + +ImageGEOINFO getImageINFO(QString in_path) +{ + std::shared_ptr df = OpenDataset(in_path); + int width = df->GetRasterXSize(); + int heigh = df->GetRasterYSize(); + int band_num = df->GetRasterCount(); + ImageGEOINFO result; + result.width = width; + result.height = heigh; + result.bandnum = band_num; + + return result; +} + +GDALDataType getGDALDataType(QString fileptah) +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDataset* rasterDataset = (GDALDataset*)(GDALOpen( + fileptah.toUtf8().constData(), GA_ReadOnly)); // 锟斤拷只斤拷式锟斤拷取斤拷影锟斤拷 + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + + GDALClose((GDALDatasetH)rasterDataset); + omp_unset_lock(&lock); // 锟酵放伙拷斤拷 + omp_destroy_lock(&lock); // 劫伙拷斤拷 + + return gdal_datatype; +} + +gdalImage::gdalImage() +{ + this->height = 0; + this->width = 0; + this->data_band_ids = 1; + this->start_row = 0; + this->start_col = 0; +} + +/// +/// 斤拷图饺∮帮拷锟?1锟?7 +/// +/// +gdalImage::gdalImage(const QString& raster_path) +{ + omp_lock_t lock; + omp_init_lock(&lock); // 锟斤拷始斤拷斤拷 + omp_set_lock(&lock); // 锟斤拷没斤拷锟?1锟?7 + this->img_path = raster_path; + + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); // 注绞斤拷斤拷锟?1锟?7 + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + // 锟斤拷DEM影锟斤拷 + GDALDataset* rasterDataset = (GDALDataset*)(GDALOpen( + raster_path.toUtf8().constData(), GA_ReadOnly)); // 锟斤拷只斤拷式锟斤拷取斤拷影锟斤拷 + this->width = rasterDataset->GetRasterXSize(); + this->height = rasterDataset->GetRasterYSize(); + this->band_num = rasterDataset->GetRasterCount(); + + double* gt = new double[6]; + // 锟斤拷梅斤拷锟斤拷 + rasterDataset->GetGeoTransform(gt); + this->gt = Eigen::MatrixXd(2, 3); + this->gt << gt[0], gt[1], gt[2], gt[3], gt[4], gt[5]; + + this->projection = rasterDataset->GetProjectionRef(); + // 斤拷斤拷 + // double* inv_gt = new double[6];; + // GDALInvGeoTransform(gt, inv_gt); // 斤拷斤拷 + // 斤拷投影 + GDALFlushCache((GDALDatasetH)rasterDataset); + GDALClose((GDALDatasetH)rasterDataset); + rasterDataset = NULL; // 指矫匡拷 + this->InitInv_gt(); + delete[] gt; + ////GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + omp_unset_lock(&lock); // 锟酵放伙拷斤拷 + omp_destroy_lock(&lock); // 劫伙拷斤拷 +} + +gdalImage::~gdalImage() {} + +void gdalImage::setHeight(int height) +{ + this->height = height; +} + +void gdalImage::setWidth(int width) +{ + this->width = width; +} + +void gdalImage::setTranslationMatrix(Eigen::MatrixXd gt) +{ + this->gt = gt; +} + +void gdalImage::setData(Eigen::MatrixXd, int data_band_ids) +{ + this->data = data; + this->data_band_ids = data_band_ids; +} + +Eigen::MatrixXd gdalImage::getData(int start_row, int start_col, int rows_count, int cols_count, + int band_ids = 1) +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDataset* rasterDataset = (GDALDataset*)(GDALOpen( + this->img_path.toUtf8().constData(), GA_ReadOnly)); // 锟斤拷只斤拷式锟斤拷取斤拷影锟斤拷 + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* demBand = rasterDataset->GetRasterBand(band_ids); + + rows_count = start_row + rows_count <= this->height ? rows_count : this->height - start_row; + cols_count = start_col + cols_count <= this->width ? cols_count : this->width - start_col; + + Eigen::MatrixXd datamatrix(rows_count, cols_count); + + if(gdal_datatype == GDT_Byte) { + char* temp = new char[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + rows_count, gdal_datatype, 0, 0); + for(int i = 0; i < rows_count; i++) { + for(int j = 0; j < cols_count; j++) { + datamatrix(i, j) = temp[i * cols_count + j]; + } + } + delete[] temp; + } else if(gdal_datatype == GDT_UInt16) { + unsigned short* temp = new unsigned short[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + rows_count, gdal_datatype, 0, 0); + for(int i = 0; i < rows_count; i++) { + for(int j = 0; j < cols_count; j++) { + datamatrix(i, j) = temp[i * cols_count + j]; + } + } + delete[] temp; + } else if(gdal_datatype == GDT_Int16) { + short* temp = new short[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + rows_count, gdal_datatype, 0, 0); + for(int i = 0; i < rows_count; i++) { + for(int j = 0; j < cols_count; j++) { + datamatrix(i, j) = temp[i * cols_count + j]; + } + } + delete[] temp; + } else if(gdal_datatype == GDT_UInt32) { + unsigned int* temp = new unsigned int[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + rows_count, gdal_datatype, 0, 0); + for(int i = 0; i < rows_count; i++) { + for(int j = 0; j < cols_count; j++) { + datamatrix(i, j) = temp[i * cols_count + j]; + } + } + delete[] temp; + } else if(gdal_datatype == GDT_Int32) { + int* temp = new int[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + rows_count, gdal_datatype, 0, 0); + for(int i = 0; i < rows_count; i++) { + for(int j = 0; j < cols_count; j++) { + datamatrix(i, j) = temp[i * cols_count + j]; + } + } + delete[] temp; + } + // else if (gdal_datatype == GDT_UInt64) { + // unsigned long* temp = new unsigned long[rows_count * cols_count]; + // demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + //rows_count, gdal_datatype, 0, 0); for (int i = 0; i < rows_count; i++) { for (int j = 0; j < + //cols_count; j++) { datamatrix(i, j) = temp[i * cols_count + j]; + // } + // } + // delete[] temp; + // } + // else if (gdal_datatype == GDT_Int64) { + // long* temp = new long[rows_count * cols_count]; + // demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + //rows_count, gdal_datatype, 0, 0); for (int i = 0; i < rows_count; i++) { for (int j = 0; j < + //cols_count; j++) { datamatrix(i, j) = temp[i * cols_count + j]; + // } + // } + // delete[] temp; + // } + else if(gdal_datatype == GDT_Float32) { + float* temp = new float[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + rows_count, gdal_datatype, 0, 0); + + for(int i = 0; i < rows_count; i++) { + for(int j = 0; j < cols_count; j++) { + datamatrix(i, j) = temp[i * cols_count + j]; + } + } + delete[] temp; + } else if(gdal_datatype == GDT_Float64) { + double* temp = new double[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + rows_count, gdal_datatype, 0, 0); + + for(int i = 0; i < rows_count; i++) { + for(int j = 0; j < cols_count; j++) { + datamatrix(i, j) = temp[i * cols_count + j]; + } + } + delete[] temp; + } else { + } + GDALClose((GDALDatasetH)rasterDataset); + omp_unset_lock(&lock); // 锟酵放伙拷斤拷 + omp_destroy_lock(&lock); // 劫伙拷斤拷 + // GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + return datamatrix; +} + +Eigen::MatrixXd gdalImage::getGeoTranslation() +{ + return this->gt; +} + +GDALDataType gdalImage::getDataType() +{ + GDALDataset* rasterDataset = + (GDALDataset*)(GDALOpen(this->img_path.toUtf8().constData(), GA_ReadOnly)); + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + return gdal_datatype; +} + +/// +/// +/// +/// +/// +/// +/// +void gdalImage::saveImage(Eigen::MatrixXd data, int start_row = 0, int start_col = 0, + int band_ids = 1) +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + if(start_row + data.rows() > this->height || start_col + data.cols() > this->width) { + QString tip = u8"file path: " + this->img_path+" image size :( "+QString::number(this->height)+" , "+ QString::number(this->width)+" ) "+" input size ("+ QString::number(start_row + data.rows())+", "+ QString::number(start_col + data.cols())+") "; + qDebug() << tip; + throw std::exception(tip.toUtf8().constData()); + } + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("GTiff"); + GDALDataset* poDstDS = nullptr; + if(exists_test(this->img_path)) { + poDstDS = (GDALDataset*)(GDALOpen(this->img_path.toUtf8().constData(), GA_Update)); + } else { + poDstDS = poDriver->Create(this->img_path.toUtf8().constData(), this->width, this->height, + this->band_num, GDT_Float32, NULL); // 斤拷锟斤拷 + poDstDS->SetProjection(this->projection.toUtf8().constData()); + + double gt_ptr[6]; + for(int i = 0; i < this->gt.rows(); i++) { + for(int j = 0; j < this->gt.cols(); j++) { + gt_ptr[i * 3 + j] = this->gt(i, j); + } + } + poDstDS->SetGeoTransform(gt_ptr); + //delete gt_ptr; + } + + int datarows = data.rows(); + int datacols = data.cols(); + + float* databuffer = + new float[datarows * datacols]; // (float*)malloc(datarows * datacols * sizeof(float)); + + for(int i = 0; i < datarows; i++) { + for(int j = 0; j < datacols; j++) { + float temp = float(data(i, j)); + databuffer[i * datacols + j] = temp; + } + } + // poDstDS->RasterIO(GF_Write,start_col, start_row, datacols, datarows, databuffer, datacols, + // datarows, GDT_Float32,band_ids, num,0,0,0); + poDstDS->GetRasterBand(band_ids)->RasterIO(GF_Write, start_col, start_row, datacols, datarows, + databuffer, datacols, datarows, GDT_Float32, 0, 0); + + GDALFlushCache(poDstDS); + GDALClose((GDALDatasetH)poDstDS); + // GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + delete[] databuffer; + omp_unset_lock(&lock); // 锟酵放伙拷斤拷 + omp_destroy_lock(&lock); // 劫伙拷斤拷 +} + +void gdalImage::saveImage() +{ + this->saveImage(this->data, this->start_row, this->start_col, this->data_band_ids); +} + +void gdalImage::setNoDataValue(double nodatavalue = -9999, int band_ids = 1) +{ + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); // 注绞斤拷斤拷锟?1锟?7 + // GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("GTiff"); + GDALDataset* poDstDS = (GDALDataset*)(GDALOpen(img_path.toUtf8().constData(), GA_Update)); + poDstDS->GetRasterBand(band_ids)->SetNoDataValue(nodatavalue); + GDALFlushCache((GDALDatasetH)poDstDS); + GDALClose((GDALDatasetH)poDstDS); +} + +int gdalImage::InitInv_gt() +{ + // 1 lon lat = x + // 1 lon lat = y + Eigen::MatrixXd temp = Eigen::MatrixXd::Zero(2, 3); + this->inv_gt = temp; + double a = this->gt(0, 0); + double b = this->gt(0, 1); + double c = this->gt(0, 2); + double d = this->gt(1, 0); + double e = this->gt(1, 1); + double f = this->gt(1, 2); + double g = 1; + double det_gt = b * f - c * e; + if(det_gt == 0) { + return 0; + } + this->inv_gt(0, 0) = (c * d - a * f) / det_gt; // 2 + this->inv_gt(0, 1) = f / det_gt; // lon + this->inv_gt(0, 2) = -c / det_gt; // lat + this->inv_gt(1, 0) = (a * e - b * d) / det_gt; // 1 + this->inv_gt(1, 1) = -e / det_gt; // lon + this->inv_gt(1, 2) = b / det_gt; // lat + return 1; +} + +Landpoint gdalImage::getRow_Col(double lon, double lat) +{ + Landpoint p{ 0, 0, 0 }; + p.lon = this->inv_gt(0, 0) + lon * this->inv_gt(0, 1) + lat * this->inv_gt(0, 2); // x + p.lat = this->inv_gt(1, 0) + lon * this->inv_gt(1, 1) + lat * this->inv_gt(1, 2); // y + return p; +} + +Landpoint gdalImage::getLandPoint(double row, double col, double ati = 0) +{ + Landpoint p{ 0, 0, 0 }; + p.lon = this->gt(0, 0) + col * this->gt(0, 1) + row * this->gt(0, 2); // x + p.lat = this->gt(1, 0) + col * this->gt(1, 1) + row * this->gt(1, 2); // y + p.ati = ati; + return p; +} + +void gdalImage::getLandPoint(double row, double col, double ati, Landpoint& Lp) +{ + Lp.lon = this->gt(0, 0) + col * this->gt(0, 1) + row * this->gt(0, 2); // x + Lp.lat = this->gt(1, 0) + col * this->gt(1, 1) + row * this->gt(1, 2); // y + Lp.ati = ati; + +} + + + + + + +double gdalImage::mean(int bandids) +{ + double mean_value = 0; + double count = this->height * this->width; + int line_invert = 100; + int start_ids = 0; + do { + Eigen::MatrixXd sar_a = this->getData(start_ids, 0, line_invert, this->width, bandids); + mean_value = mean_value + sar_a.sum() / count; + start_ids = start_ids + line_invert; + } while(start_ids < this->height); + return mean_value; +} + +double gdalImage::BandmaxValue(int bandids) +{ + double max_value = 0; + bool state_max = true; + int line_invert = 100; + int start_ids = 0; + double temp_max = 0; + do { + Eigen::MatrixXd sar_a = this->getData(start_ids, 0, line_invert, this->width, bandids); + if(state_max) { + state_max = false; + max_value = sar_a.maxCoeff(); + } else { + temp_max = sar_a.maxCoeff(); + if(max_value < temp_max) { + max_value = temp_max; + } + } + start_ids = start_ids + line_invert; + } while(start_ids < this->height); + return max_value; +} + +double gdalImage::BandminValue(int bandids) +{ + double min_value = 0; + bool state_min = true; + int line_invert = 100; + int start_ids = 0; + double temp_min = 0; + do { + Eigen::MatrixXd sar_a = this->getData(start_ids, 0, line_invert, this->width, bandids); + if(state_min) { + state_min = false; + min_value = sar_a.minCoeff(); + } else { + temp_min = sar_a.minCoeff(); + if(min_value < temp_min) { + min_value = temp_min; + } + } + start_ids = start_ids + line_invert; + } while(start_ids < this->height); + return min_value; +} + +GDALRPCInfo gdalImage::getRPC() +{ + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO"); + CPLSetConfigOption("GDAL_DATA", "./data"); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); // 注斤拷锟斤拷 + // 斤拷锟斤拷 + GDALDataset* pDS = (GDALDataset*)GDALOpen(this->img_path.toUtf8().constData(), GA_ReadOnly); + // 锟斤拷元斤拷锟叫伙拷取RPC锟斤拷息 + char** papszRPC = pDS->GetMetadata("RPC"); + + // 斤拷取锟斤拷RPC锟斤拷息斤拷山峁癸拷锟?1锟?7 + GDALRPCInfo oInfo; + GDALExtractRPCInfo(papszRPC, &oInfo); + + GDALClose((GDALDatasetH)pDS); + + return oInfo; +} + +Eigen::MatrixXd gdalImage::getLandPoint(Eigen::MatrixXd points) +{ + if(points.cols() != 3) { + throw new std::exception("the size of points is equit 3!!!"); + } + + Eigen::MatrixXd result(points.rows(), 3); + result.col(2) = points.col(2); // 锟竭筹拷 + points.col(2) = points.col(2).array() * 0 + 1; + points.col(0).swap(points.col(2)); // 斤拷 + Eigen::MatrixXd gts(3, 2); + gts.col(0) = this->gt.row(0); + gts.col(1) = this->gt.row(1); + + result.block(0, 0, points.rows(), 2) = points * gts; + return result; +} + +Eigen::MatrixXd gdalImage::getHist(int bandids) +{ + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); // 注绞斤拷斤拷锟?1锟?7 + // 锟斤拷DEM影锟斤拷 + GDALDataset* rasterDataset = (GDALDataset*)(GDALOpen( + this->img_path.toUtf8().constData(), GA_ReadOnly)); // 锟斤拷只斤拷式锟斤拷取斤拷影锟斤拷 + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* xBand = rasterDataset->GetRasterBand(bandids); + + double dfMin = this->BandminValue(bandids); + double dfMax = this->BandmaxValue(bandids); + int count = int((dfMax - dfMin) / 0.01); + count = count > 255 ? count : 255; + GUIntBig* panHistogram = new GUIntBig[count]; + xBand->GetHistogram(dfMin, dfMax, count, panHistogram, TRUE, FALSE, NULL, NULL); + Eigen::MatrixXd result(count, 2); + double delta = (dfMax - dfMin) / count; + for(int i = 0; i < count; i++) { + result(i, 0) = dfMin + i * delta; + result(i, 1) = double(panHistogram[i]); + } + delete[] panHistogram; + GDALClose((GDALDatasetH)rasterDataset); + return result; +} + +gdalImage CreategdalImage(const QString& img_path, int height, int width, int band_num, + Eigen::MatrixXd gt, QString projection, bool need_gt, bool overwrite) +{ + if(exists_test(img_path.toUtf8().constData())) { + if(overwrite) { + gdalImage result_img(img_path); + return result_img; + } else { + throw "file has exist!!!"; + exit(1); + } + } + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); // 注锟斤拷锟绞斤拷锟斤拷锟斤拷锟?1锟?7 + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("GTiff"); + GDALDataset* poDstDS = poDriver->Create(img_path.toUtf8().constData(), width, height, band_num, + GDT_Float32, NULL); // 锟斤拷锟斤拷锟斤拷 + if(need_gt) { + poDstDS->SetProjection(projection.toUtf8().constData()); + double gt_ptr[6] = { 0 }; + for(int i = 0; i < gt.rows(); i++) { + for(int j = 0; j < gt.cols(); j++) { + gt_ptr[i * 3 + j] = gt(i, j); + } + } + poDstDS->SetGeoTransform(gt_ptr); + } + for(int i = 1; i <= band_num; i++) { + poDstDS->GetRasterBand(i)->SetNoDataValue(-9999); + } + GDALFlushCache((GDALDatasetH)poDstDS); + GDALClose((GDALDatasetH)poDstDS); + ////GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + gdalImage result_img(img_path); + return result_img; +} + +gdalImageComplex CreategdalImageComplex(const QString& img_path, int height, int width, + int band_num, Eigen::MatrixXd gt, QString projection, + bool need_gt, bool overwrite) +{ + if(exists_test(img_path.toUtf8().constData())) { + if(overwrite) { + gdalImageComplex result_img(img_path); + return result_img; + } else { + throw "file has exist!!!"; + exit(1); + } + } + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + GDALDataset* poDstDS = poDriver->Create(img_path.toUtf8().constData(), width, height, band_num, + GDT_CFloat64, NULL); + if(need_gt) { + poDstDS->SetProjection(projection.toUtf8().constData()); + + // 锟斤拷锟斤拷转锟斤拷锟斤拷锟斤拷 + double gt_ptr[6] = { 0 }; + for(int i = 0; i < gt.rows(); i++) { + for(int j = 0; j < gt.cols(); j++) { + gt_ptr[i * 3 + j] = gt(i, j); + } + } + poDstDS->SetGeoTransform(gt_ptr); + } + //for(int i = 1; i <= band_num; i++) { + // poDstDS->GetRasterBand(i)->SetNoDataValue(0); // 回波部分 + //} + GDALFlushCache((GDALDatasetH)poDstDS); + GDALClose((GDALDatasetH)poDstDS); + ////GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + gdalImageComplex result_img(img_path); + return result_img; +} + +gdalImageComplex CreateEchoComplex(const QString& img_path, int height, int width, int band_num) +{ + // 创建图像 + Eigen::MatrixXd gt = Eigen::MatrixXd::Zero(2, 3); + //Xgeo = GeoTransform[0] + Xpixel * GeoTransform[1] + Ypixel * GeoTransform[2] + //Ygeo = GeoTransform[3] + Xpixel * GeoTransform[4] + Ypixel * GeoTransform[5] + // X + gt(0, 0) = 0; gt(0, 2) = 1; gt(0, 2) = 0; + gt(1, 0) = 0; gt(1, 1) = 0; gt(1, 2) = 1; + // Y + QString projection = ""; + gdalImageComplex echodata = CreategdalImageComplex(img_path, height, width, 1, gt, projection, false, true); + return echodata; + +} + +int ResampleGDAL(const char* pszSrcFile, const char* pszOutFile, double* gt, int new_width, + int new_height, GDALResampleAlg eResample) +{ + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO"); + GDALDataset* pDSrc = (GDALDataset*)GDALOpen(pszSrcFile, GA_ReadOnly); + if(pDSrc == NULL) { + return -1; + } + + GDALDriver* pDriver = GetGDALDriverManager()->GetDriverByName("GTiff"); + if(pDriver == NULL) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + return -2; + } + int width = pDSrc->GetRasterXSize(); + int height = pDSrc->GetRasterYSize(); + int nBandCount = pDSrc->GetRasterCount(); + GDALDataType dataType = pDSrc->GetRasterBand(1)->GetRasterDataType(); + + char* pszSrcWKT = NULL; + pszSrcWKT = const_cast(pDSrc->GetProjectionRef()); + + // 锟斤拷锟矫伙拷锟酵队帮拷锟斤拷锟轿?拷锟斤拷锟揭伙拷锟?1锟?7 + if(strlen(pszSrcWKT) <= 0) { + OGRSpatialReference oSRS; + oSRS.importFromEPSG(4326); + // oSRS.SetUTM(50, true); //锟斤拷锟斤拷锟斤拷 锟斤拷锟斤拷120锟斤拷 + // oSRS.SetWellKnownGeogCS("WGS84"); + oSRS.exportToWkt(&pszSrcWKT); + } + qDebug() << "GDALCreateGenImgProjTransformer " << Qt::endl; + void* hTransformArg; + hTransformArg = GDALCreateGenImgProjTransformer((GDALDatasetH)pDSrc, pszSrcWKT, NULL, pszSrcWKT, + FALSE, 0.0, 1); + qDebug() << "no proj " << Qt::endl; + //(没锟斤拷投影锟斤拷影锟斤拷锟斤拷锟斤拷锟斤拷卟锟酵?拷锟?1锟?7) + if(hTransformArg == NULL) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + return -3; + } + qDebug() << "has proj " << Qt::endl; + double dGeoTrans[6] = { 0 }; + int nNewWidth = 0, nNewHeight = 0; + if(GDALSuggestedWarpOutput((GDALDatasetH)pDSrc, GDALGenImgProjTransform, hTransformArg, + dGeoTrans, &nNewWidth, &nNewHeight) + != CE_None) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + return -3; + } + + // GDALDestroyGenImgProjTransformer(hTransformArg); + + GDALDataset* pDDst = + pDriver->Create(pszOutFile, new_width, new_height, nBandCount, dataType, NULL); + if(pDDst == NULL) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + return -2; + } + + pDDst->SetProjection(pszSrcWKT); + pDDst->SetGeoTransform(gt); + + GDALWarpOptions* psWo = GDALCreateWarpOptions(); + + // psWo->papszWarpOptions = CSLDuplicate(NULL); + psWo->eWorkingDataType = dataType; + psWo->eResampleAlg = eResample; + + psWo->hSrcDS = (GDALDatasetH)pDSrc; + psWo->hDstDS = (GDALDatasetH)pDDst; + qDebug() << "GDALCreateGenImgProjTransformer" << Qt::endl; + psWo->pfnTransformer = GDALGenImgProjTransform; + psWo->pTransformerArg = GDALCreateGenImgProjTransformer( + (GDALDatasetH)pDSrc, pszSrcWKT, (GDALDatasetH)pDDst, pszSrcWKT, FALSE, 0.0, 1); + ; + + qDebug() << "GDALCreateGenImgProjTransformer has created" << Qt::endl; + psWo->nBandCount = nBandCount; + psWo->panSrcBands = (int*)CPLMalloc(nBandCount * sizeof(int)); + psWo->panDstBands = (int*)CPLMalloc(nBandCount * sizeof(int)); + for(int i = 0; i < nBandCount; i++) { + psWo->panSrcBands[i] = i + 1; + psWo->panDstBands[i] = i + 1; + } + + GDALWarpOperation oWo; + if(oWo.Initialize(psWo) != CE_None) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + GDALClose((GDALDatasetH)(GDALDatasetH)pDDst); + return -3; + } + qDebug() << "ChunkAndWarpImage:" << new_width << "," << new_height << Qt::endl; + oWo.ChunkAndWarpMulti(0, 0, new_width, new_height); + GDALFlushCache(pDDst); + qDebug() << "ChunkAndWarpImage over" << Qt::endl; + // GDALDestroyGenImgProjTransformer(psWo->pTransformerArg); + // GDALDestroyWarpOptions(psWo); + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + GDALClose((GDALDatasetH)(GDALDatasetH)pDDst); + ////GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + return 0; +} + +int ResampleGDALs(const char* pszSrcFile, int band_ids, GDALRIOResampleAlg eResample) +{ + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO"); + GDALDataset* pDSrc = (GDALDataset*)GDALOpen(pszSrcFile, GA_Update); + if(pDSrc == NULL) { + return -1; + } + + GDALDataType gdal_datatype = pDSrc->GetRasterBand(1)->GetRasterDataType(); + + GDALRasterBand* demBand = pDSrc->GetRasterBand(band_ids); + + int width = pDSrc->GetRasterXSize(); + int height = pDSrc->GetRasterYSize(); + int start_col = 0, start_row = 0, rows_count = 0, cols_count; + + int row_delta = int(120000000 / width); + + GDALRasterIOExtraArg psExtraArg; + INIT_RASTERIO_EXTRA_ARG(psExtraArg); + psExtraArg.eResampleAlg = eResample; + + do { + rows_count = start_row + row_delta < height ? row_delta : height - start_row; + cols_count = width; + + if(gdal_datatype == GDALDataType::GDT_UInt16) { + unsigned short* temp = new unsigned short[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, + cols_count, rows_count, gdal_datatype, 0, 0); + demBand->RasterIO(GF_Write, start_col, start_row, cols_count, rows_count, temp, + cols_count, rows_count, gdal_datatype, 0, 0, &psExtraArg); + delete[] temp; + + } else if(gdal_datatype == GDALDataType::GDT_Int16) { + short* temp = new short[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, + cols_count, rows_count, gdal_datatype, 0, 0); + demBand->RasterIO(GF_Write, start_col, start_row, cols_count, rows_count, temp, + cols_count, rows_count, gdal_datatype, 0, 0, &psExtraArg); + delete[] temp; + } else if(gdal_datatype == GDALDataType::GDT_Float32) { + float* temp = new float[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, + cols_count, rows_count, gdal_datatype, 0, 0); + demBand->RasterIO(GF_Write, start_col, start_row, cols_count, rows_count, temp, + cols_count, rows_count, gdal_datatype, 0, 0, &psExtraArg); + delete[] temp; + } + start_row = start_row + rows_count; + } while(start_row < height); + GDALClose((GDALDatasetH)pDSrc); + + return 0; +} + +int alignRaster(QString inputPath, QString referencePath, QString outputPath, GDALResampleAlg eResample) +{ + + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + + GDALDataset* pDSrc = (GDALDataset*)GDALOpen(inputPath.toUtf8().constData(), GA_ReadOnly); + if (pDSrc == NULL) { + return -1; + } + + GDALDataset* tDSrc = (GDALDataset*)GDALOpen(referencePath.toUtf8().constData(), GA_ReadOnly); + if (tDSrc == NULL) { + return -1; + } + + long new_width = tDSrc->GetRasterXSize(); + long new_height = tDSrc->GetRasterYSize(); + + GDALDriver* pDriver = GetGDALDriverManager()->GetDriverByName("GTiff"); + if (pDriver == NULL) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + return -2; + } + int width = pDSrc->GetRasterXSize(); + int height = pDSrc->GetRasterYSize(); + int nBandCount = pDSrc->GetRasterCount(); + GDALDataType dataType = pDSrc->GetRasterBand(1)->GetRasterDataType(); + + char* pszSrcWKT = NULL; + pszSrcWKT = const_cast(pDSrc->GetProjectionRef()); + + if (strlen(pszSrcWKT) <= 0) { + OGRSpatialReference oSRS; + oSRS.importFromEPSG(4326); + oSRS.exportToWkt(&pszSrcWKT); + } + qDebug() << "GDALCreateGenImgProjTransformer " << Qt::endl; + void* hTransformArg; + hTransformArg = GDALCreateGenImgProjTransformer((GDALDatasetH)pDSrc, pszSrcWKT, NULL, pszSrcWKT,FALSE, 0.0, 1); + qDebug() << "no proj " << Qt::endl; + if (hTransformArg == NULL) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + return -3; + } + qDebug() << "has proj " << Qt::endl; + std::shared_ptr dGeoTrans(new double[6], delArrPtr); + int nNewWidth = 0, + nNewHeight = 0; + + if (GDALSuggestedWarpOutput((GDALDatasetH)pDSrc, GDALGenImgProjTransform, hTransformArg,dGeoTrans.get(), &nNewWidth, &nNewHeight) != CE_None) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + return -3; + } + + GDALDataset* pDDst = pDriver->Create(outputPath.toUtf8().constData(), new_width, new_height, pDSrc->GetRasterCount(), dataType, NULL); + if (pDDst == NULL) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + return -2; + } + + std::shared_ptr gt(new double[6], delArrPtr); + tDSrc->GetGeoTransform(gt.get()); + pDDst->SetProjection(pszSrcWKT); + pDDst->SetGeoTransform(gt.get()); + + GDALWarpOptions* psWo = GDALCreateWarpOptions(); + + psWo->eWorkingDataType = dataType; + psWo->eResampleAlg = eResample; + + psWo->hSrcDS = (GDALDatasetH)pDSrc; + psWo->hDstDS = (GDALDatasetH)pDDst; + + psWo->pfnTransformer = GDALGenImgProjTransform; + psWo->pTransformerArg = GDALCreateGenImgProjTransformer((GDALDatasetH)pDSrc, pszSrcWKT, (GDALDatasetH)pDDst, pszSrcWKT, FALSE, 0.0, 1); + + + qDebug() << "GDALCreateGenImgProjTransformer has created" << Qt::endl; + psWo->nBandCount = nBandCount; + psWo->panSrcBands = (int*)CPLMalloc(nBandCount * sizeof(int)); + psWo->panDstBands = (int*)CPLMalloc(nBandCount * sizeof(int)); + for (int i = 0; i < nBandCount; i++) { + psWo->panSrcBands[i] = i + 1; + psWo->panDstBands[i] = i + 1; + } + + GDALWarpOperation oWo; + if (oWo.Initialize(psWo) != CE_None) { + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + GDALClose((GDALDatasetH)(GDALDatasetH)pDDst); + GDALClose((GDALDatasetH)(GDALDatasetH)tDSrc); + return -3; + } + qDebug() << "ChunkAndWarpImage:" << new_width << "," << new_height << Qt::endl; + oWo.ChunkAndWarpMulti(0, 0, new_width, new_height); + GDALFlushCache(pDDst); + qDebug() << "ChunkAndWarpImage over" << Qt::endl; + GDALClose((GDALDatasetH)(GDALDatasetH)pDSrc); + GDALClose((GDALDatasetH)(GDALDatasetH)pDDst); + GDALClose((GDALDatasetH)(GDALDatasetH)tDSrc); + return 0; + + +} + + + + +int saveMatrixXcd2TiFF(Eigen::MatrixXcd data, QString out_tiff_path) +{ + int rows = data.rows(); + int cols = data.cols(); + + Eigen::MatrixXd gt = Eigen::MatrixXd::Zero(2, 3); + + gdalImage image_tiff = + CreategdalImage(out_tiff_path, rows, cols, 2, gt, "", false, true); // 注意这里保留仿真结果 + // 保存二进制文件 + Eigen::MatrixXd real_img = data.array().real(); + Eigen::MatrixXd imag_img = data.array().imag(); + image_tiff.saveImage(real_img, 0, 0, 1); + image_tiff.saveImage(imag_img, 0, 0, 2); + return -1; +} + +gdalImageComplex::gdalImageComplex(const QString& raster_path) +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + this->img_path = raster_path; + + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDataset* rasterDataset = (GDALDataset*)(GDALOpen( + raster_path.toUtf8().constData(), GA_ReadOnly)); // 锟斤拷只斤拷式锟斤拷取斤拷影锟斤拷 + this->width = rasterDataset->GetRasterXSize(); + this->height = rasterDataset->GetRasterYSize(); + this->band_num = rasterDataset->GetRasterCount(); + + double* gt = new double[6]; + rasterDataset->GetGeoTransform(gt); + this->gt = Eigen::MatrixXd(2, 3); + this->gt << gt[0], gt[1], gt[2], gt[3], gt[4], gt[5]; + + double a = this->gt(0, 0); + double b = this->gt(0, 1); + double c = this->gt(0, 2); + double d = this->gt(1, 0); + double e = this->gt(1, 1); + double f = this->gt(1, 2); + + this->projection = rasterDataset->GetProjectionRef(); + + // 斤拷投影 + GDALFlushCache((GDALDatasetH)rasterDataset); + GDALClose((GDALDatasetH)rasterDataset); + rasterDataset = NULL; // 指矫匡拷 + this->InitInv_gt(); + delete[] gt; + ////GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + omp_unset_lock(&lock); // 锟酵放伙拷斤拷 + omp_destroy_lock(&lock); // 劫伙拷斤拷 +} + +gdalImageComplex::~gdalImageComplex() {} + +void gdalImageComplex::setData(Eigen::MatrixXcd data) +{ + this->data = data; +} + +void gdalImageComplex::saveImage(Eigen::MatrixXcd data, int start_row, int start_col, int band_ids) +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + if(start_row + data.rows() > this->height || start_col + data.cols() > this->width) { + QString tip = u8"file path: " + this->img_path; + qDebug() << tip; + throw std::exception(tip.toUtf8().constData()); + } + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + GDALDataset* poDstDS = nullptr; + if(exists_test(this->img_path)) { + poDstDS = (GDALDataset*)(GDALOpen(this->img_path.toUtf8().constData(), GA_Update)); + } else { + poDstDS = poDriver->Create(this->img_path.toUtf8().constData(), this->width, this->height, + this->band_num, GDT_CFloat64, NULL); // 斤拷锟斤拷 + poDstDS->SetProjection(this->projection.toUtf8().constData()); + + double gt_ptr[6]; + for(int i = 0; i < this->gt.rows(); i++) { + for(int j = 0; j < this->gt.cols(); j++) { + gt_ptr[i * 3 + j] = this->gt(i, j); + } + } + poDstDS->SetGeoTransform(gt_ptr); + //delete[] gt_ptr; + } + + int datarows = data.rows(); + int datacols = data.cols(); + + double* databuffer = new double[data.size() * 2]; + for(int i = 0; i < data.rows(); i++) { + for(int j = 0; j < data.cols(); j++) { + databuffer[i * data.cols() * 2 + j * 2] = data(i, j).real(); + databuffer[i * data.cols() * 2 + j * 2 + 1] = data(i, j).imag(); + } + } + + // poDstDS->RasterIO(GF_Write,start_col, start_row, datacols, datarows, databuffer, datacols, + // datarows, GDT_Float32,band_ids, num,0,0,0); + poDstDS->GetRasterBand(band_ids)->RasterIO(GF_Write, start_col, start_row, datacols, datarows, + databuffer, datacols, datarows, GDT_CFloat64, 0, 0); + GDALFlushCache(poDstDS); + delete databuffer; + GDALClose((GDALDatasetH)poDstDS); + // GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // +} + +Eigen::MatrixXcd gdalImageComplex::getDataComplex(int start_row, int start_col, int rows_count, + int cols_count, int band_ids) +{ + GDALDataset* poDataset; + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + + // 打开TIFF文件 + poDataset = (GDALDataset*)GDALOpen(this->img_path.toUtf8().constData(), GA_ReadOnly); + if(poDataset == nullptr) { + QMessageBox::warning(nullptr, u8"错误", u8"无法打开文件：" + this->img_path); + qDebug() << u8"无法打开文件：" + this->img_path; + } + + // 获取数据集的第一个波段 + GDALRasterBand* poBand; + poBand = poDataset->GetRasterBand(1); + + // 读取波段信息，假设是复数类型 + int nXSize = poBand->GetXSize(); + int nYSize = poBand->GetYSize(); + + double* databuffer = new double[nXSize * nYSize * 2]; + poBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, databuffer, cols_count, + rows_count, GDT_CFloat64, 0, 0); + GDALClose((GDALDatasetH)poDataset); + Eigen::MatrixXcd rasterData(nYSize, nXSize); // 使用Eigen的MatrixXcd + for(size_t i = 0; i < nYSize; i++) { + for(size_t j = 0; j < nXSize; j++) { + rasterData(i, j) = std::complex(databuffer[i * nXSize * 2 + j * 2], + databuffer[i * nXSize * 2 + j * 2 + 1]); + } + } + + delete[] databuffer; + return rasterData; +} + +void gdalImageComplex::saveImage() +{ + this->saveImage(this->data, this->start_row, this->start_col, this->data_band_ids); +} +void gdalImageComplex::savePreViewImage() +{ + qDebug()<<"void gdalImageComplex::savePreViewImage()"; + Eigen::MatrixXd data_abs = Eigen::MatrixXd::Zero(this->height, this->width); + data_abs = (this->data.array().real().pow(2) + this->data.array().imag().pow(2)) + .array() + .log10()*10.0; // 计算振幅 + + double min_abs = data_abs.minCoeff(); // 最大值 + double max_abs = data_abs.maxCoeff(); // 最小值 + double delta = (max_abs - min_abs) / 1000; // 1000分位档 + Eigen::MatrixX data_idx = + ((data_abs.array() - min_abs).array() / delta).array().floor().cast(); + + // 初始化 + double hist[1001]; + for(size_t i = 0; i < 1001; i++) { + hist[i] = 0; // 初始化 + } + for(size_t i = 0; i < this->height; i++) { + for(size_t j = 0; j < this->width; j++) { + hist[data_idx(i, j)]++; + } + } + + // 统计 + size_t count = this->height * this->width; + double precent = 0; + size_t curCount = 0; + double pre2 = 0; + bool findprec_2 = true; + double pre98 = 0; + bool findprec_98 = true; + for(size_t i = 0; i < 1001; i++) { + precent = precent + hist[i]; + if(findprec_2 & precent / count > 0.02) { + pre2 = i * delta + min_abs; + findprec_2 = false; + } + if(findprec_98 & precent / count > 0.98) { + pre98 = (i-1) * delta + min_abs; + findprec_98 = false; + } + } + // 拉伸 + delta = (pre98-pre2)/200; + data_idx= + ((data_abs.array() - pre2).array() / delta).array().floor().cast(); + + for(size_t i = 0; i < this->height; i++) { + for(size_t j = 0; j < this->width; j++) { + if(data_idx(i,j)<0){ + data_idx(i,j)=0; + } + else if(data_idx(i,j)>255){ + data_idx(i,j)=255; + }else{ + + } + } + } + + // 赋值 + QString filePath = this->img_path; + QFile file(filePath); + QFileInfo fileInfo(file); + + QString directory = fileInfo.absolutePath(); + qDebug() << "文件所在目录：" << directory; + QString baseName = fileInfo.completeBaseName(); + qDebug() << "无后缀文件名：" << baseName; + + // 创建文件路径 + QString previewImagePath = JoinPath(directory, baseName+"_preview.png"); + cv::Mat img(this->height, this->width, CV_8U ,cv::Scalar(0)); + + for(size_t i = 0; i < this->height; i++) { + for(size_t j = 0; j < this->width; j++) { + img.at(i,j)= (uchar)(data_idx(i,j)); + } + } + //std::string outimgpath=previewImagePath.toUtf8().data(); + cv::imwrite(previewImagePath.toUtf8().data(), img); +} + + + + + + + + diff --git a/ImageOperatorBase.h b/ImageOperatorBase.h new file mode 100644 index 0000000..7d87945 --- /dev/null +++ b/ImageOperatorBase.h @@ -0,0 +1,352 @@ +#pragma once +/** + * 影像数据操作项 + * 所有数据相关的操作，都是基于ENVI的数据格式 + * 影像读写操作时基于 GDAL 库 + * **/ + + + +#ifndef IMAGEOPERATORBASE_H +#define IMAGEOPERATORBASE_H +#include "BaseConstVariable.h" +#include "GeoOperator.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include // for CPLMalloc() + + +struct ImageGEOINFO { + int width; + int height; + int bandnum; +}; + + +enum GDALREADARRCOPYMETHOD { + MEMCPYMETHOD, // 直接拷贝 + VARIABLEMETHOD // 变量赋值 + +}; + + +// 判断是否需要输出为DLL +#define DLLOUT +// 文件打开 // 当指令销毁时，调用GDALClose 销毁类型 +std::shared_ptr OpenDataset(const QString& in_path, GDALAccess rwmode= GA_ReadOnly); +void CloseDataset(GDALDataset* ptr); + +// 数据格式转换 +int TIFF2ENVI(QString in_tiff_path,QString out_envi_path); +int ENVI2TIFF(QString in_envi_path,QString out_tiff_path); + +// 保存影像数据 --直接保存 ENVI 文件 + +int CreateDataset(QString new_file_path, int height, int width, int band_num, double* gt, QString projection, GDALDataType gdal_dtype, bool need_gt); // 创建文件 + +int saveDataset(QString new_file_path, int start_line, int start_cols, int band_ids, int datacols, int datarows, void* databuffer); + +// 根据限制条件估算分块大小 +int block_num_pre_memory(int width, int height, GDALDataType gdal_dtype,double memey_size); + +// 将结果转换为复数 或者 实数 +Eigen::Matrix ReadComplexMatrixData(int start_line,int width, int line_num, std::shared_ptr rasterDataset, GDALDataType gdal_datatype); + +Eigen::Matrix ReadMatrixDoubleData(int start_line, int width, int line_num, std::shared_ptr rasterDataset, GDALDataType gdal_datatype,int band_idx); + +Eigen::MatrixXd getGeoTranslationArray(QString in_path); +ImageGEOINFO getImageINFO(QString in_path); + +GDALDataType getGDALDataType(QString fileptah); + + +struct DemBox { + double min_lat; //纬度 + double min_lon;//经度 + double max_lat;//纬度 + double max_lon;//经度 +}; + + + + + +/// +/// gdalImage图像操作类 +/// + +class gdalImage +{ + +public: // 方法 + gdalImage(); + gdalImage(const QString& raster_path); + ~gdalImage(); + virtual void setHeight(int); + virtual void setWidth(int); + virtual void setTranslationMatrix(Eigen::MatrixXd gt); + virtual void setData(Eigen::MatrixXd,int data_band_ids=1); + virtual Eigen::MatrixXd getData(int start_row, int start_col, int rows_count, int cols_count, int band_ids); + virtual Eigen::MatrixXd getGeoTranslation(); + virtual GDALDataType getDataType(); + virtual void saveImage(Eigen::MatrixXd, int start_row, int start_col, int band_ids); + virtual void saveImage(); + virtual void setNoDataValue(double nodatavalue, int band_ids); + virtual int InitInv_gt(); + virtual Landpoint getRow_Col(double lon, double lat); + virtual Landpoint getLandPoint(double i, double j, double ati); + + virtual void getLandPoint(double i, double j, double ati, Landpoint& Lp); + + virtual double mean(int bandids = 1); + double BandmaxValue(int bandids = 1); + double BandminValue(int bandids = 1); + virtual GDALRPCInfo getRPC(); + virtual Eigen::MatrixXd getLandPoint(Eigen::MatrixXd points); + + virtual Eigen::MatrixXd getHist(int bandids); + + +public: + QString img_path; // 图像文件 + int height; // 高 + int width; // 宽 + int band_num;// 波段数 + int start_row;// + int start_col;// + int data_band_ids; + Eigen::MatrixXd gt; // 变换矩阵 + Eigen::MatrixXd inv_gt; // 逆变换矩阵 + Eigen::MatrixXd data; + QString projection; +}; + + + +/// +/// gdalImage图像操作类 +/// +class gdalImageComplex:public gdalImage +{ + +public: // 方法 + gdalImageComplex(const QString& raster_path); + ~gdalImageComplex(); + void setData(Eigen::MatrixXcd); + void saveImage(Eigen::MatrixXcd data, int start_row, int start_col, int band_ids); + Eigen::MatrixXcd getDataComplex(int start_row, int start_col, int rows_count, int cols_count, int band_ids); + + void saveImage() override; + void savePreViewImage(); +public: + Eigen::MatrixXcd data; +}; + +// 创建影像 +gdalImage CreategdalImage(const QString& img_path, int height, int width, int band_num, Eigen::MatrixXd gt, QString projection, bool need_gt = true, bool overwrite = false); + +gdalImageComplex CreategdalImageComplex(const QString& img_path, int height, int width, int band_num, Eigen::MatrixXd gt, QString projection, bool need_gt = true, bool overwrite = false); + +gdalImageComplex CreateEchoComplex(const QString& img_path, int height, int width, int band_num); + + + +int ResampleGDAL(const char* pszSrcFile, const char* pszOutFile, double* gt, int new_width, int new_height, GDALResampleAlg eResample); + +int ResampleGDALs(const char* pszSrcFile, int band_ids, GDALRIOResampleAlg eResample = GRIORA_Bilinear); + +int alignRaster(QString inputPath, QString referencePath, QString outputPath, GDALResampleAlg eResample); + +//--------------------- 保存文博 ------------------------------- + +int saveMatrixXcd2TiFF(Eigen::MatrixXcd data, QString out_tiff_path); + +//---------------------------------------------------- + + +template +std::shared_ptr readDataArr(gdalImage &img,int start_row, int start_col, int rows_count, int cols_count, int band_ids, GDALREADARRCOPYMETHOD method); + + +#ifndef DLLOUT + + + + + +#else +//#define DllExport __declspec( dllexport ) +//double __declspec(dllexport) ProcessMGCMathX_MGC(int Xbetaidx, int Xbwidx, double XTao, double satH, char* sigma_path, char* output_path, +// double p1_x, double p1_y, double p2_x, double p2_y, double p3_x, double p3_y, double p4_x, double p4_y) + +#endif + +#endif + +template +std::shared_ptr readDataArr(gdalImage& imgds, int start_row, int start_col, int rows_count, int cols_count, int band_ids, GDALREADARRCOPYMETHOD method) +{ + std::shared_ptr result = nullptr; + + + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDataset* rasterDataset = (GDALDataset*)(GDALOpen(imgds.img_path.toUtf8().constData(), GA_ReadOnly)); // 锟斤拷只斤拷式锟斤拷取斤拷影锟斤拷 + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* demBand = rasterDataset->GetRasterBand(band_ids); + + rows_count = start_row + rows_count <= imgds.height ? rows_count : imgds.height - start_row; + cols_count = start_col + cols_count <= imgds.width ? cols_count : imgds.width - start_col; + + Eigen::MatrixXd datamatrix(rows_count, cols_count); + + if (gdal_datatype == GDT_Byte) { + char* temp = new char[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); + result=std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { + std::memcpy(result.get(), temp, rows_count * cols_count); + } + else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { + long count = rows_count * cols_count; + for (long i = 0; i < count; i++) { + result.get()[i] = T(temp[i]); + } + } + + + delete[] temp; + } + else if (gdal_datatype == GDT_UInt16) { + unsigned short* temp = new unsigned short[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); + result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { + std::memcpy(result.get(), temp, rows_count * cols_count); + } + else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { + long count = rows_count * cols_count; + for (long i = 0; i < count; i++) { + result.get()[i] = T(temp[i]); + } + } + delete[] temp; + } + else if (gdal_datatype == GDT_Int16) { + short* temp = new short[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); + result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { + std::memcpy(result.get(), temp, rows_count * cols_count); + } + else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { + long count = rows_count * cols_count; + for (long i = 0; i < count; i++) { + result.get()[i] = T(temp[i]); + } + } + delete[] temp; + } + else if (gdal_datatype == GDT_UInt32) { + unsigned int* temp = new unsigned int[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); + result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { + std::memcpy(result.get(), temp, rows_count * cols_count); + } + else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { + long count = rows_count * cols_count; + for (long i = 0; i < count; i++) { + result.get()[i] = T(temp[i]); + } + } + delete[] temp; + } + else if (gdal_datatype == GDT_Int32) { + int* temp = new int[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); + result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { + std::memcpy(result.get(), temp, rows_count * cols_count); + } + else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { + long count = rows_count * cols_count; + for (long i = 0; i < count; i++) { + result.get()[i] = T(temp[i]); + } + } + delete[] temp; + } + // else if (gdal_datatype == GDT_UInt64) { + // unsigned long* temp = new unsigned long[rows_count * cols_count]; + // demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + //rows_count, gdal_datatype, 0, 0); for (int i = 0; i < rows_count; i++) { for (int j = 0; j < + //cols_count; j++) { datamatrix(i, j) = temp[i * cols_count + j]; + // } + // } + // delete[] temp; + // } + // else if (gdal_datatype == GDT_Int64) { + // long* temp = new long[rows_count * cols_count]; + // demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, + //rows_count, gdal_datatype, 0, 0); for (int i = 0; i < rows_count; i++) { for (int j = 0; j < + //cols_count; j++) { datamatrix(i, j) = temp[i * cols_count + j]; + // } + // } + // delete[] temp; + // } + else if (gdal_datatype == GDT_Float32) { + float* temp = new float[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); + + result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { + std::memcpy(result.get(), temp, rows_count * cols_count); + } + else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { + long count = rows_count * cols_count; + for (long i = 0; i < count; i++) { + result.get()[i] = T(temp[i]); + } + } + delete[] temp; + } + else if (gdal_datatype == GDT_Float64) { + double* temp = new double[rows_count * cols_count]; + demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); + + result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { + std::memcpy(result.get(), temp, rows_count * cols_count); + } + else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { + long count = rows_count * cols_count; + for (long i = 0; i < count; i++) { + result.get()[i] = T(temp[i]); + } + } + delete[] temp; + } + else { + } + GDALClose((GDALDatasetH)rasterDataset); + omp_unset_lock(&lock); // 锟酵放伙拷斤拷 + omp_destroy_lock(&lock); // 劫伙拷斤拷 + // GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + return result; +} + diff --git a/ImageShowDialogClass.cpp b/ImageShowDialogClass.cpp new file mode 100644 index 0000000..8206c12 --- /dev/null +++ b/ImageShowDialogClass.cpp @@ -0,0 +1,495 @@ +#include "stdafx.h" +#include "ImageShowDialogClass.h" +#include "ui_ImageShowDialogClass.h" +ImageShowDialogClass ::ImageShowDialogClass(QWidget *parent) + : QDialog(parent) +{ + ui=new Ui::ImageShowDialogClass; + ui->setupUi(this); + ui->m_plot->setInteractions(QCP::iRangeDrag | QCP::iRangeZoom); + connect(this->ui->m_plot, SIGNAL(mouseMove(QMouseEvent*)), this, SLOT(updateCursor(QMouseEvent*))); + this->graphclass=LAMPDATASHOWCLASS::NOWINDOWS; + + this->menubar = new QMenuBar(this); + QMenu* windowsMenu = this->menubar->addMenu(u8"数据"); + QAction* action_cursor_enable = windowsMenu->addAction(u8"打开游标"); + QObject::connect(action_cursor_enable,SIGNAL(triggered()),this,SLOT(on_action_cursor_enable_trigged())); + this->tracer=new QCPItemTracer(this->ui->m_plot); + this->m_plot = ui->m_plot; + + this->desCursor = nullptr; + this->HlineCursor = nullptr; + this->VlineCursor = nullptr; + + this->desCursorflag=false; + this->HlineCursorflag=false; + this->VlineCursorflag = false; + + this->statusbar=new QStatusBar(this); + this->statusbar->setMaximumHeight(20); + this->statusbar->setMinimumHeight(20); + + ui->verticalLayout->setMenuBar(this->menubar); + ui->verticalLayout->addWidget(this->statusbar); + this->setLayout(ui->verticalLayout); +} + + +ImageShowDialogClass::~ImageShowDialogClass() +{} + + +void ImageShowDialogClass::on_action_cursor_enable_trigged() +{ + this->desCursor = new ImageShowCursorDesClass(this); + + this->desCursorflag = true; + for (size_t i = 0; i < this->getGraphCount(); i++) { + if (this->getGraphClass(i) == LAMPDATASHOWCLASS::LAMPColorMap) { + this->HlineCursorflag = true; + this->VlineCursorflag = true; + } + } + + //下面的代码就是设置游标的外观 + this->setMouseTracking(true); + tracer->setInterpolating(false);//禁用插值 + tracer->setPen(QPen(Qt::DashLine));//虚线游标 + tracer->setStyle(QCPItemTracer::tsPlus);//游标样式：十字星、圆圈、方框 + tracer->setBrush(Qt::red);//游标颜色 + + if (this->graphclass == LAMPDATASHOWCLASS::LAMPColorMap) { + this->HlineCursor = new ImageShowCursorLineClass(this); + this->VlineCursor = new ImageShowCursorLineClass(this); + QObject::connect(this->HlineCursor, SIGNAL(windowsClose()), this, SLOT(on_action_Hlinecursor_close_trigged())); + QObject::connect(this->VlineCursor, SIGNAL(windowsClose()), this, SLOT(on_action_VVlinecursor_close_trigged())); + QObject::connect(this->ui->m_plot->xAxis, SIGNAL(rangeChanged(QCPRange)), this->HlineCursor, SLOT(xAxisRangeChanged(QCPRange))); + QObject::connect(this->ui->m_plot->yAxis, SIGNAL(rangeChanged(QCPRange)), this->HlineCursor, SLOT(xAxisRangeChanged(QCPRange))); + this->HlineCursor->show(); + this->VlineCursor->show(); + } + this->desCursor->show(); + QObject::connect(this->desCursor, SIGNAL(windowsClose()), this, SLOT(on_action_descursor_close_trigged())); +} + +void ImageShowDialogClass::PlotLine(QVector xs, QVector dataValues,QString Name) +{ + if (dataValues.count() == 0) { return; } + QCPGraph* graph = ui->m_plot->addGraph(); + graph->setName(Name); + graph->setData(xs, dataValues); + double min_y = dataValues[0]; + double max_y = dataValues[0]; + for (size_t i = 0; i < dataValues.count(); i++) { + min_y = min_y > dataValues[i] ? dataValues[i] : min_y; + max_y = max_y < dataValues[i] ? dataValues[i] : max_y; + } + + double min_x = xs[0]; + double max_x = xs[0]; + for (size_t i = 0; i < xs.count(); i++) { + min_x = min_x > xs[i] ? xs[i] : min_x; + max_x = max_x < xs[i] ? xs[i] : max_x; + } + ui->m_plot->xAxis->setRange(min_x, max_x); + ui->m_plot->yAxis->setRange(min_y, max_y); + ui->m_plot->legend->setVisible(true); + ui->m_plot->replot(); + + + //下面的代码就是设置游标的外观 + this->setMouseTracking(true); + tracer->setInterpolating(false);//禁用插值 + tracer->setPen(QPen(Qt::DashLine));//虚线游标 + tracer->setStyle(QCPItemTracer::tsPlus);//游标样式：十字星、圆圈、方框 + tracer->setBrush(Qt::red);//游标颜色 + tracer->setGraph(graph); // 确保游标跟随 + this->m_graphMap.insert(Name, graph); + +} + +void ImageShowDialogClass::Scatter(double dx,double dy) +{ + QVectorxs(0); + QVectorys(0); + xs.push_back(dx); + ys.push_back(dy); + // 设置画笔风格 + QPen drawPen; + drawPen.setColor(Qt::red); + drawPen.setWidth(4); + + // 绘制散点 + QCPGraph* curGraph = ui->m_plot->addGraph(); + curGraph->setPen(drawPen); + curGraph->setLineStyle(QCPGraph::lsNone); + curGraph->setScatterStyle(QCPScatterStyle(QCPScatterStyle::ssCircle, 2)); + curGraph->setData(xs, ys); + + ui->m_plot->xAxis->setVisible(true); + ui->m_plot->yAxis->setVisible(true); + +} + +void ImageShowDialogClass::Scatter(QVector xs, QVector ys) +{ + // 设置画笔风格 + QPen drawPen; + drawPen.setColor(Qt::red); + drawPen.setWidth(4); + + // 绘制散点 + QCPGraph* curGraph = ui->m_plot->addGraph(); + curGraph->setPen(drawPen); + curGraph->setLineStyle(QCPGraph::lsNone); + curGraph->setScatterStyle(QCPScatterStyle(QCPScatterStyle::ssCircle, 2)); + curGraph->setData(xs, ys); + + ui->m_plot->xAxis->setVisible(true); + ui->m_plot->yAxis->setVisible(true); + +} + +void ImageShowDialogClass::load_double_MatrixX_data(Eigen::MatrixXd data, QString name) +{ + this->setWindowTitle(name); + this->graphclass = LAMPDATASHOWCLASS::LAMPColorMap; + int ny = data.rows(); // 行数 + int nx = data.cols(); // 列数 + // 创建 Color Map + QCPColorMap* colorMap = new QCPColorMap(ui->m_plot->xAxis, ui->m_plot->yAxis); + colorMap->data()->setSize(nx, ny); // 设置 Color Map 的大小 + colorMap->data()->setRange(QCPRange(0, nx), QCPRange(0, ny)); // 设置坐标轴的范围 + + // 填充数据 + for (int xIndex = 0; xIndex < nx; ++xIndex) { + for (int yIndex = 0; yIndex < ny; ++yIndex) { + double magnitude = data(yIndex, xIndex); // 或者使用 std::arg() 获取相位 + colorMap->data()->setCell(xIndex, yIndex,magnitude); + } + } + + colorMap->setGradient(QCPColorGradient::gpJet); + colorMap->rescaleDataRange(true); + ui->m_plot->rescaleAxes(); + ui->m_plot->replot(); + +} + +void ImageShowDialogClass::load_double_MatrixX_data(Eigen::MatrixXd X, Eigen::MatrixXd Y, Eigen::MatrixXd data, QString name) +{ + this->graphclass = LAMPDATASHOWCLASS::LAMPColorMap; + int nx = data.cols(); // 行数 + int ny = data.rows(); // 列数 + // 创建 Color Map + ui->m_plot->xAxis->setRange(X(0, 0), X(0, nx - 1)); + ui->m_plot->yAxis->setRange(Y(0, 0), Y(ny - 1, 0)); + QCPColorMap* colorMap = new QCPColorMap(ui->m_plot->xAxis, ui->m_plot->yAxis); + colorMap->data()->setSize(ny, nx); // 设置 Color Map 的大小 + colorMap->data()->setRange(QCPRange(X(0,0), X(0,nx-1)), QCPRange(Y(0,0),Y(ny-1,0))); // 设置坐标轴的范围 + // 填充数据 + for (int xIndex = 0; xIndex < nx; ++xIndex) { + for (int yIndex = 0; yIndex < ny; ++yIndex) { + double magnitude = data(yIndex, xIndex); // 或者使用 std::arg() 获取相位 + colorMap->data()->setCell(yIndex, xIndex, magnitude); + } + } + + // add a color scale: + QCPColorScale* m_colorScale = new QCPColorScale(m_plot); + m_plot->plotLayout()->addElement(0, 1, m_colorScale); + m_colorScale->setType(QCPAxis::atRight); + colorMap->setColorScale(m_colorScale); + + colorMap->setGradient(QCPColorGradient::gpJet); + colorMap->rescaleDataRange(true); + ui->m_plot->rescaleAxes(); + ui->m_plot->replot(); +} + +void ImageShowDialogClass::remove_Data(QString name) +{ + for(size_t i=0;im_plot->graphCount();i++){ + if (this->m_plot->graph(i)->name() == name) { + this->m_plot->removeGraph(i); + this->m_plot->replot(); + return; + } + } +} + + + +LAMPDATASHOWCLASS ImageShowDialogClass::getGraphClass(size_t i) +{ + if (this->m_plot->graphCount() == 0 || i > this->m_plot->graphCount()) + { + return LAMPDATASHOWCLASS::NOWINDOWS; + } + QCPAbstractPlottable* plottable = this->ui->m_plot->plottable(i); + if (dynamic_cast(plottable)) { + return LAMPDATASHOWCLASS::LAMPGraph; + } + else if (dynamic_cast(plottable)) { + return LAMPDATASHOWCLASS::LAMPScatterStyle; + } + else if (dynamic_cast(plottable)) { + return LAMPDATASHOWCLASS::LAMPColorMap; + } + else { + return LAMPDATASHOWCLASS::NOWINDOWS; + } +} + + +size_t ImageShowDialogClass::getGraphCount() +{ + return this->m_plot->graphCount(); +} + +void ImageShowDialogClass::on_action_descursor_close_trigged() +{ + this->desCursorflag = false; +} + +void ImageShowDialogClass::on_action_Hlinecursor_close_trigged() +{ + this->HlineCursorflag = false; +} + +void ImageShowDialogClass::on_action_VVlinecursor_close_trigged() +{ + this->VlineCursorflag = false; +} + +void ImageShowDialogClass::updateCursor(QMouseEvent *event) +{ + if (this->desCursorflag &&this->HlineCursorflag &&this->VlineCursorflag) { + //下面的代码就是设置游标的外观 + this->setMouseTracking(false); + } + else { + // 准备获取数据 + QPoint pos = event->pos(); + double x=this->m_plot->xAxis->pixelToCoord(pos.x()); // 将鼠标位置映射到图表坐标系中 + double y = this->m_plot->yAxis->pixelToCoord(pos.y()); + this->statusbar->showMessage(u8"X: "+QString::number(x,'f', 6)+" y: "+QString::number(y, 'f', 6)); + if (this->desCursorflag) { + if (this->graphclass == LAMPDATASHOWCLASS::LAMPColorMap) { + QCPColorMap* colorMap = dynamic_cast(this->ui->m_plot->plottable(0)); + double dataValue = colorMap->data()->data(x, y); + this->desCursor->updateCursorContent(u8"X: " + QString::number(x, 'f', 6) + " y: " + QString::number(y, 'f', 6) +u8"\n" +u8"DataValue: "+QString::number(dataValue)); + + + double xMin = this->m_plot->xAxis->range().lower; + double xMax = this->m_plot->xAxis->range().upper; + double yMin = this->m_plot->yAxis->range().lower; + double yMax = this->m_plot->yAxis->range().upper; + + long XCount = 1000; + double xdelte = (xMax - xMin) / (XCount-1); + double ydelte = (yMax - yMin) / (XCount-1); + + qDebug() << xMin << "," << xMax << " " << yMin << "," << yMax << " " << XCount<<" "< Hx(XCount); + QVector Hy(XCount); + QVector Vx(XCount); + QVector Vy(XCount); + + for (long i = 0; i < XCount; i++) + { + Hx[i] = i * xdelte + xMin; + Hy[i] = colorMap->data()->data(Hx[i], y); + Vx[i] = i * ydelte + yMin; + Vy[i] = colorMap->data()->data(x, Vx[i]); + } + + this->HlineCursor->UpdatePlotLine(Hx, Hy, u8"水平"); // 水平 + this->VlineCursor->UpdatePlotLine(Vx, Vy, u8"垂直"); // 垂直 + + } + else if (this->graphclass == LAMPDATASHOWCLASS::LAMPGraph) { + + } + else if (this->graphclass == LAMPDATASHOWCLASS::LAMPScatterStyle) { + + } + else { + this->desCursor->updateCursorContent(u8"无法识别图像类型"); + } + } + if (this->HlineCursorflag) { + + + } + if (this->VlineCursorflag) { + + } + } +} + + +ImageShowCursorDesClass::ImageShowCursorDesClass(QWidget* parent) +{ + // 创建 QLabel + label = new QLabel(this); + QFont font; + font.setPointSize(20); // 设置字体大小为20 + label->setFont(font); + // 创建布局 + QVBoxLayout* layout = new QVBoxLayout(this); + layout->addWidget(label); + + // 设置布局 + setLayout(layout); + this->resize(100, 50); +} + +ImageShowCursorDesClass::~ImageShowCursorDesClass() +{ +} + +void ImageShowCursorDesClass::closeEvent(QCloseEvent* event) +{ + QDialog::closeEvent(event); +} + + +void ImageShowCursorDesClass::updateCursorContent(QString content) { + + this->label->setText(content); + +} + + +ImageShowCursorLineClass::ImageShowCursorLineClass(QWidget* parent):QDialog(parent) +{ + this->menubar = new QMenuBar(this); + //QMenu* toolMenu = this->menubar->addMenu(u8"工具"); + //QAction* action_add_compare_line = toolMenu->addAction(u8"添加对比数据"); + //QObject::connect(action_add_compare_line, SIGNAL(triggered()), this, SLOT(load_new_compare_line())); + + this->tracerEnable = false; + this->tracer = nullptr; + this->tracerXText = nullptr; + this->tracerYText = nullptr; + // 创建 QLabel + this->customPlot = new QCustomPlot(this); + + // 创建布局 + QVBoxLayout* layout = new QVBoxLayout(this); + layout->addWidget(this->customPlot); + + // 设置布局 + setLayout(layout); + this->tracer = new QCPItemTracer(this->customPlot); + this->tracerYText = new QCPItemText(this->customPlot); + this->tracerXText = new QCPItemText(this->customPlot); + QObject::connect(this->customPlot, SIGNAL(mouseMove(QMouseEvent*)), this, SLOT(updateCursor(QMouseEvent*))); + this->customPlot->legend->setVisible(true); + this->setMinimumWidth(600); + this->setMinimumHeight(400); +} + +ImageShowCursorLineClass::~ImageShowCursorLineClass() +{ + +} + +void ImageShowCursorLineClass::loadPlotLine(QVector xs, QVector dataValues,QString Name) +{ + if (dataValues.count() == 0) { return; } + QCPGraph* graph = this->customPlot->addGraph(); + graph->setName(Name); + graph->setData(xs, dataValues); + double min_y = dataValues[0]; + double max_y = dataValues[0]; + for (size_t i = 0; i < dataValues.count(); i++) { + min_y = min_y > dataValues[i] ? dataValues[i] : min_y; + max_y = max_y < dataValues[i] ? dataValues[i] : max_y; + } + + double min_x = xs[0]; + double max_x = xs[0]; + for (size_t i = 0; i < xs.count(); i++) { + min_x = min_x > xs[i] ? xs[i] : min_x; + max_x = max_x < xs[i] ? xs[i] : max_x; + } + customPlot->xAxis->setRange(min_x, max_x); + customPlot->yAxis->setRange(min_y, max_y); + customPlot->legend->setVisible(true); + customPlot->replot(); + + + //下面的代码就是设置游标的外观 + this->setMouseTracking(true); + tracer->setInterpolating(false);//禁用插值 + tracer->setPen(QPen(Qt::DashLine));//虚线游标 + tracer->setStyle(QCPItemTracer::tsPlus);//游标样式：十字星、圆圈、方框 + tracer->setBrush(Qt::red);//游标颜色 + tracer->setGraph(graph); // 确保游标跟随 + this->m_graphMap.insert(Name, graph); + this->TrackName = Name; + +} + +void ImageShowCursorLineClass::UpdatePlotLine(QVector xs, QVector dataValues, QString Name) +{ + this->customPlot->clearGraphs(); + this->loadPlotLine(xs, dataValues, Name); +} + +void ImageShowCursorLineClass::updateCursor(QMouseEvent* event) +{ + if (tracerEnable)//游标使能 + { + + double x = this->customPlot->xAxis->pixelToCoord(event->pos().x());//鼠标点的像素坐标转plot坐标 + + tracer->setGraphKey(x);//设置游标的X值（这就是游标随动的关键代码） + double traceX = tracer->position->key(); + double traceY = tracer->position->value(); + + tracerXText->setText(QDateTime::fromMSecsSinceEpoch(traceX * 1000.0).toString("hh:mm:ss.zzz"));//游标文本框，指示游标的X值 + tracerYText->setText(QString::number(traceY));//游标文本框，指示游标的Y值 + + ////计算游标X值对应的所有曲线的Y值 + //for (int i = 0; i < this->customPlot->graphCount(); i++) + //{ + // QCPGraph* graph = dynamic_cast(this->customPlot->plottable(i)); + // //搜索左边离traceX最近的key对应的点，详情参考findBegin函数的帮助 + // QCPDataContainer::const_iterator coorPoint = graph->data().data()->findBegin(traceX, true);//true代表向左搜索 + // qDebug() << QString("graph%1对应的Y值是").arg(i) << coorPoint->value; + //} + } + + +} + +void ImageShowCursorLineClass::closeEvent(QCloseEvent* event) +{ + QDialog::closeEvent(event); +} + +void ImageShowCursorLineClass::xAxisRangeChanged(QCPRange range) { + this->customPlot->xAxis->setRange(range); + customPlot->replot(); +} + +void ImageShowCursorLineClass::yAxisRangeChanged(QCPRange range) { + this->customPlot->yAxis->setRange(range); + customPlot->replot(); +} + +void ImageShowCursorLineClass::on_SwichTracerGraph_Name() +{ +} + +void ImageShowCursorLineClass::load_new_compare_line() +{ + + + +} + \ No newline at end of file diff --git a/ImageShowDialogClass.h b/ImageShowDialogClass.h new file mode 100644 index 0000000..4c08e08 --- /dev/null +++ b/ImageShowDialogClass.h @@ -0,0 +1,154 @@ + +#pragma once + +#ifndef IMAGESHOWDIALOGCLASS_H +#define IMAGESHOWDIALOGCLASS_H +#include "BaseConstVariable.h" +#include +#include +#include "qcustomplot.h" +#include +#include + + +//=========================== +// 定义需要绘制的图像的类型 +//=========================== +enum LAMPDATASHOWCLASS { + LAMPGraph, + LAMPColorMap, + LAMPScatterStyle, + NOWINDOWS +}; + +namespace Ui { + class ImageShowDialogClass; +}; + + + + +//=========================== +// 构建游标类型 +// 1. 单纯的描述游标，主要用来展示坐标，还有当前数据信息 +// 2. 区域性描述游标，通过线，等用来展示某一个区域的信息 +//=========================== +class ImageShowCursorDesClass : public QDialog +{ + Q_OBJECT +public: + QLabel* label; + + +public: + ImageShowCursorDesClass(QWidget* parent = nullptr); + ~ImageShowCursorDesClass(); +public slots: + void updateCursorContent(QString content); +signals: + bool windowsClose(); +protected: + void closeEvent(QCloseEvent* event) override; + +}; + + +class ImageShowCursorLineClass :public QDialog { + Q_OBJECT +public: + QMenuBar* menubar; + QMenu* DataList; + QMenu* RemoveList; + + QCustomPlot* customPlot; + bool tracerEnable; + + QCPItemTracer* tracer; //游标 + QCPItemText* tracerYText; //图标标签 + QCPItemText* tracerXText; //游标标签 + QMap m_graphMap; + QMap data_action_map; + QString TrackName; + +public: + ImageShowCursorLineClass(QWidget* parent = nullptr); + ~ImageShowCursorLineClass(); + +public: + void loadPlotLine(QVector xs, QVector dataValues, QString Name); + void UpdatePlotLine(QVector xs, QVector dataValues, QString Name); + void updateCursor(QMouseEvent* event); + + +public slots: + void xAxisRangeChanged(QCPRange); + void yAxisRangeChanged(QCPRange); + void on_SwichTracerGraph_Name(); + void load_new_compare_line(); + +signals: + bool windowsClose(); +protected: + void closeEvent(QCloseEvent* event) override; +}; + + + + +//=========================== +// 构建数据展示窗口 +// 1. 单纯的描述游标，主要用来展示坐标，还有当前数据信息 +// 2. 区域性描述游标，通过线，等用来展示某一个区域的信息 +//=========================== +class ImageShowDialogClass : public QDialog +{ + Q_OBJECT +private: + LAMPDATASHOWCLASS graphclass; +public: + QMenuBar* menubar; + QStatusBar* statusbar; + ImageShowCursorDesClass* desCursor; // 描述游标 + ImageShowCursorLineClass* HlineCursor;// H 游标 + ImageShowCursorLineClass* VlineCursor;// V 游标 + + bool desCursorflag; + bool HlineCursorflag; + bool VlineCursorflag; + + QCustomPlot* m_plot; + QCPItemTracer* tracer; //游标 + +public: + QMap m_graphMap; +public: + ImageShowDialogClass(QWidget *parent = nullptr); + ~ImageShowDialogClass(); +public: + + void PlotLine(QVector xs, QVector ys, QString Name); + + void Scatter(double dx, double dy); + void Scatter(QVector xs, QVector ys); + void load_double_MatrixX_data(Eigen::MatrixXd data, QString name); + void load_double_MatrixX_data(Eigen::MatrixXd X,Eigen::MatrixXd Y,Eigen::MatrixXd data, QString name); + void remove_Data( QString name); + + LAMPDATASHOWCLASS getGraphClass(size_t i = 0); + size_t getGraphCount(); +private: + Ui::ImageShowDialogClass* ui; + +public slots: + void updateCursor(QMouseEvent* event); + + +public slots: // cursor + void on_action_cursor_enable_trigged(); + void on_action_descursor_close_trigged(); + void on_action_Hlinecursor_close_trigged(); + void on_action_VVlinecursor_close_trigged(); +}; + + +#endif diff --git a/ImageShowDialogClass.ui b/ImageShowDialogClass.ui new file mode 100644 index 0000000..271d5bc --- /dev/null +++ b/ImageShowDialogClass.ui @@ -0,0 +1,33 @@ + + + ImageShowDialogClass + + + + 0 + 0 + 600 + 400 + + + + 复数展示 + + + + + + + + + + + QCustomPlot + QWidget +
qcustomplot.h
+ 1 + + + + + diff --git a/LogInfoCls.cpp b/LogInfoCls.cpp new file mode 100644 index 0000000..47cba8f --- /dev/null +++ b/LogInfoCls.cpp @@ -0,0 +1,228 @@ +#include "stdafx.h" +#include "LogInfoCls.h" + +std::string errorCode2errInfo(ErrorCode e) +{ + switch (e) + { + _CASE_STR(SUCCESS); + _CASE_STR(VIRTUALABSTRACT); + _CASE_STR(FILENOFOUND); + _CASE_STR(OrbitNodeNotEnough); + _CASE_STR(XYDataPointNotEqual); + _CASE_STR(FILEOPENFAIL ); + _CASE_STR(XMLPARSEFAIL ); + _CASE_STR(XMLNOTFOUNDElEMENT ); + _CASE_STR(FILEPATHISEMPTY); + _CASE_STR(FOLDER_NOT_EXIST); + _CASE_STR(FILE_NOT_EXIST); + _CASE_STR(FIND_ID_ERROR); + _CASE_STR(INSERT_ID_ERROR); + //GSL 1xx + _CASE_STR(Error_GSL_FAILURE ); + _CASE_STR(Error_GSL_CONTINUE ); + _CASE_STR(Error_GSL_EDOM ); + _CASE_STR(Error_GSL_ERANGE ); + _CASE_STR(Error_GSL_EFAULT ); + _CASE_STR(Error_GSL_EINVAL ); + _CASE_STR(Error_GSL_EFAILED ); + _CASE_STR(Error_GSL_EFACTOR ); + _CASE_STR(Error_GSL_ESANITY ); + _CASE_STR(Error_GSL_ENOMEM ); + _CASE_STR(Error_GSL_EBADFUNC ); + _CASE_STR(Error_GSL_ERUNAWAY ); + _CASE_STR(Error_GSL_EMAXITER ); + _CASE_STR(Error_GSL_EZERODIV ); + _CASE_STR(Error_GSL_EBADTOL ); + _CASE_STR(Error_GSL_ETOL ); + _CASE_STR(Error_GSL_EUNDRFLW ); + _CASE_STR(Error_GSL_EOVRFLW ); + _CASE_STR(Error_GSL_ELOSS ); + _CASE_STR(Error_GSL_EROUND ); + _CASE_STR(Error_GSL_EBADLEN ); + _CASE_STR(Error_GSL_ENOTSQR ); + _CASE_STR(Error_GSL_ESING ); + _CASE_STR(Error_GSL_EDIVERGE ); + _CASE_STR(Error_GSL_EUNSUP ); + _CASE_STR(Error_GSL_EUNIMPL ); + _CASE_STR(Error_GSL_ECACHE ); + _CASE_STR(Error_GSL_ETABLE ); + _CASE_STR(Error_GSL_ENOPROG ); + _CASE_STR(Error_GSL_ENOPROGJ ); + _CASE_STR(Error_GSL_ETOLF ); + _CASE_STR(Error_GSL_ETOLX ); + _CASE_STR(Error_GSL_ETOLG ); + _CASE_STR(Error_GSL_EOF ); + // RTPC + _CASE_STR(RTPC_PARAMSISEMPTY); + _CASE_STR(ECHO_L0DATA_NOTOPEN); + _CASE_STR(ECHO_L0DATA_ROW_COL_NOEQUAL); + _CASE_STR(ECHO_L0DATA_PRFIDXOUTRANGE); + _CASE_STR(ECHO_L0DATA_GPSFILEFORMATERROR); + _CASE_STR(ECHO_L0DATA_ECHOFILEFORMATERROR); + _CASE_STR(ECHO_L0DATA_ECHOFILENOTOPEN); + _CASE_STR(ECHO_L0DATA_GPSFILEFNOTOPEN); + _CASE_STR(ECHO_L0DATA_XMLFILENOTOPEN); + _CASE_STR(OUTOFRANGE); + _CASE_STR(ECHO_L0DATA_XMLNAMEERROR); + // + _CASE_STR(TBP_L0OPENFAIL); + + // + _CASE_STR(IMAGE_L1DATA_XMLNAMEERROR); + _CASE_STR(IMAGE_L1DATA_XMLNAMEOPENERROR); + + + default: + break; + } + return "UNKNOW_EVENT!"; +} + +ErrorCode GSLState2ErrorCode(int gslState) +{ + + switch (gslState) + { + case 0: + return ErrorCode::SUCCESS; + break; + case -1: + return ErrorCode::Error_GSL_FAILURE; + break; + case -2: + return ErrorCode::Error_GSL_CONTINUE; + break; + case 1: + return ErrorCode::Error_GSL_EDOM; // sqrt(-1) + break; + + case 2: + return ErrorCode::Error_GSL_ERANGE; // Χ exp(1e100) + break; + + case 3: + return ErrorCode::Error_GSL_EFAULT; // Чָ + break; + + case 4: + return ErrorCode::Error_GSL_EINVAL; // ûṩЧ + break; + + case 5: + return ErrorCode::Error_GSL_EFAILED; // ͨʧ + break; + + case 6: + return ErrorCode::Error_GSL_EFACTOR; // ӷֽʧ + break; + + case 7: + return ErrorCode::Error_GSL_ESANITY; // ǼʧܣͨӦ÷ + break; + + case 8: + return ErrorCode::Error_GSL_ENOMEM; // ڴʧ + break; + + case 9: + return ErrorCode::Error_GSL_EBADFUNC; // ûṩĺ + break; + + case 10: + return ErrorCode::Error_GSL_ERUNAWAY; // ʧ + break; + + case 11: + return ErrorCode::Error_GSL_EMAXITER; // + break; + + case 12: + return ErrorCode::Error_GSL_EZERODIV; // Խ + break; + + case 13: + return ErrorCode::Error_GSL_EBADTOL; // û̶ָЧ + break; + + case 14: + return ErrorCode::Error_GSL_ETOL; // δܴﵽ̶ָ + break; + + case 15: + return ErrorCode::Error_GSL_EUNDRFLW; // + break; + + case 16: + return ErrorCode::Error_GSL_EOVRFLW; // + break; + + case 17: + return ErrorCode::Error_GSL_ELOSS; // ʧ + break; + + case 18: + return ErrorCode::Error_GSL_EROUND; // ʧ + break; + + case 19: + return ErrorCode::Error_GSL_EBADLEN; // Ȳƥ + break; + + case 20: + return ErrorCode::Error_GSL_ENOTSQR; // Ƿ + break; + + case 21: + return ErrorCode::Error_GSL_ESING; // ⵽Ե + break; + + case 22: + return ErrorCode::Error_GSL_EDIVERGE; // ֻɢ + break; + + case 23: + return ErrorCode::Error_GSL_EUNSUP; // ԲӲ֧ + break; + + case 24: + return ErrorCode::Error_GSL_EUNIMPL; // δʵ + break; + + case 25: + return ErrorCode::Error_GSL_ECACHE; // + break; + + case 26: + return ErrorCode::Error_GSL_ETABLE; // + break; + + case 27: + return ErrorCode::Error_GSL_ENOPROG; // δܳȡýչ + break; + + case 28: + return ErrorCode::Error_GSL_ENOPROGJ; // ſɱδܸƽ + break; + + case 29: + return ErrorCode::Error_GSL_ETOLF; // ޷ﵽ F ̶ָ + break; + + case 30: + return ErrorCode::Error_GSL_ETOLX; // ޷ﵽ X ̶ָ + break; + + case 31: + return ErrorCode::Error_GSL_ETOLG; // ޷ﵽݶȵ̶ָ + break; + + case 32: + return ErrorCode::Error_GSL_EOF; // ļ + break; + + default: + return ErrorCode::Error_GSL_FAILURE; // δ֪󣬷һʧ + break; + } +} diff --git a/LogInfoCls.h b/LogInfoCls.h new file mode 100644 index 0000000..c439b9c --- /dev/null +++ b/LogInfoCls.h @@ -0,0 +1,97 @@ +#pragma once +/*****************************************************************//** + * \file LogInfoCls.h + * \brief ö࣬ԼشϢ + * + * \author + * \date October 2024 + *********************************************************************/ + +#include + +// 任 +#define _CASE_STR(x) case x : return #x; break; + + +enum ErrorCode { + SUCCESS = -1,// ִгɹ + VIRTUALABSTRACT = -2,// virtual abstract function not implement + FILENOFOUND = 1, + OrbitNodeNotEnough = 2, + XYDataPointNotEqual = 3, + FILEOPENFAIL = 4, + XMLPARSEFAIL = 5, + XMLNOTFOUNDElEMENT = 6, + FILEPATHISEMPTY = 7, + FOLDER_NOT_EXIST = 8, + FILE_NOT_EXIST = 9, + FIND_ID_ERROR = 10, + INSERT_ID_ERROR = 11, + //GSL 1xx + Error_GSL_FAILURE = -101, + Error_GSL_CONTINUE = -102, /* iteration has not converged */ + Error_GSL_EDOM = 101, /* input domain error, e.g sqrt(-1) */ + Error_GSL_ERANGE = 102, /* output range error, e.g. exp(1e100) */ + Error_GSL_EFAULT = 103, /* invalid pointer */ + Error_GSL_EINVAL = 104, /* invalid argument supplied by user */ + Error_GSL_EFAILED = 105, /* generic failure */ + Error_GSL_EFACTOR = 106, /* factorization failed */ + Error_GSL_ESANITY = 107, /* sanity check failed - shouldn't happen */ + Error_GSL_ENOMEM = 108, /* malloc failed */ + Error_GSL_EBADFUNC = 109, /* problem with user-supplied function */ + Error_GSL_ERUNAWAY = 110, /* iterative process is out of control */ + Error_GSL_EMAXITER = 111, /* exceeded max number of iterations */ + Error_GSL_EZERODIV = 112, /* tried to divide by zero */ + Error_GSL_EBADTOL = 113, /* user specified an invalid tolerance */ + Error_GSL_ETOL = 114, /* failed to reach the specified tolerance */ + Error_GSL_EUNDRFLW = 115, /* underflow */ + Error_GSL_EOVRFLW = 116, /* overflow */ + Error_GSL_ELOSS = 117, /* loss of accuracy */ + Error_GSL_EROUND = 118, /* failed because of roundoff error */ + Error_GSL_EBADLEN = 119, /* matrix, vector lengths are not conformant */ + Error_GSL_ENOTSQR = 120, /* matrix not square */ + Error_GSL_ESING = 121, /* apparent singularity detected */ + Error_GSL_EDIVERGE = 122, /* integral or series is divergent */ + Error_GSL_EUNSUP = 123, /* requested feature is not supported by the hardware */ + Error_GSL_EUNIMPL = 124, /* requested feature not (yet) implemented */ + Error_GSL_ECACHE = 125, /* cache limit exceeded */ + Error_GSL_ETABLE = 126, /* table limit exceeded */ + Error_GSL_ENOPROG = 127, /* iteration is not making progress towards solution */ + Error_GSL_ENOPROGJ = 128, /* jacobian evaluations are not improving the solution */ + Error_GSL_ETOLF = 129, /* cannot reach the specified tolerance in F */ + Error_GSL_ETOLX = 130, /* cannot reach the specified tolerance in X */ + Error_GSL_ETOLG = 131, /* cannot reach the specified tolerance in gradient */ + Error_GSL_EOF = 132, /* end of file */ + + // RTPC + RTPC_PARAMSISEMPTY = 201, + // L0 + ECHO_L0DATA_NOTOPEN = 202, + ECHO_L0DATA_ROW_COL_NOEQUAL = 203, // л + ECHO_L0DATA_PRFIDXOUTRANGE = 204, // PRF Χ + ECHO_L0DATA_GPSFILEFORMATERROR = 205, // GPSļ + ECHO_L0DATA_ECHOFILEFORMATERROR = 206, // زļʽ + ECHO_L0DATA_ECHOFILENOTOPEN = 207, // زļ򿪴 + ECHO_L0DATA_GPSFILEFNOTOPEN = 208, // GPSļ򿪴 + ECHO_L0DATA_XMLFILENOTOPEN = 209, // xmlļ򿪴 + OUTOFRANGE = 210, // Χ + ECHO_L0DATA_XMLNAMEERROR = 211, // Χ + // BP + TBP_L0OPENFAIL = 301, // 0ļ򿪴 + + // L1ͼ + IMAGE_L1DATA_XMLNAMEERROR = 401, + IMAGE_L1DATA_XMLNAMEOPENERROR = 402, + IMAGE_L1DATA_XMLNAMEPARASEERROR = 403, + +}; + + + +std::string errorCode2errInfo(ErrorCode code); + + +ErrorCode GSLState2ErrorCode(int gslState); + + + diff --git a/SARSimulationImageL1.cpp b/SARSimulationImageL1.cpp new file mode 100644 index 0000000..13ffcc3 --- /dev/null +++ b/SARSimulationImageL1.cpp @@ -0,0 +1,515 @@ +#include "stdafx.h" +#include "SARSimulationImageL1.h" +#include "LogInfoCls.h" +#include +#include +#include +#include + +SARSimulationImageL1Dataset::SARSimulationImageL1Dataset() +{ +} + +SARSimulationImageL1Dataset::~SARSimulationImageL1Dataset() +{ +} + +QString SARSimulationImageL1Dataset::getoutFolderPath() +{ + return outFolderPath; +} + +QString SARSimulationImageL1Dataset::getDatesetName() +{ + return L1DatasetName; +} + +QString SARSimulationImageL1Dataset::getxmlfileName() +{ + return xmlfileName; +} + +QString SARSimulationImageL1Dataset::getxmlFilePath() +{ + return xmlFilePath; +} + +QString SARSimulationImageL1Dataset::getImageRasterName() +{ + return ImageRasterName; +} + +QString SARSimulationImageL1Dataset::getImageRasterPath() +{ + return ImageRasterPath; +} + +QString SARSimulationImageL1Dataset::getGPSPointFilename() +{ + return GPSPointFilename; +} + +QString SARSimulationImageL1Dataset::getGPSPointFilePath() +{ + return GPSPointFilePath; +} + +ErrorCode SARSimulationImageL1Dataset::OpenOrNew(QString folder, QString filename, long inrowCount, long incolCount) +{ + qDebug() << "--------------Image Raster L1 Data OpenOrNew ---------------------------------------"; + qDebug() << "Folder: " << folder; + qDebug() << "Filename: " << filename; + QDir dir(folder); + if (dir.exists() == false) + { + dir.mkpath(folder); + } + else {} + + if (dir.exists() == false) { + return ErrorCode::FOLDER_NOT_EXIST; + } + else {} + QString filePath = dir.filePath(filename); // ļ· + this->rowCount = inrowCount; + this->colCount = incolCount; + + this->outFolderPath = folder; + this->L1DatasetName = filename; + + this->xmlfileName = filename + ".xml"; + this->GPSPointFilename = filename + ".gpspos.data"; + this->ImageRasterName = filename + ".bin"; + + this->xmlFilePath = dir.filePath(this->xmlfileName); + this->GPSPointFilePath = dir.filePath(this->GPSPointFilename); + this->ImageRasterPath = dir.filePath(this->ImageRasterName); + + if (QFile(this->xmlFilePath).exists()) + { + this->loadFromXml(); + } + else + { + this->saveToXml(); + } + + if (QFile(this->GPSPointFilePath).exists() == false) { + // ļ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + std::shared_ptr poDstDS(poDriver->Create(this->GPSPointFilePath.toUtf8().constData(), 16, rowCount, 1, GDT_Float64, NULL)); + GDALFlushCache((GDALDatasetH)poDstDS.get()); + poDstDS.reset(); + omp_unset_lock(&lock); + omp_destroy_lock(&lock); + } + + if (QFile(this->ImageRasterPath).exists() == false) { + + // ļ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + std::shared_ptr poDstDS(poDriver->Create(this->ImageRasterPath.toUtf8().constData(), colCount, rowCount, 1, GDT_CFloat64, NULL)); + GDALFlushCache((GDALDatasetH)poDstDS.get()); + poDstDS.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + } + + return ErrorCode::SUCCESS; +} + +ErrorCode SARSimulationImageL1Dataset::Open(QString folder, QString filename) +{ + + QDir dir(folder); + if (dir.exists() == false) + { + return ErrorCode::FOLDER_NOT_EXIST; + } + else {} + + if (dir.exists() == false) { + return ErrorCode::FOLDER_NOT_EXIST; + } + else {} + QString filePath = dir.filePath(filename); // ļ· + + this->outFolderPath = folder; + this->L1DatasetName = filename; + + this->xmlfileName = filename + ".xml"; + this->GPSPointFilename = filename + ".gpspos.data"; + this->ImageRasterName = filename + ".bin"; + + this->xmlFilePath = dir.filePath(this->xmlfileName); + this->GPSPointFilePath = dir.filePath(this->GPSPointFilename); + this->ImageRasterPath = dir.filePath(this->ImageRasterName); + + if (QFile(this->xmlFilePath).exists()) + { + this->loadFromXml(); + } + else + { + return ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR; + } +} + +ErrorCode SARSimulationImageL1Dataset::Open(QString xmlPath) +{ + QFileInfo fileInfo(xmlPath); + QString fileName = fileInfo.fileName(); // ȡļ + QString fileSuffix = fileInfo.suffix(); // ȡ׺ + QString fileNameWithoutSuffix = fileInfo.baseName(); // ȡ׺ļ + QString directoryPath = fileInfo.path(); // ȡļ· + if (fileSuffix.toLower() == "xml" || fileSuffix.toLower() == ".xml") { + return this->Open(directoryPath, fileNameWithoutSuffix); + } + else { + return ErrorCode::IMAGE_L1DATA_XMLNAMEERROR; + } + return ErrorCode::SUCCESS; + +} + +void SARSimulationImageL1Dataset::saveToXml() +{ + QFile file(this->xmlFilePath); + if (!file.open(QIODevice::WriteOnly | QIODevice::Text)) { + qDebug() << "Cannot open file for writing:" << file.errorString(); + return; + } + + QXmlStreamWriter xmlWriter(&file); + xmlWriter.setAutoFormatting(true); + + xmlWriter.writeStartDocument(); + xmlWriter.writeStartElement("Parameters"); + + xmlWriter.writeTextElement("RowCount", QString::number(this->rowCount)); + xmlWriter.writeTextElement("ColCount", QString::number(this->colCount)); + xmlWriter.writeTextElement("Rnear", QString::number(this->Rnear)); + xmlWriter.writeTextElement("Rfar", QString::number(this->Rfar)); + xmlWriter.writeTextElement("Rref", QString::number(this->Rref)); + xmlWriter.writeTextElement("CenterFreq", QString::number(this->centerFreq)); + xmlWriter.writeTextElement("Fs", QString::number(this->Fs)); + xmlWriter.writeTextElement("CenterAngle", QString::number(this->CenterAngle)); + xmlWriter.writeTextElement("LookSide", this->LookSide); + + xmlWriter.writeEndElement(); // Parameters + xmlWriter.writeEndDocument(); + + file.close(); +} + +ErrorCode SARSimulationImageL1Dataset::loadFromXml() +{ + QFile file(this->xmlFilePath); + if (!file.open(QIODevice::ReadOnly | QIODevice::Text)) { + qDebug() << "Cannot open file for reading:" << file.errorString(); + return ErrorCode::IMAGE_L1DATA_XMLNAMEOPENERROR; + } + + QXmlStreamReader xmlReader(&file); + while (!xmlReader.atEnd() && !xmlReader.hasError()) { + QXmlStreamReader::TokenType token = xmlReader.readNext(); + if (token == QXmlStreamReader::StartDocument) { + continue; + } + if (token == QXmlStreamReader::StartElement) { + if (xmlReader.name() == "Parameters") { + continue; + } + else if (xmlReader.name() == "RowCount") { + this->rowCount = xmlReader.readElementText().toLong(); + } + else if (xmlReader.name() == "ColCount") { + this->colCount = xmlReader.readElementText().toLong(); + } + else if (xmlReader.name() == "Rnear") { + this->Rnear = xmlReader.readElementText().toDouble(); + } + else if (xmlReader.name() == "Rfar") { + this->Rfar = xmlReader.readElementText().toDouble(); + } + else if (xmlReader.name() == "Rref") { + this->Rref = xmlReader.readElementText().toDouble(); + } + else if (xmlReader.name() == "CenterFreq") { + this->centerFreq = xmlReader.readElementText().toDouble(); + } + else if (xmlReader.name() == "Fs") { + this->Fs = xmlReader.readElementText().toDouble(); + } + else if (xmlReader.name() == "CenterAngle") { + this->CenterAngle = xmlReader.readElementText().toDouble(); + } + else if (xmlReader.name() == "LookSide") { + this->LookSide = xmlReader.readElementText(); + } + } + } + + if (xmlReader.hasError()) { + qDebug() << "XML error:" << xmlReader.errorString(); + return ErrorCode::IMAGE_L1DATA_XMLNAMEPARASEERROR; + } + + file.close(); + return ErrorCode::SUCCESS; +} + +std::shared_ptr SARSimulationImageL1Dataset::getAntPos() +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + // ȡļ + std::shared_ptr rasterDataset = OpenDataset(this->GPSPointFilePath, GDALAccess::GA_ReadOnly); + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* demBand = rasterDataset->GetRasterBand(1); + + long width = rasterDataset->GetRasterXSize(); + long height = rasterDataset->GetRasterYSize(); + long band_num = rasterDataset->GetRasterCount(); + + std::shared_ptr temp = nullptr; + + if (gdal_datatype == GDT_Float64) { + temp = std::shared_ptr(new double[this->rowCount * 16], delArrPtr); + demBand->RasterIO(GF_Read, 0, 0, 10, this->rowCount, temp.get(), 16, this->rowCount, gdal_datatype, 0, 0); + } + else { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_GPSFILEFORMATERROR)); + } + rasterDataset.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return temp; +} + +ErrorCode SARSimulationImageL1Dataset::saveAntPos(std::shared_ptr ptr) +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + // ȡļ + std::shared_ptr rasterDataset = OpenDataset(this->GPSPointFilePath, GDALAccess::GA_Update); + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* demBand = rasterDataset->GetRasterBand(1); + + long width = rasterDataset->GetRasterXSize(); + long height = rasterDataset->GetRasterYSize(); + long band_num = rasterDataset->GetRasterCount(); + + if (gdal_datatype == GDT_Float64) { + demBand->RasterIO(GF_Write, 0, 0, 16, this->rowCount, ptr.get(), 16, this->rowCount, gdal_datatype, 0, 0); + } + else { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_GPSFILEFORMATERROR)); + omp_unset_lock(&lock); + omp_destroy_lock(&lock); + return ErrorCode::ECHO_L0DATA_GPSFILEFORMATERROR; + } + rasterDataset.reset(); + omp_unset_lock(&lock); + omp_destroy_lock(&lock); + return ErrorCode::SUCCESS; +} + +std::shared_ptr> SARSimulationImageL1Dataset::getImageRaster() +{ + return this->getImageRaster(0, this->rowCount); +} + +ErrorCode SARSimulationImageL1Dataset::saveImageRaster(std::shared_ptr> echoPtr, long startPRF, long PRFLen) +{ + if (!(startPRF < this->rowCount)) { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_PRFIDXOUTRANGE)); + return ErrorCode::ECHO_L0DATA_PRFIDXOUTRANGE; + } + else {} + + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + std::shared_ptr rasterDataset = OpenDataset(this->ImageRasterPath, GDALAccess::GA_Update); + + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* poBand = rasterDataset->GetRasterBand(1); + if (NULL == poBand || nullptr == poBand) { + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + return ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR; + } + else {} + long width = rasterDataset->GetRasterXSize(); + long height = rasterDataset->GetRasterYSize(); + long band_num = rasterDataset->GetRasterCount(); + + if (height != this->rowCount || width != this->colCount) { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR; + } + else { + if (gdal_datatype == GDT_CFloat64) { + poBand->RasterIO(GF_Write, 0, startPRF, width, PRFLen, echoPtr.get(), width, PRFLen, GDT_CFloat64, 0, 0); + GDALFlushCache((GDALDatasetH)rasterDataset.get()); + } + else { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + } + } + rasterDataset.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return ErrorCode::SUCCESS; +} + +std::shared_ptr> SARSimulationImageL1Dataset::getImageRaster(long startPRF, long PRFLen) +{ + if (!(startPRF < this->rowCount)) { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_PRFIDXOUTRANGE)) << startPRF << " " << this->rowCount; + return nullptr; + } + else {} + + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + std::shared_ptr rasterDataset = OpenDataset(this->ImageRasterPath, GDALAccess::GA_ReadOnly); + + + GDALDataType gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType(); + GDALRasterBand* poBand = rasterDataset->GetRasterBand(1); + if (NULL == poBand || nullptr == poBand) { + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + return nullptr; + } + else {} + long width = rasterDataset->GetRasterXSize(); + long height = rasterDataset->GetRasterYSize(); + long band_num = rasterDataset->GetRasterCount(); + std::shared_ptr> temp = nullptr; + if (height != this->rowCount || width != this->colCount) { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + } + else { + if (gdal_datatype == GDT_CFloat64) { + temp = std::shared_ptr>(new std::complex[PRFLen * width], delArrPtr); + poBand->RasterIO(GF_Read, 0, startPRF, width, PRFLen, temp.get(), width, PRFLen, GDT_CFloat64, 0, 0); + GDALFlushCache((GDALDatasetH)rasterDataset.get()); + } + else { + qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR)); + } + } + + rasterDataset.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + return temp; +} + +Eigen::MatrixXcd SARSimulationImageL1Dataset::getImageRasterMatrix() +{ + std::shared_ptr> data = this->getImageRaster(); + Eigen::MatrixXcd dataimg = Eigen::MatrixXcd::Zero(this->rowCount, this->colCount); + for (long i = 0; i < this->rowCount; i++) { + for (long j = 0; j < this->colCount; j++) { + dataimg(i, j) = data.get()[i*colCount+j]; + } + } + return Eigen::MatrixXcd(); +} + +ErrorCode SARSimulationImageL1Dataset::saveImageRaster(Eigen::MatrixXcd& dataimg, long startPRF) +{ + std::shared_ptr> data(new std::complex[dataimg.rows()* dataimg.cols()]); + long dataimgrows = dataimg.rows(); + long dataimgcols = dataimg.cols(); + for (long i = 0; i < dataimgrows; i++) { + for (long j = 0; j < dataimgcols; j++) { + data.get()[i * dataimgcols + j] = dataimg(i, j); + } + } + this->saveImageRaster(data, startPRF,dataimgrows); + return ErrorCode(); +} + + + + +// Getter Setter ʵ +long SARSimulationImageL1Dataset::getrowCount() { return this->rowCount; } +void SARSimulationImageL1Dataset::setrowCount(long rowCountin) { this->rowCount = rowCountin; } + +long SARSimulationImageL1Dataset::getcolCount() { return this->colCount; } +void SARSimulationImageL1Dataset::setcolCount(long pulsePointsin) { this->colCount = pulsePointsin; } + + +double SARSimulationImageL1Dataset::getNearRange() { return this->Rnear; } +void SARSimulationImageL1Dataset::setNearRange(double nearRange) { this->Rnear = nearRange; } + +double SARSimulationImageL1Dataset::getFarRange() { return this->Rfar; } +void SARSimulationImageL1Dataset::setFarRange(double farRange) { this->Rfar = farRange; } + +double SARSimulationImageL1Dataset::getRefRange() +{ + return this->Rref; +} + +void SARSimulationImageL1Dataset::setRefRange(double refRange) +{ + this->Rref = refRange; +} + +double SARSimulationImageL1Dataset::getCenterFreq() { return this->centerFreq; } +void SARSimulationImageL1Dataset::setCenterFreq(double freq) { this->centerFreq = freq; } + +double SARSimulationImageL1Dataset::getFs() { return this->Fs; } +void SARSimulationImageL1Dataset::setFs(double samplingFreq) { this->Fs = samplingFreq; } + + + +double SARSimulationImageL1Dataset::getCenterAngle() +{ + return this->CenterAngle; +} + +void SARSimulationImageL1Dataset::setCenterAngle(double angle) +{ + this->CenterAngle = angle; +} + +QString SARSimulationImageL1Dataset::getLookSide() +{ + return this->LookSide; +} + +void SARSimulationImageL1Dataset::setLookSide(QString looksides) +{ + this->LookSide = looksides; +} diff --git a/SARSimulationImageL1.h b/SARSimulationImageL1.h new file mode 100644 index 0000000..e95d0b9 --- /dev/null +++ b/SARSimulationImageL1.h @@ -0,0 +1,105 @@ +#pragma once + +#include "BaseConstVariable.h" + +#include "BaseTool.h" +#include "ImageOperatorBase.h" +#include "GeoOperator.h" +#include "FileOperator.h" +#include "LogInfoCls.h" + + + +class SARSimulationImageL1Dataset +{ +public: + SARSimulationImageL1Dataset(); + ~SARSimulationImageL1Dataset(); + +private : + QString outFolderPath; + QString L1DatasetName; + QString xmlfileName; + QString xmlFilePath; + QString ImageRasterName; + QString ImageRasterPath; + QString GPSPointFilename; + QString GPSPointFilePath; + +public: + QString getoutFolderPath(); + QString getDatesetName(); + QString getxmlfileName(); + QString getxmlFilePath(); + QString getImageRasterName(); + QString getImageRasterPath(); + QString getGPSPointFilename(); // GPS + QString getGPSPointFilePath(); + + +public: + ErrorCode OpenOrNew(QString folder, QString filename, long rowCount, long colCount); + ErrorCode Open(QString folderPath, QString Name); + ErrorCode Open(QString xmlPath); +public: + void saveToXml(); + ErrorCode loadFromXml(); + std::shared_ptr getAntPos(); + ErrorCode saveAntPos(std::shared_ptr ptr); // עΣգдǰǷȷ + + std::shared_ptr> getImageRaster(); + ErrorCode saveImageRaster(std::shared_ptr> echoPtr, long startPRF, long PRFLen); + std::shared_ptr> getImageRaster(long startPRF, long PRFLen); + + Eigen::MatrixXcd getImageRasterMatrix(); + ErrorCode saveImageRaster(Eigen::MatrixXcd& data, long startPRF); + +private: // xmlв + + long rowCount; + long colCount; + + double Rnear; + double Rfar; + double Rref; + + double centerFreq; + double Fs; + + double CenterAngle; + QString LookSide; + +public: + long getrowCount(); + void setrowCount(long rowCount); + + long getcolCount(); + void setcolCount(long pulsePoints); + + double getNearRange(); + void setNearRange(double nearRange); + + double getFarRange(); + void setFarRange(double farRange); + + double getRefRange(); + void setRefRange(double refRange); + + double getCenterFreq(); + void setCenterFreq(double freq); + + double getFs(); + void setFs(double samplingFreq); + + double getCenterAngle(); + void setCenterAngle(double angle); + + QString getLookSide(); + void setLookSide(QString lookside); + + + + + +}; + diff --git a/stdafx.cpp b/stdafx.cpp new file mode 100644 index 0000000..fd4f341 --- /dev/null +++ b/stdafx.cpp @@ -0,0 +1 @@ +#include "stdafx.h" diff --git a/stdafx.h b/stdafx.h new file mode 100644 index 0000000..e69de29