01 数据类型和变量命名
数据类型:
常见的数据类型有浮点(32Bit)和存储位数更多的双精度浮点型(64Bit)、整型(16Bit)和长整型(32Bit),字符串型(与python类似并无单字符的数据类型),至于其它数据类型较不常见不予介绍;
变量命名(在IDL命令行生成):
由于IDL为解释型语言,因此无需进行变量类型的声明,且变量之间没有不可逾越的鸿沟,例如整型变量可以很轻易的重新赋值为字符串型变量,这在C语言和CPlusPlus中是不允许的。
1. IDL> ; 输出浮点数 2. IDL> num1 = 2.0 ; 只需要该数包含小数即默认创建浮点型变量 3. IDL> num2 = 3. ; 即使如此语法上也是支持的 4. IDL> num3 = 3.1415926 ; 正常的创建一个浮点数 5. IDL> print, num1 ; 输出num1 6. 2.00000 7. IDL> print, num2 ; 输出num2 8. 3.00000 9. IDL> print, num3 ; 输出num3 10. 3.14159 11. 12. IDL> ; 创建双精度浮点型变量 13. IDL> num1 = 3.2D ; 只需要在正常的小数末端添加D 14. IDL> num2 = 3.d ; 其实在小数末端添加d也是可以的 15. IDL> num3 = 3.1415926D ; 创建一个正常的浮点数,或许你已经注意到此前创建的num3并未全部显示,这是存储位数有限对于浮点数而言 16. IDL> print, num1 ; 输出num1 17. 3.2000000 18. IDL> print, num2 ; 输出num2 19. 3.0000000 20. IDL> print, num3 ; 输出num3,这次可以正常显示小数部分的全部了,因为双精度浮点数的存储位数更多了 21. 3.1415926 22. 23. IDL> ; 创建整型变量 24. IDL> num2 = 123S ; 你在整数末端添加S或者s指代该整数为整型变量 25. IDL> print, num2 ; 输出num2 26. 314 27. 28. IDL> ; 创建长整型变量 29. IDL> num1 = 123L ; 在正常的整数末端添加L或者l即可 30. IDL> print, num1 ; 输出num1 31. 123 32. 33. IDL> ; 若什么类型也不指定那么IDL依据大小进行数据类型的选择 34. IDL> num1 = 123 35. IDL> num2 = 123456789 36. IDL> num3 = 123456789123 37. IDL> help, num1 ; help函数获取当前变量的数据类型等基本信息 38. NUM1 INT = 123 39. IDL> help, num2 40. NUM2 LONG = 123456789 41. IDL> help, num3 42. NUM3 LONG64 = 123456789123 43. 44. IDL> ; 创建字符串 45. IDL> info = "hello, world" 46. IDL> print, info 47. hello, world 48. IDL> help, info 49. INFO STRING = 'hello, world'
02 类型转换
这里我们说明一些常见的类型转换,例如最最常见的整型转化为浮点型,浮点型转化为整型(这里又涉及小数部分向上向下取整部分);
1. IDL> ; 我们创建一个整型变量和浮点型变量 2. IDL> num1 = 1234 3. IDL> num2 = 12.34 4. IDL> ; 将整型转化为浮点型有以下几种方法 5. IDL> num1 = float(num1) ; 使用FLOAT函数,注意IDL不区分大小写(不管是变量还是函数) 6. IDL> num1 7. 1234.0000 8. IDL> num1 = 1234 ; 为了保证实验可靠性这里我们重新赋值 9. IDL> num1 = num1 * 1.0 ; 很简单,乘以1.0或者除以1.0都可,但是记得赋值回给num1 10. IDL> num1 11. 1234.0000 12. 13. IDL> ; 下面我们继续将浮点数转化为整型 14. IDL> num2 = fix(num2) ; 使用FIX函数可将浮点数转化为整数,注意其是向下取整 15. IDL> num2 16. 12 17. IDL> num2 = 12.34 ; 重新赋值继续其它函数的测验 18. IDL> num2 = ceil(num2) ; 使用ceil函数向上取整 19. IDL> num2 20. 13 21. IDL> num2 = 12.34 ; 重新赋值继续其它函数的测验 22. IDL> num2 = round(num2) ; 使用round函四舍五入取整 23. IDL> num2 24. 12
03 创建数组
这里仅列举常见数组的创建方法,例如整型数组(即数组中每一个数均为整型),浮点数组,索引数组(有序数字组成)
1. IDL> ; 创建整型数组 2. IDL> arr = intarr(3, 4) ; 表示创建3列4行的数组,默认初始化为0 3. IDL> arr 4. 0 0 0 5. 0 0 0 6. 0 0 0 7. 0 0 0 8. IDL> arr = intarr(3) ; 刚刚是生成二维数组,现在是生成一维数组 9. IDL> arr 10. 0 0 0 11. IDL> arr = intarr(2, 3, 4) ; 生成2行3列4个通道的三维数组 12. IDL> arr 13. 0 0 14. 0 0 15. 0 0 16. 17. 0 0 18. 0 0 19. 0 0 20. 21. 0 0 22. 0 0 23. 0 0 24. 25. 0 0 26. 0 0 27. 0 0 28. IDL> ; 至于创建浮点数组使用fltarr()函数即可,与上文类似,这里不再重复 29. IDL> ; 创建索引数组 30. IDL> arr = indgen(3, 4) ; 创建3行4列的整型顺序数组,默认从0开始 31. IDL> arr 32. 0 1 2 33. 3 4 5 34. 6 7 8 35. 9 10 11 36. IDL> arr = indgen(3, 4, /double, start=123) ; 还可以指定开始数字,以及数组元素的类型 37. IDL> help, arr 38. ARR DOUBLE = Array[3, 4] 39. IDL> arr 40. 123.00000000000000 124.00000000000000 125.00000000000000 41. 126.00000000000000 127.00000000000000 128.00000000000000 42. 129.00000000000000 130.00000000000000 131.00000000000000 43. 132.00000000000000 133.00000000000000 134.00000000000000 44.
04 读取数组和数组运算
4.1 读取数组
5. IDL> ; 创建一个整型数组 6. IDL> arr = indgen(3, 4) 7. IDL> ; 先查看以下吧 8. IDL> arr 9. 0 1 2 10. 3 4 5 11. 6 7 8 12. 9 10 11 13. IDL> ; 获取第2行第3列的数字 14. IDL> arr[1, 2] ; 首先是用[]进行元素的选择,且行列的索引从0开始,且以后我所说的第n行第m列也是指代索引为0开始的行列号 15. 7 16. IDL> ; 获取第2行第3列的数字(索引从0开始的行列) 17. IDL> ; 获取第2列第3行的数字(索引从0开始的行列) 18. IDL> arr[1, 2] ; 首先是用[]进行元素的选择,且行列的索引从0开始 19. 7 20. IDL> ; 获取第2列所有的数字 21. IDL> arr[2, *] ; 用*星号代替该行 22. 2 23. 5 24. 8 25. 11 26. IDL> ; 获取第1至2列,第2至3行的所有数字 27. IDL> arr[1:2, 2:3] ; 使用:表示从哪到哪, 且均为闭区间 28. 7 8 29. 10 11 30. IDL> ; 取第1列第2行,第2列第3行的两个数字 31. IDL> arr[[1, 2], [2, 3]] 32. 7 11 33. IDL> ; 上述涉及的取法可能比较少,需要在[]中嵌套[],两个[]对应起来 34. IDL> ; 很自然地我们想到是否可以使用多个 单个索引进行取值?当然可以 35. IDL> arr[[5, 6, 7, 8, 9]] ; 但是需要传入地是一个数组,这里也是后面一个例子的核心 36. 5 6 7 8 9 37.
4.2 数组运算
注意:数组大家很容易会和线性代数中的向量和矩阵等联系起来,这本没有任何问题,但是在IDL中的基础运算例如×,它与线性代数中的×不一样,IDL只是两个数组对应元素的加减乘除,而不是转置取逆等。
所以只要没有特殊说明,数组的运算均为对应位置的像元的算术运算。
加减乘除不说大家都明白,这里说一下其它较不常见的运算符(例如幂、取余);
1. IDL> ; 创建一个整型数组 2. IDL> arr1 = indgen(2, 3) 3. IDL> arr2 = indgen(2, 3, start=3) 4. IDL> arr1 5. 0 1 6. 2 3 7. 4 5 8. IDL> arr2 9. 3 4 10. 5 6 11. 7 8 12. IDL> arr1 ^ arr2 13. 0 1 14. 32 729 15. 16384 -2591 16. IDL> ; 可以发现依旧是对应像元的幂运算,即第一个数组arr1像元的m次方(m为第二数组arr2对应像元值)
05 求取多幅TIFF影像对应像元的最大值-输出为Geotiff文件
既然都学习了这么多基础知识,那么我们开始学习一下吧,关于IF和For循环语法、读取TIFF等相关内容你都可以在该pro文件中得到了解,为方便了解,这里尽量写的详细方便理解。
pro max_tiff ; 输入路径 ==> 有一点注意,既可以\亦可以/,但是在读取HDF文件时尽量使用/(因为HDF文件内组的路径只允许/),所以建议使用/ in_path = "D:/task/RemoteSencingImageProcessing_HuiChen/IDL/experiment_1/data1_tiff" ; 输出路径 ==> 还有一点注意,尽量避免中文路径 out_path = "D:/task/RemoteSencingImageProcessing_HuiChen/IDL\experiment_1/OutputData\max.tif" ; 获取以npp开头,后缀为tiff的文件,*表示任意个字符 files_path = FILE_SEARCH(in_path, "npp*.tif", count=files_amount) ; 返回字符串数组,每一个字符串都是检索到的TIFF文件的完整路径 ; 上述file_search函数中变量files_amount接收参数count返回的检索到的文件个数 ; 或许你会使用n_elements()函数获取files_path得到文件个数,但是不是特别推荐,一旦没有检索到任何文件个数那么files_path为空字符串,此时函数返回1 ==> 导致后面函数出现错误 ; 但是如果使用参数count返回那么files_amount得到的值为0,这就是区别 ; 所以我们应当添加一个if判断,如果没有检索TIFF文件并进行接下来的代码 ==> 当然相信就算没有这个if判断你也很快能发现没有检索到任何TIFF文件这个问题 if files_amount eq 0 then begin ; IF语法:if condition then begin,最后以endif结束 return ; 表示立刻退出函数或过程,从而不执行后续的所有代码 endif ; 读取tif文件 file_data = READ_TIFF(files_path[0], GEOTIFF=file_geotiff) ; 参数geotiff返回一个关于该TIFF文件投影信息的结构体,赋值给变量file_geotiff ; 读取第一个TIFF文件的size大小 file_size = size(file_data) ; 对于当前矩阵(二维),返回5个元素的数组,分别是[维度数, 列数, 行数, 矩阵元素的数据类型编码, 矩阵元素的总数] ; 初始化一个浮点矩阵,用于存储所有TIFF文件中对应像元中最大的那个值 max_arr = FLTARR(file_size[1], file_size[2]) ; 读取数组时从索引0开始,所以列数和行数的读取应该是索引1, 2 ; 循环处理每一个TIFF文件 for file_i = 0, files_amount - 1 do begin ; for循环的语法: for start(开始数字), end(结尾数字,包含其) do begin, 末尾是endfor ; 获取当前循环下的TIFF文件路径 file_path = files_path[file_i] ; 读取当前循环下的TIFF文件 file_data = READ_TIFF(file_path) ; 传入TIFF文件的路径获取该TIFF文件的栅格阵列 ; 以下内容实际上是判断当前TIFF文件每一个像元是否大于max_trr,若大于则替换否则不替换,但是逐像元的进行二重循环太过麻烦所以使用矩阵运算-会有少许难理解 ; 判断当前TIFF的file_data是否大于max_arr中对应像元值,若大于,返回其索引(索引有多个,因此以数组形式返回) position = WHERE(file_data gt max_arr) ; where函数十分好用,返回数组中非0的所有数的索引(以数组形式) ; 但是这些索引都不是大家想的行列号读取的那种二维索引,而是一维的。 ; 这里有一个bug,如果矩阵file_data中没有满足条件(file_data gt max_arr)那么该函数会返回-1(非数组形式); ; 或许你会这么写:if position eq -1 then continue ==> 总体思路没有问题,但是如果position是数组那么比较会出错 if position[0] eq -1 then continue ; 若返回数组那么数组中第一个数不会是-1而是正数,若是返回-1那么其是整型,获取第一个值其实就是它本身; ; 当然你也可以去判断position的维度 ; if SIZE(position, /N_DIMENSIONS) eq 1 then continue ; 将大于的像元值赋值到max_arr而不大于的不赋值 max_arr[position] = file_data[position] endfor ; 输出为TIFF文件 ==> write_tiff第一个参数为输出TIFF文件的路径, 第二个参数为输出TIFF文件的栅格阵列即矩阵或者说是数组 WRITE_TIFF, out_path, max_arr, GEOTIFF=file_geotiff, /float ; 记住输出数据类型为float,否则默认为整型, 另外将之前得到的file_geotiff传到参geotiff参数中 ; 也就是输出的tiff文件的投影信息就是输入的tiff文件的投影信息。 end
06 求取多幅TIFF影像对应像元的平均值-输出为Geotiff文件
这一个是求取平均值,由于实验没有要求,这里不详细讲述其中代码部分.
pro aver_tiff ; 输入路径 in_path = "D:\task\RemoteSencingImageProcessing_HuiChen\IDL\experiment_1\data1_tiff" out_path = "D:\task\RemoteSencingImageProcessing_HuiChen\IDL\experiment_1\OutputData\aver.tif" ; 获取in_path中所有的tiff文件 files_path = FILE_SEARCH(in_path, "*.tif", COUNT=files_amount) ; 存储求和的容器 data_box = fltarr(1035, 934) frequency_box = intarr(1035, 934) ; 进入循环对每一个tif文件进行处理 for i = 0, files_amount - 1 do begin ; 获取当前循环下的TIFF文件路径 file_path = files_path[i] ; 读取该TIFF文件 file_data = READ_TIFF(file_path, GEOTIFF=file_geotiff) ; 将该数据累加到box中 data_box += (file_data gt 0) * file_data frequency_box += (file_data gt 0) endfor ; 求取平均值 frequency_box = (frequency_box eq 0) + (frequency_box ne 0) * frequency_box mean_box = data_box / frequency_box ; 输出为Geotiff格式 WRITE_TIFF, out_path, mean_box, geotiff=file_geotiff, /float end