From 26fd7c20aec67c8dfa106df45593559aac24a4c9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=E9=99=88=E5=A2=9E=E8=BE=89?= <3045316072@qq.com> Date: Fri, 15 Nov 2024 17:25:05 +0800 Subject: [PATCH] =?UTF-8?q?=E6=9B=B4=E6=96=B0=E5=90=88=E5=B9=B6=E4=BB=A3?= =?UTF-8?q?=E7=A0=81?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- ImageOperatorBase.cpp | 218 +++++++++++++++++++++++++++--------------- ImageOperatorBase.h | 25 +++-- 2 files changed, 154 insertions(+), 89 deletions(-) diff --git a/ImageOperatorBase.cpp b/ImageOperatorBase.cpp index 5ee933c..c8b2ed2 100644 --- a/ImageOperatorBase.cpp +++ b/ImageOperatorBase.cpp @@ -21,8 +21,10 @@ #include #include #include - - +#include +#include +#include +#include /** * 输入数据是ENVI格式数据 @@ -569,28 +571,28 @@ Eigen::MatrixXd gdalImage::getData(int start_row, int start_col, int rows_count, } 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_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, @@ -695,28 +697,27 @@ Eigen::MatrixXi gdalImage::getDatai(int start_row, int start_col, int rows_count } 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_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, @@ -773,6 +774,7 @@ GDALDataType gdalImage::getDataType() void gdalImage::saveImage(Eigen::MatrixXd data, int start_row = 0, int start_col = 0, int band_ids = 1) { + GDALDataType datetype = this->getDataType(); omp_lock_t lock; omp_init_lock(&lock); omp_set_lock(&lock); @@ -789,7 +791,15 @@ void gdalImage::saveImage(Eigen::MatrixXd data, int start_row = 0, int start_col 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); // 斤拷锟斤拷 + this->band_num, datetype, NULL); // 斤拷锟斤拷 + + if (nullptr == poDstDS) { + 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()); + return; + } + poDstDS->SetProjection(this->projection.toUtf8().constData()); double gt_ptr[6]; @@ -802,6 +812,9 @@ void gdalImage::saveImage(Eigen::MatrixXd data, int start_row = 0, int start_col //delete gt_ptr; } + + + int datarows = data.rows(); int datacols = data.cols(); @@ -817,7 +830,7 @@ void gdalImage::saveImage(Eigen::MatrixXd data, int start_row = 0, int start_col // 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); + databuffer, datacols, datarows, datetype, 0, 0); GDALFlushCache(poDstDS); GDALClose((GDALDatasetH)poDstDS); @@ -829,6 +842,8 @@ void gdalImage::saveImage(Eigen::MatrixXd data, int start_row = 0, int start_col void gdalImage::saveImage(Eigen::MatrixXi data, int start_row, int start_col, int band_ids) { + + GDALDataType datetype=this->getDataType(); omp_lock_t lock; omp_init_lock(&lock); omp_set_lock(&lock); @@ -872,7 +887,7 @@ void gdalImage::saveImage(Eigen::MatrixXi data, int start_row, int start_col, in // 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_Int64, 0, 0); + databuffer, datacols, datarows, datetype, 0, 0); GDALFlushCache(poDstDS); GDALClose((GDALDatasetH)poDstDS); @@ -1544,8 +1559,10 @@ int saveMatrixXcd2TiFF(Eigen::MatrixXcd data, QString out_tiff_path) return -1; } -ErrorCode MergeRasterProcess(QVector filepaths, QString outfileptah, QString mainString, MERGEMODE mergecode, bool isENVI) +ErrorCode MergeRasterProcess(QVector filepaths, QString outfileptah, QString mainString, MERGEMODE mergecode, bool isENVI, ShowProessAbstract* dia ) { + + // 参数检查 if (!isExists(mainString)) { qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::FILENOFOUND) )<< "\t" << mainString; @@ -1588,8 +1605,6 @@ ErrorCode MergeRasterProcess(QVector filepaths, QString outfileptah, QS } } - - Eigen::MatrixXd maingt = mainimg.getGeoTranslation(); Eigen::MatrixXd rgt = Eigen::MatrixXd::Zero(2,3); RasterExtend mainExtend = mainimg.getExtend(); @@ -1620,38 +1635,52 @@ ErrorCode MergeRasterProcess(QVector filepaths, QString outfileptah, QS gdalImage resultImage = CreategdalImage(outfileptah, height, width, mainimg.band_num, rgt, EPSGCode, mainType, true, true, isENVI); QString resultMaskString = addMaskToFileName(outfileptah, QString("_MASK")); - gdalImage maskImage = CreategdalImage(resultMaskString, height, width,1 , rgt, EPSGCode, GDT_UInt16, true, true, isENVI); + gdalImage maskImage = CreategdalImage(resultMaskString, height, width,1 , rgt, EPSGCode, GDT_Int32, true, true, isENVI); // 初始化 - - // 分块计算 - long resultline = Memory1MB * 100 / 8 / resultImage.width; - resultline = resultline < 100 ? resultline : 100; // 最多100行 - + long resultline = Memory1MB * 500 / 8 / resultImage.width; + resultline = resultline < 3000 ? resultline : 3000; // 最多100行 + resultline = resultline > 0 ? resultline : 2; long bandnum = resultImage.band_num + 1; long starti = 0; long rasterCount = imgdslist.count(); - + + + + QProgressDialog progressDialog(u8"初始化影像", u8"终止", 0, resultImage.height); + progressDialog.setWindowTitle(u8"初始化影像"); + progressDialog.setWindowModality(Qt::WindowModal); + progressDialog.setAutoClose(true); + progressDialog.setValue(0); + progressDialog.setMaximum(resultImage.height); + progressDialog.setMinimum(0); + progressDialog.show(); + + for (starti = 0; starti < resultImage.height; starti = starti + resultline) { long blocklines = resultline; blocklines = starti + blocklines < resultImage.height ? blocklines : resultImage.height - starti; for (long b = 1; b < bandnum; b++) { Eigen::MatrixXd data = resultImage.getData(starti, 0, blocklines, resultImage.width, b); - Eigen::MatrixXd maskdata = maskImage.getData(starti, 0, blocklines, resultImage.width, b); + Eigen::MatrixXi maskdata = maskImage.getDatai(starti, 0, blocklines, resultImage.width, b); data = data.array() * 0; maskdata = maskdata.array() * 0; resultImage.saveImage(data, starti, 0, b); maskImage.saveImage(maskdata, starti, 0, b); } + if (nullptr != dia) { + dia->showProcess(starti * 1.0 / resultImage.height, u8"初始化影像数据"); + } + progressDialog.setValue(starti+ blocklines); } - + progressDialog.close(); switch (mergecode) { case MERGE_GEOCODING: - return MergeRasterInGeoCoding(imgdslist, resultImage, maskImage); + return MergeRasterInGeoCoding(imgdslist, resultImage, maskImage, dia); default: break; } @@ -1660,7 +1689,7 @@ ErrorCode MergeRasterProcess(QVector filepaths, QString outfileptah, QS return ErrorCode::SUCCESS; } -ErrorCode MergeRasterInGeoCoding(QVector imgdslist, gdalImage resultimg, gdalImage maskimg) +ErrorCode MergeRasterInGeoCoding(QVector imgdslist, gdalImage resultimg, gdalImage maskimg, ShowProessAbstract* dia) { omp_set_num_threads(Paral_num_thread); // 逐点合并计算 @@ -1680,12 +1709,22 @@ ErrorCode MergeRasterInGeoCoding(QVector imgdslist, gdalImage resulti long processNumber = 0; + QProgressDialog progressDialog(u8"合并影像", u8"终止", 0, resultimg.height); + progressDialog.setWindowTitle(u8"合并影像"); + progressDialog.setWindowModality(Qt::WindowModal); + progressDialog.setAutoClose(true); + progressDialog.setValue(0); + progressDialog.setMaximum(resultimg.height); + progressDialog.setMinimum(0); + progressDialog.show(); + + + + omp_lock_t lock; omp_init_lock(&lock); - - -#pragma omp parallel for +//#pragma omp parallel for for (starti = 0; starti < resultimg.height; starti = starti + resultline) { long blocklines = resultline; blocklines = starti + blocklines < resultimg.height ? blocklines : resultimg.height - starti; @@ -1700,9 +1739,9 @@ ErrorCode MergeRasterInGeoCoding(QVector imgdslist, gdalImage resulti long minRid = imgdslist[ir].height; long maxRid = 0; - - Eigen::MatrixXd ridlist = Eigen::MatrixXd::Zero(blocklines, resultimg.width); - Eigen::MatrixXd cidlist = Eigen::MatrixXd::Zero(blocklines, resultimg.width); + Eigen::MatrixXd ridlist = resultimg.getData(starti, 0, blocklines, resultimg.width, 1); + ridlist = ridlist.array() * 0; + Eigen::MatrixXd cidlist = ridlist.array() * 0; for (long i = 0; i < blocklines; i++) {// 行号 rid = starti + i; @@ -1764,6 +1803,14 @@ ErrorCode MergeRasterInGeoCoding(QVector imgdslist, gdalImage resulti r0 = ridtemp - std::floor(ridtemp); c0 = cidtemp - std::floor(cidtemp); + lastr = std::floor(ridtemp); + nextr = std::ceil(ridtemp); + lastc = std::floor(cidtemp); + lastc = std::ceil(cidtemp); + + lastr = lastr - minRid; + nextr = nextr - minRid; + p0 = Landpoint{ c0,r0,0 }; p11 = Landpoint{ 0,0,data(lastr,lastc) }; p21 = Landpoint{ 0,1,data(nextr,lastc) }; @@ -1778,16 +1825,22 @@ ErrorCode MergeRasterInGeoCoding(QVector imgdslist, gdalImage resulti maskimg.saveImage(resultmask, starti, 0, b); } } - omp_set_lock(&lock); processNumber = processNumber + blocklines; std::cout << "\rprocess bar:\t" << processNumber * 100.0 / resultimg.height << " % " << "\t\t\t"; - omp_unset_lock(&lock); // 锟酵放伙拷斤拷 + if (nullptr != dia) { + dia->showProcess(processNumber * 1.0 / resultimg.height, u8"合并图像"); + } + if (progressDialog.maximum() >= processNumber) { + processNumber = progressDialog.maximum() - 1; + } + progressDialog.setValue(processNumber); + omp_unset_lock(&lock); } - omp_destroy_lock(&lock); // 劫伙拷斤拷 + omp_destroy_lock(&lock); std::cout << std::endl; + progressDialog.setWindowTitle(u8"影像掩膜"); - // 全体处理 for (starti = 0; starti < resultimg.height; starti = starti + resultline) { long blocklines = resultline; blocklines = starti + blocklines < resultimg.height ? blocklines : resultimg.height - starti; @@ -1798,6 +1851,7 @@ ErrorCode MergeRasterInGeoCoding(QVector imgdslist, gdalImage resulti for (long i = 0; i < data.rows(); i++) { for (long j = 0; j < data.cols(); j++) { if (maskdata(i, j) == 0) { + data(i, j) = -9999; continue; } data(i, j) = data(i, j) / maskdata(i, j); @@ -1807,9 +1861,10 @@ ErrorCode MergeRasterInGeoCoding(QVector imgdslist, gdalImage resulti resultimg.saveImage(data, starti, 0, b); maskimg.saveImage(maskdata, starti, 0, b); } + progressDialog.setValue(starti + blocklines); } - - + resultimg.setNoDataValue(-9999); + progressDialog.close(); return ErrorCode::SUCCESS; } @@ -1862,6 +1917,9 @@ void gdalImageComplex::setData(Eigen::MatrixXcd 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); @@ -2183,6 +2241,8 @@ long GetEPSGFromRasterFile(QString filepath) return -1; } + std::cout << pszProjection << std::endl; + long epsgCode = atoi(oSRS.GetAuthorityCode(nullptr)); // 获取EPSG代码 if (epsgCode != 0) { @@ -2193,6 +2253,7 @@ long GetEPSGFromRasterFile(QString filepath) else { qDebug() << "EPSG code could not be determined from the spatial reference."; GDALClose(poDataset); + return -1; } } @@ -2303,6 +2364,11 @@ CoordinateSystemType getCoordinateSystemTypeByEPSGCode(long epsg_code) } } +void ShowProessAbstract::showProcess(double precent, QString tip) +{ +} - +void ShowProessAbstract::showToolInfo(QString tip) +{ +} diff --git a/ImageOperatorBase.h b/ImageOperatorBase.h index 24efe21..0e95baf 100644 --- a/ImageOperatorBase.h +++ b/ImageOperatorBase.h @@ -74,6 +74,15 @@ enum GDALREADARRCOPYMETHOD { +class ShowProessAbstract { +public: + virtual void showProcess(double precent,QString tip); + virtual void showToolInfo( QString tip) ; + + +}; + + /// 根据经纬度获取 /// EPSG代码，根据经纬度返回对应投影坐标系统，其中如果在中华人民共和国境内，默认使用 @@ -97,14 +106,7 @@ long GetEPSGFromRasterFile(QString filepath); std::shared_ptr GetCenterPointInRaster(QString filepath); CoordinateSystemType getCoordinateSystemTypeByEPSGCode(long EPSGCODE); - - - - - - - - + // 文件打开 // 当指令销毁时，调用GDALClose 销毁类型 std::shared_ptr OpenDataset(const QString& in_path, GDALAccess rwmode= GA_ReadOnly); @@ -248,13 +250,10 @@ enum MERGEMODE }; -ErrorCode MergeRasterProcess(QVector filepath, QString outfileptah, QString mainString, MERGEMODE mergecode = MERGEMODE::MERGE_GEOCODING, bool isENVI = false); - - -ErrorCode MergeRasterInGeoCoding(QVector inimgs, gdalImage resultimg,gdalImage maskimg); - +ErrorCode MergeRasterProcess(QVector filepath, QString outfileptah, QString mainString, MERGEMODE mergecode = MERGEMODE::MERGE_GEOCODING, bool isENVI = false, ShowProessAbstract* dia=nullptr); +ErrorCode MergeRasterInGeoCoding(QVector inimgs, gdalImage resultimg,gdalImage maskimg, ShowProessAbstract* dia = nullptr); template std::shared_ptr readDataArr(gdalImage &img,int start_row, int start_col, int rows_count, int cols_count, int band_ids, GDALREADARRCOPYMETHOD method);