C/C++教程

C++ GDAL——在原始影像中根据云检图生成抠掉云的原始影像

本文主要是介绍C++ GDAL——在原始影像中根据云检图生成抠掉云的原始影像,对大家解决编程问题具有一定的参考价值,需要的程序猿们随着小编来一起学习吧!
#include <fstream>
#include <iostream>
#include <Windows.h>
#include "shlwapi.h"
#include <algorithm>
#include <direct.h> 
#include <imagehlp.h>
#pragma comment(lib,"imagehlp.lib")
#pragma comment(lib, "shlwapi.lib")
#include "CreatePyramid.h"
#include "gdal_priv.h"
#pragma warning(disable:4996)

#define RAM_SIZE  100

using namespace std;

int main(int argc, char* argv[])
{

    if (argc != 2)
    {
        cout << "argument is false" << endl;
        return 1;
    }
    string strError = argv[1]; int posError = strError.find_last_of("\\"); strError.erase(posError + 1); strError += "Error.log";
    std::ofstream finerror(strError, ios::app);
    FILE*pErrFile;
    if ((pErrFile = fopen(strError.c_str(), "r")) == NULL)
    {

        printf("Open Tsk File failed:\n%s\n", strError);

        return 1;

    }
    string strTsk = argv[1];
    std::ifstream fin(strTsk);
    if (!fin)
    {
        finerror << "ERROR: " << argv[1] << endl;
        std::cout << "Can not open ini file!  " << strTsk << std::endl;
        return 1;

    }


    string Argument; vector<string>strArg;
    while (getline(fin, Argument))
    {
        strArg.push_back(Argument);
    }
    string resultDom = strArg[2];
    int poResult = resultDom.find_last_of("\\");
    resultDom.erase(poResult+1);
    const char*resultDomPath = resultDom.c_str();
    if (!PathIsDirectory(resultDomPath))
    {
        MakeSureDirectoryPathExists((PCSTR)resultDomPath);
    }
    string TskLog = argv[1]; int posLog = TskLog.find_last_of("."); TskLog.erase(posLog); TskLog += ".log";
    ofstream LogWrite(TskLog);
    FILE*pLogFile;
    if ((pLogFile = fopen(TskLog.c_str(), "r")) == NULL)
    {

        printf("Open Tsk File failed:\n%s\n", TskLog);
        finerror << "ERROR: " << strArg [0]<< endl;
        return 1;

    }
    GDALDataset* poSrcDS;
    const char*strDom = strArg[0].c_str();
    GDALAllRegister();
    poSrcDS = (GDALDataset*)GDALOpen(strDom, GA_ReadOnly);
    if (poSrcDS == NULL)
    {
        LogWrite << strDom<<"  :no exit or errror" << endl;
        cout << "Can't Open the File and Image!" << endl;
        return EXIT_FAILURE;
    }
    LogWrite << strDom << "  Open succeed" << endl;
    int iBandCount = poSrcDS->GetRasterCount();
    int iSrcWidth = poSrcDS->GetRasterXSize();
    int iSrcHeight = poSrcDS->GetRasterYSize();

    const char* pszFormat = "GTiff";
    GDALDriver *poDriver = GetGDALDriverManager()->GetDriverByName(pszFormat);
    if (poDriver == NULL)
    {
        finerror << "ERROR: " << strArg[0] << endl;
        GDALClose((GDALDatasetH)poSrcDS);
        return 1;
    }
    double dGeoTransform[6];
    poSrcDS->GetGeoTransform(dGeoTransform);  
    const char* pszProj = poSrcDS->GetProjectionRef();


    GDALDataset* poSrcDS1;
    poSrcDS1 = (GDALDataset*)GDALOpen(strArg[1].c_str(), GA_ReadOnly);
    if (poSrcDS1 == NULL)
    {
        finerror << "ERROR: " << strArg[0] << endl;
        LogWrite << strArg[1]<< "  :no exit or errror" << endl;
        cout << "Can't Open the File and Image!" << endl;
        return EXIT_FAILURE;
    }
    LogWrite << strArg[1] << "  Open succeed" << endl;

    int iBandCount1 = poSrcDS1->GetRasterCount();
    int iSrcWidth1 = poSrcDS1->GetRasterXSize();
    int iSrcHeight1 = poSrcDS1->GetRasterYSize();


    const char* pszFormat1 = "GTiff";
    GDALDriver *poDriver1 = GetGDALDriverManager()->GetDriverByName(pszFormat1);
    if (poDriver1 == NULL)
    {
        finerror << "ERROR: " << strArg[0] << endl;
        GDALClose((GDALDatasetH)poSrcDS1);
        return 1;
    }
    double dGeoTransform1[6];
    poSrcDS1->GetGeoTransform(dGeoTransform1); 
    const char* pszProj1 = poSrcDS1->GetProjectionRef();


    GDALDataset* poDstDS;
    poDstDS = poDriver->Create(strArg[2].c_str(), iSrcWidth, iSrcHeight, iBandCount, GDT_UInt16, NULL);
    int nStepSize = (RAM_SIZE * 1024 * 1024) / (iSrcWidth * iBandCount);
    int nStepNum = iSrcHeight / nStepSize;
    if (iSrcHeight%nStepSize)nStepNum++;
    int pBandMap1[] = { 1 };
    int pBandMap[] = { 1, 2, 3, 4 };
    double nodata;
    nodata = poSrcDS->GetRasterBand(1)->GetNoDataValue();
    poDstDS->GetRasterBand(1)->SetNoDataValue(nodata);
    poDstDS->SetGeoTransform(dGeoTransform);
    poDstDS->SetProjection(pszProj);
    cout << "Nocld DOM Images are generating:" << strArg[2].c_str() << endl;
    for (int k = 0; k < nStepNum; k++)
    {
        cout << double(k*100 / nStepNum) << "%" << endl;
        int ybeg = max(0, min(nStepSize*k, iSrcHeight - 1));
        int yend = max(0, min(nStepSize*(k + 1), iSrcHeight));
        WORD *pImg = new WORD[(yend - ybeg)*iSrcWidth*iBandCount];       
        memset(pImg, 0, (yend - ybeg)*iSrcWidth*iBandCount*sizeof(WORD));
        BYTE *pImg1 = new BYTE[(yend - ybeg)*iSrcWidth*iBandCount1];
        memset(pImg1, 0, (yend - ybeg)*iSrcWidth*iBandCount1*sizeof(BYTE));
        poSrcDS->RasterIO(GF_Read, 0, ybeg, iSrcWidth, (yend - ybeg), pImg, iSrcWidth, (yend - ybeg), GDT_UInt16, iBandCount, pBandMap, 0, 0, 0);

        poSrcDS1->RasterIO(GF_Read, 0, ybeg, iSrcWidth1, (yend - ybeg), pImg1, iSrcWidth1, (yend - ybeg), GDT_Byte, iBandCount1, pBandMap1, 0, 0, 0);
        for (int j = 0; j < (yend - ybeg)*iSrcWidth*iBandCount1; ++j)
        {
            if (pImg1[j] == 255)
            {
                pImg[j] = 0;
                pImg[j + iSrcWidth*(yend - ybeg) * 1] = 0;
                pImg[j + iSrcWidth*(yend - ybeg) * 2] = 0;
                pImg[j + iSrcWidth*(yend - ybeg) * 3] = 0;
            }
        }
        poDstDS->RasterIO(GF_Write, 0, ybeg, iSrcWidth, (yend - ybeg), pImg, iSrcWidth, (yend - ybeg), GDT_UInt16, iBandCount, pBandMap, 0, 0, 0);
        delete[]pImg; pImg = NULL;
        delete[]pImg1; pImg1 = NULL;
    }
    cout << endl;
    GDALClose((GDALDatasetH)poSrcDS);
    GDALClose((GDALDatasetH)poSrcDS1); 
    GDALClose((GDALDatasetH)poDstDS);
    cout << "\nCreating Pyramids:" << endl;

    CConsoleProcess *pProgress = new CConsoleProcess();

    bool f = CreatePyramids(strArg[2].c_str(), pProgress);

    if (f == true)
    {
        cout << "***********************CreatePyramids Successed...***********************\n" << endl;
        LogWrite << "CreatePyramids Successed" << endl;
    }

    else
    {
        finerror << "ERROR: " << strArg[0] << endl;
        LogWrite << "CreatePyramids Failed" << endl;
        cout << "***********************CreatePyramids Failed...***********************\n" << endl;
    }
    delete pProgress;
    LogWrite.close();
    fclose(pLogFile);
    finerror.close();
    fclose(pErrFile);
    return 0;
}

/**
 * 在原始影像中根据云检图生成抠掉云的原始影像
 * 程序主要依赖的方法是gdal中rasterIO函数对tif影像的分块读写,使用性较强.
 */

 

这篇关于C++ GDAL——在原始影像中根据云检图生成抠掉云的原始影像的文章就介绍到这儿,希望我们推荐的文章对大家有所帮助,也希望大家多多支持为之网!