GDAL求取LAI叶面积指数

using System;
using System.Collections.Generic;
using System.ComponentModel;
using System.Data;
using System.Drawing;
using System.Linq;
using System.Text;
using System.Windows.Forms;
using OSGeo.GDAL;
using System.IO;
using System.Diagnostics;

namespace Yaan_AppSysWinForm
{
    public partial class LaiForm : Form
    {
        public LaiForm()
        {
            InitializeComponent();
        }

        string fileNameOnly = "";

        private void btnInputImage_Click(object sender, EventArgs e)
        {
            OpenFileDialog ofd = new OpenFileDialog();
            if (ofd.ShowDialog() == DialogResult.OK)
            {
                tbInputImage.Text = ofd.FileName;
            }
            int index = tbInputImage.Text.LastIndexOf("\\");
            fileNameOnly = tbInputImage.Text.Substring(index + 1);
            
        }

        private void btnInputClassified_Click(object sender, EventArgs e)
        {
            OpenFileDialog ofd = new OpenFileDialog();
            if (ofd.ShowDialog() == DialogResult.OK)
            {
                tbInputClassified.Text = ofd.FileName;
            }
        }

        private void btnInputTxt_Click(object sender, EventArgs e)
        {
            OpenFileDialog ofd = new OpenFileDialog();
            if (ofd.ShowDialog() == DialogResult.OK)
            {
                tbInputTxt.Text = ofd.FileName;
            }
        }

        private void btnOutput_Click(object sender, EventArgs e)
        {
            string saveFolder = "";
            FolderBrowserDialog fbd = new FolderBrowserDialog();
            if (fbd.ShowDialog() == DialogResult.OK)
            {
                saveFolder = fbd.SelectedPath;
            }

            int index = fileNameOnly.LastIndexOf(".");
            tbOutput.Text = saveFolder + "\\" + "lai_result.tif";
        }

        private void btnCancel_Click(object sender, EventArgs e)
        {
            //MessageBox.Show("取消");
            this.Close();
        }

        private void btnExecute_Click(object sender, EventArgs e)
        {
            Gdal.AllRegister();
            Gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8", "YES");

            if (tbInputImage.Text == "")
            {
                MessageBox.Show("影像路径不能为空");
                return;
            }

            Dataset dsSrc = Gdal.Open(tbInputImage.Text, Access.GA_ReadOnly);//原影像

            int imgWidth = dsSrc.RasterXSize;    //获得影像宽
            int imgHeight = dsSrc.RasterYSize;   //获得影像高
            int bandNum = dsSrc.RasterCount;     //获得波段数量

            //原影像的灰度值
            double[] dSrcArray = new double[imgHeight * imgWidth];

            Band bandSrc = dsSrc.GetRasterBand(1);

            //将原图像的像元值读取到dSrcArray数组中
            bandSrc.ReadRaster(0, 0, imgWidth, imgHeight, dSrcArray, imgWidth, imgHeight, 0, 0);

            #region 验证原影像的灰度值数据--已注释
            //for (int ii = 0; ii < dSrcArray.Length; ++ii)
            //{
            //    if (ii % 40 == 0)
            //    {
            //        Debug.WriteLine("\n");
            //    } 
            //    Debug.Write(dSrcArray[ii].ToString() + "  ");
            //}
            #endregion

            //分类影像
            Dataset dsClassified = Gdal.Open(tbInputClassified.Text, Access.GA_ReadOnly);
            int classWidth = dsClassified.RasterXSize;//获得分类影像宽度
            int classHeight = dsClassified.RasterYSize;//获得分类影像高度
            int[] dClassifiedArray = new int[classWidth * classHeight];//存储分类影像的灰度值
            Band bandClassified = dsClassified.GetRasterBand(1);

            //将原图像的像元值读取到dClassifiedArray数组中
            bandClassified.ReadRaster(0, 0, classWidth, classHeight, dClassifiedArray, classWidth, classHeight, 0, 0);
            
            
            #region  验证分类影像的灰度值数据--已注释
            //for (int ii = 0; ii < dClassifiedArray.Length; ++ii)
            //{
            //    if (ii % 40 == 0)
            //    {
            //        Debug.WriteLine("\n");
            //    }
            //    Debug.Write(dClassifiedArray[ii].ToString() + "  ");
            //}
            #endregion


            //构建位图
            //Bitmap bmp = new Bitmap(classWidth, classHeight, System.Drawing.Imaging.PixelFormat.Format24bppRgb);
            //Bitmap b = new Bitmap(400, 400, System.Drawing.Imaging.PixelFormat.Format16bppGrayScale);


            #region 读取查找表数据,类型：typeOfShuiDao，typeOfYuMi，typeOfHuaSheng，typeOfLinDi，typeOfCaoDi，typeOfOther；回归系数：aShuiDao, aYuMi, aHuaSheng, aLinDi, aCaoDi, aOther，bShuiDao, bYuMi, bHuaSheng, bLinDi, bCaoDi, bOther

            /**************读取查找表****************/
            string[] lines = File.ReadAllLines(tbInputTxt.Text,Encoding.Default);
            int typeOfShuiDao = 1111;
            int typeOfYuMi = 1131;
            int typeOfHuaSheng = 1132;
            int typeOfLinDi = 131;
            int typeOfCaoDi = 141;
            int typeOfOther = 0;
            double aShuiDao, aYuMi, aHuaSheng, aLinDi, aCaoDi, aOther = 0;
            double bShuiDao, bYuMi, bHuaSheng, bLinDi, bCaoDi, bOther = 0;
           
            string[] temp1=lines[1].Split(new char[]{' ','\t'}, StringSplitOptions.RemoveEmptyEntries);
            aShuiDao = Convert.ToDouble(temp1[2]);
            bShuiDao = Convert.ToDouble(temp1[3]);

            string[] temp2 = lines[2].Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
            aYuMi = Convert.ToDouble(temp2[2]);
            bYuMi = Convert.ToDouble(temp2[3]);

            string[] temp3 = lines[3].Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
            aHuaSheng = Convert.ToDouble(temp3[2]);
            bHuaSheng = Convert.ToDouble(temp3[3]);

            string[] temp4 = lines[4].Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
            aLinDi = Convert.ToDouble(temp4[2]);
            bLinDi = Convert.ToDouble(temp4[3]);

            string[] temp5 = lines[5].Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
            aCaoDi = Convert.ToDouble(temp5[2]);
            bCaoDi = Convert.ToDouble(temp5[3]);

            string[] temp6 = lines[6].Split(new char[] { ' ', '\t' }, StringSplitOptions.RemoveEmptyEntries);
            aOther = Convert.ToDouble(temp6[2]);
            bOther = Convert.ToDouble(temp6[3]);

            Debug.WriteLine(aShuiDao.ToString() + "   " + bShuiDao.ToString());
            Debug.WriteLine(aYuMi.ToString() + "   " + bYuMi.ToString());
            Debug.WriteLine(aHuaSheng.ToString() + "   " + bHuaSheng.ToString());
            Debug.WriteLine(aLinDi.ToString() + "   " + bLinDi.ToString());
            Debug.WriteLine(aCaoDi.ToString() + "   " + bCaoDi.ToString());
            Debug.WriteLine(aOther.ToString() + "   " + bOther.ToString());

            #endregion

            double[] resultArray = new double[imgHeight * imgWidth];//将该数组的数据写入新建的空白tiff中

            /*************创建TIFF***************/
            //获取指定格式的驱动，用于创建图像
            string pszFormat = "GTiff";
            //GTIiff创建为tif图像，HFA创建为Erdas的img格式，ENVI创建为ENVI的hdr文件

            Driver oDriver = Gdal.GetDriverByName(pszFormat);
            if (oDriver == null)
            {
                Console.WriteLine("格式{0}不支持Create()方法。\n", pszFormat);
                return;
            }

            /********************新建tiff文件，注意数据格式********************/
            //Dataset oDs = oDriver.Create(tbOutput.Text, imgWidth, imgHeight, 1, DataType.GDT_Byte, null);
            
       <span style="font-size:14px;color:#FF0000;">   <strong>  //支持double数据写入tiff
            Dataset oDs = oDriver.Create(tbOutput.Text, imgWidth, imgHeight, 1, DataType.GDT_Float64, null);</strong></span>
            if (oDs == null)
            {
                MessageBox.Show("创建图像{0}失败。\n", tbOutput.Text);
                return;
            }
            Band pBand = oDs.GetRasterBand(1);
            

            for (int i = 0; i < imgWidth; i++)
            {
                for (int j = 0; j < imgHeight; j++)
                {
                    //分类处理
                    switch (dClassifiedArray[i + j * imgWidth])
                    {
                        case 1111:
                            if ((dSrcArray[i + j * imgWidth] * aShuiDao + bShuiDao )<= 0)
                            {
                                resultArray[i + j * imgWidth] = 0;
                            }
                            else
                            {
                                resultArray[i + j * imgWidth] = dSrcArray[i + j * imgWidth] * aShuiDao + bShuiDao;
                            }
                            break;
                        case 1131:
                            if ((dSrcArray[i + j * imgWidth]*aYuMi+bYuMi)<= 0)
                            {
                                resultArray[i + j * imgWidth] = 0;
                            }
                            else
                            {
                                resultArray[i + j * imgWidth] = dSrcArray[i + j * imgWidth] * aYuMi + bYuMi;
                            }
                            break;
                        case 1132:
                            if ((dSrcArray[i + j * imgWidth] * aHuaSheng + bHuaSheng) <= 0)
                            {
                                resultArray[i + j * imgWidth] = 0;
                            }
                            else
                            {
                                resultArray[i + j * imgWidth] = dSrcArray[i + j * imgWidth] * aHuaSheng + bHuaSheng;
                            }
                            break;
                        case 131:
                            if ((dSrcArray[i + j * imgWidth] * aLinDi + bLinDi) <= 0)
                            {
                                resultArray[i + j * imgWidth] = 0;
                            }
                            else
                            {
                                resultArray[i + j * imgWidth] = dSrcArray[i + j * imgWidth] * aLinDi + bLinDi;
                            }
                            break;
                        case 141:
                            if ((dSrcArray[i + j * imgWidth] * aCaoDi + bCaoDi) <= 0)
                            {
                                resultArray[i + j * imgWidth] = 0;
                            }
                            else
                            {
                                resultArray[i + j * imgWidth] = dSrcArray[i + j * imgWidth] * aCaoDi + bCaoDi;
                            }
                            break;
                        case 0:
                            resultArray[i + j * imgWidth] = dSrcArray[i + j * imgWidth] * aOther + bOther;
                            break;
                        default:
                            resultArray[i + j * imgWidth] = dSrcArray[i + j * imgWidth] * aOther + bOther;
                            break;
                    }
                }
            }

            Debug.WriteLine("Max:" + resultArray.Max().ToString());
            Debug.WriteLine("Min:" + resultArray.Min().ToString());

            pBand.WriteRaster(0, 0, imgWidth, imgHeight, resultArray, imgWidth, imgHeight, 0, 0);
            oDs.FlushCache();
            oDs.Dispose();
            MessageBox.Show("数据集生成成功!\n所在目录：" + tbOutput.Text);
        }
    }
}

成功！