最近一直在做大气方面相关遥感数据的影像处理,由于大气遥感数据的幅宽一般较大,原始数据一般会自带经纬度文件,手动一景一景处理效率极低,故编写IDL代码进行处理。那么这篇文章就带各位在面对不同大气遥感数据时候,如何快速对该数据的处理流程有个大致了解。
其实在我这篇博客里面https://mp.csdn.net/mp_blog/creation/editor/107604176 已经对自带经纬度文件的遥感数据处理在ENVI里面的流程大致叙述了一遍。如果有还不知道的小伙伴,可以翻看一下。主要就是:
首先你得知道需要处理数据的格式,一般有两种格式–NC文件和HDF文件,知道了数据格式,就知道如何在ENVI里面打开数据,如下图:
因为哨兵5P数据格式为NetCDF-4格式,所以使用 Open as ---- Generic Formats ---- HDF5/NetCDF4方式打开nc数据。接下来就是寻找经纬度文件,一般就找 latitude or longitude or lat or lon 这几个关键字就能找到,像S5P里面就是 latitude和longitude。
那么在IDL代码里面是如何打开经纬度数据呢?
先上代码:
file = 'G:\CH4\S5P_OFFL_L2__CH4____20200529T043518_20200529T061648_13603_01_010302_20200530T210438.nc'
latitude= e.openraster(file,DATASET_NAME='/PRODUCT/latitude')
longitude= e.openraster(file,DATASET_NAME='/PRODUCT/longitude')
这里主要用到了 openraster这个函数里面的一个参数---- DATASET_NAME,先看一下官方解释:
这里可以去看一下NC格式文件和HDF文件时如何存储数据的。可以简单这么理解:把NC数据看成一个大文件夹,这个大文件夹下面可能会有若干级子文件夹,DATASET_NAME就是去寻找底层文件夹里面放置的数据,相当于一个相对路径。回过头看我代码:
latitude= e.openraster(file,DATASET_NAME='/PRODUCT/latitude')
已经知道DATASET_NAME里面写的是一个路径,路径由两部分组成: Path 和 Name。这里纬度文件的Path 就是下图里面 Path后面那一串,Name 就是 Path上面有个Name后面那一串。经度数据和目标数据都是同理打开。我这里处理的是哨兵5P——甲烷数据,目标数据集就为methane_mixing_ratio
FILE = 'G:\CH4\S5P_OFFL_L2__CH4____20200529T043518_20200529T061648_13603_01_010302_20200530T210438.nc'
raster=e.openraster(file,DATASET_NAME='/PRODUCT/methane_mixing_ratio') ;打开目标数据集
现在已经找到了经纬度文件,那么构建GLT文件就很简单啦,直接上代码:
latitude= e.openraster(file,DATASET_NAME='/PRODUCT/latitude')
longitude= e.openraster(file,DATASET_NAME='/PRODUCT/longitude')
out_name=e.gettemporaryfilename() ;选择创建临时文件输出路径
lat_id=ENVIRASTERTOFID(latitude)
lon_id=ENVIRASTERTOFID(longitude)
envi_file_query,lat_id,nb=nb_lat
y_pos=lindgen(nb_lat)
envi_file_query,lon_id,nb=nb_lon
x_pos=lindgen(nb_lon) ;提取需要用到的参数
input_prj=ENVI_PROJ_CREATE(/geographic,datum='WGS-84')
output_prj=ENVI_PROJ_CREATE(/geographic,datum='WGS-84') ;设置输入输出的坐标系为 WGS-84坐标系
ENVI_DOIT,'ENVI_GLT_DOIT',dims=dims,I_PROJ=input_prj,O_PROJ=output_prj,pixel_size=0.08113,$
OUT_NAME=out_name,ROTATION=0,X_FID=lon_id,$
X_POS=x_pos[0], Y_FID=lat_id, Y_POS=y_pos[0]
这里需要注意一点就是在ENVI_DOIT里面有个参数是pixel_size,这个参数是需要自己去指定的,那么怎么指定呢?我一般是手动在ENVI里面构建GLT文件时,在设置GLT参数时,会自己有个默认值。也可以自己算:公式如下:这个有点不准,还是建议手动
1 ° 111.11 K M = 卫 星 分 辨 率 ( 十 进 制 度 ) 卫 星 分 辨 率 ( K M ) \frac{1°}{111.11KM}=\frac{卫星分辨率(十进制度)}{卫星分辨率(KM)} 111.11KM1°=卫星分辨率(KM)卫星分辨率(十进制度)
这里主要就是找到目标raster,从数据集中找目标raster和找经纬度文件是一样的方法,这里就不赘述了。校正部分代码如下:
glt_raster = e.openraster(out_name) ;打开GLT文件
glt_fid=envirastertofid(glt_raster)
envi_file_query,glt_fid,dims=dims
raster=e.openraster(file[i],DATASET_NAME='/PRODUCT/methane_mixing_ratio') ;打开目标数据集
raster_id=envirastertofid(raster)
envi_file_query,raster_id,ns=ns,nb=nb
pos=indgen(nb)
envi_doit,'ENVI_GEOREF_FROM_GLT_DOIT',FID=raster_id,GLT_DIMS=dims,GLT_FID=glt_fid,out_name=out_path,pos=pos
一开始我在主pro里面(就是从头执行到底不插入任何函数)进行上述所有操作,但是频繁给我报错,如下图:
总结下来就是,在for循环里面直接去构建GLT文件这种方法在处理第二景的时候就会报错,有时候提示GLT文件不存在,有时候提示不能连接到指定文件。
然后我便开始找原因,我觉得,可能是第一景文件的GLT文件还存在,导致后续创建GLT文件时,IDL认为GLT文件已存在,所以在处理第二景的时候, ENVI_DOIT,‘ENVI_GLT_DOIT’, 这一步直接给你跳过了,故在后续处理第二景的时候,读取GLT文件的时候读取不到。发现了问题,那就解决问题。按照常理想,你既然以为这个GLT文件存在,那我就在你用完这个GLT文件之后,也就是第一景处理完之后,把这个GLT文件关闭应该就可以了,但是还是不能解决问题。。。附上会报错的代码
pro CH4_deal
compile_opt idl2
e=envi(/headless)
time = systime(1)
path = 'G:\CH4\Oridata'
file = file_search(path,'*.nc',count=n)
for i=0,n-1 do begin
name = file_basename(file[i],'.nc4')
print,name
out_path='G:\CH4\result\'+name+'_result.dat'
; 重复处理数据不进行处理
if file_test(out_path) eq 1 then begin
print,'已处理'
endif else begin
;构建GLT文件
latitude= e.openraster(file[i],DATASET_NAME='/PRODUCT/latitude')
longitude= e.openraster(file[i],DATASET_NAME='/PRODUCT/longitude')
out_name=e.gettemporaryfilename()
lat_id=ENVIRASTERTOFID(latitude)
lon_id=ENVIRASTERTOFID(longitude)
envi_file_query,lat_id,nb=nb_lat
y_pos=lindgen(nb_lat)
envi_file_query,lon_id,nb=nb_lon
x_pos=lindgen(nb_lon)
input_prj=ENVI_PROJ_CREATE(/geographic,datum='WGS-84')
output_prj=ENVI_PROJ_CREATE(/geographic,datum='WGS-84')
ENVI_DOIT,'ENVI_GLT_DOIT',dims=dims,I_PROJ=input_prj,O_PROJ=output_prj,pixel_size=0.08113,$
OUT_NAME=out_name,ROTATION=0,X_FID=lon_id,$
X_POS=x_pos[0], Y_FID=lat_id, Y_POS=y_pos[0]
; 开始进行几何校正
glt_raster = e.openraster(out_name)
glt_fid=envirastertofid(glt_raster)
envi_file_query,glt_fid,dims=dims
raster=e.openraster(file[i],DATASET_NAME='/PRODUCT/methane_mixing_ratio')
raster_id=envirastertofid(raster)
envi_file_query,raster_id,ns=ns,nb=nb
pos=indgen(nb)
envi_doit,'ENVI_GEOREF_FROM_GLT_DOIT',FID=raster_id,GLT_DIMS=dims,GLT_FID=glt_fid,out_name=out_path,pos=pos
all_time = systime(2)
print,'共花费时间:',(all_time-time)/60
endelse
endfor
end
然后我就尝试把构建GLT文件的代码封装成了一个函数,结果可以执行,哦了,问题解决~,上代码:
pro CH4_deal
compile_opt idl2
e=envi(/headless)
time = systime(1)
path = 'G:\CH4\Oridata'
file = file_search(path,'*.nc',count=n)
for i=0,n-1 do begin
name = file_basename(file[i],'.nc4')
print,name
out_path='G:\CH4\result\'+name+'_result.dat'
if file_test(out_path) eq 1 then begin
print,'已处理'
endif else begin
out_name = build_glt(file[i],name)
glt_raster = e.openraster(out_name)
glt_fid=envirastertofid(glt_raster)
envi_file_query,glt_fid,dims=dims
raster=e.openraster(file[i],DATASET_NAME='/PRODUCT/methane_mixing_ratio')
raster_id=envirastertofid(raster)
envi_file_query,raster_id,ns=ns,nb=nb
pos=indgen(nb)
envi_doit,'ENVI_GEOREF_FROM_GLT_DOIT',FID=raster_id,GLT_DIMS=dims,GLT_FID=glt_fid,out_name=out_path,pos=pos
all_time = systime(2)
print,'共花费时间:',(all_time-time)/60
endelse
endfor
end
function build_glt,file,name
compile_opt idl2
e=envi(/headless)
latitude= e.openraster(file,DATASET_NAME='/PRODUCT/latitude')
longitude= e.openraster(file,DATASET_NAME='/PRODUCT/longitude')
out_name=e.gettemporaryfilename()
lat_id=ENVIRASTERTOFID(latitude)
lon_id=ENVIRASTERTOFID(longitude)
envi_file_query,lat_id,nb=nb_lat
y_pos=lindgen(nb_lat)
envi_file_query,lon_id,nb=nb_lon
x_pos=lindgen(nb_lon)
input_prj=ENVI_PROJ_CREATE(/geographic,datum='WGS-84')
output_prj=ENVI_PROJ_CREATE(/geographic,datum='WGS-84')
ENVI_DOIT,'ENVI_GLT_DOIT',dims=dims,I_PROJ=input_prj,O_PROJ=output_prj,pixel_size=0.081130,$
OUT_NAME=out_name,ROTATION=0,X_FID=lon_id,$
X_POS=x_pos[0], Y_FID=lat_id, Y_POS=y_pos[0]
return,out_name
end
如有错误,欢迎指正!
文章浏览阅读1k次。通过使用ajax方法跨域请求是浏览器所不允许的,浏览器出于安全考虑是禁止的。警告信息如下:不过jQuery对跨域问题也有解决方案,使用jsonp的方式解决,方法如下:$.ajax({ async:false, url: 'http://www.mysite.com/demo.do', // 跨域URL ty..._nginx不停的xhr
文章浏览阅读2k次。关于在 Oracle 中配置 extproc 以访问 ST_Geometry,也就是我们所说的 使用空间SQL 的方法,官方文档链接如下。http://desktop.arcgis.com/zh-cn/arcmap/latest/manage-data/gdbs-in-oracle/configure-oracle-extproc.htm其实简单总结一下,主要就分为以下几个步骤。..._extproc
文章浏览阅读1.5w次。linux下没有上面的两个函数,需要使用函数 mbstowcs和wcstombsmbstowcs将多字节编码转换为宽字节编码wcstombs将宽字节编码转换为多字节编码这两个函数,转换过程中受到系统编码类型的影响,需要通过设置来设定转换前和转换后的编码类型。通过函数setlocale进行系统编码的设置。linux下输入命名locale -a查看系统支持的编码_linux c++ gbk->utf8
文章浏览阅读750次。今天准备从生产库向测试库进行数据导入,结果在imp导入的时候遇到“ IMP-00009:导出文件异常结束” 错误,google一下,发现可能有如下原因导致imp的数据太大,没有写buffer和commit两个数据库字符集不同从低版本exp的dmp文件,向高版本imp导出的dmp文件出错传输dmp文件时,文件损坏解决办法:imp时指定..._imp-00009导出文件异常结束
文章浏览阅读143次。当下是一个大数据的时代,各个行业都离不开数据的支持。因此,网络爬虫就应运而生。网络爬虫当下最为火热的是Python,Python开发爬虫相对简单,而且功能库相当完善,力压众多开发语言。本次教程我们爬取前程无忧的招聘信息来分析Python程序员需要掌握那些编程技术。首先在谷歌浏览器打开前程无忧的首页,按F12打开浏览器的开发者工具。浏览器开发者工具是用于捕捉网站的请求信息,通过分析请求信息可以了解请..._初级python程序员能力要求
文章浏览阅读7.6k次,点赞2次,收藏6次。@Service标注的bean,类名:ABDemoService查看源码后发现,原来是经过一个特殊处理:当类的名字是以两个或以上的大写字母开头的话,bean的名字会与类名保持一致public class AnnotationBeanNameGenerator implements BeanNameGenerator { private static final String C..._@service beanname
文章浏览阅读6.9w次,点赞73次,收藏463次。1.前序创建#include<stdio.h>#include<string.h>#include<stdlib.h>#include<malloc.h>#include<iostream>#include<stack>#include<queue>using namespace std;typed_二叉树的建立
文章浏览阅读7.1k次。在Asp.net上使用Excel导出功能,如果文件名出现中文,便会以乱码视之。 解决方法: fileName = HttpUtility.UrlEncode(fileName, System.Text.Encoding.UTF8);_asp.net utf8 导出中文字符乱码
文章浏览阅读2.1k次,点赞4次,收藏23次。第一次实验 词法分析实验报告设计思想词法分析的主要任务是根据文法的词汇表以及对应约定的编码进行一定的识别,找出文件中所有的合法的单词,并给出一定的信息作为最后的结果,用于后续语法分析程序的使用;本实验针对 PL/0 语言 的文法、词汇表编写一个词法分析程序,对于每个单词根据词汇表输出: (单词种类, 单词的值) 二元对。词汇表:种别编码单词符号助记符0beginb..._对pl/0作以下修改扩充。增加单词
文章浏览阅读773次。我在使用adb.exe时遇到了麻烦.我想使用与bash相同的adb.exe shell提示符,所以我决定更改默认的bash二进制文件(当然二进制文件是交叉编译的,一切都很完美)更改bash二进制文件遵循以下顺序> adb remount> adb push bash / system / bin /> adb shell> cd / system / bin> chm..._adb shell mv 权限
文章浏览阅读6.8k次,点赞12次,收藏125次。1. 单目相机标定引言相机标定已经研究多年,标定的算法可以分为基于摄影测量的标定和自标定。其中,应用最为广泛的还是张正友标定法。这是一种简单灵活、高鲁棒性、低成本的相机标定算法。仅需要一台相机和一块平面标定板构建相机标定系统,在标定过程中,相机拍摄多个角度下(至少两个角度,推荐10~20个角度)的标定板图像(相机和标定板都可以移动),即可对相机的内外参数进行标定。下面介绍张氏标定法(以下也这么称呼)的原理。原理相机模型和单应矩阵相机标定,就是对相机的内外参数进行计算的过程,从而得到物体到图像的投影_相机-投影仪标定
文章浏览阅读2.2k次。文章目录Wayland 架构Wayland 渲染Wayland的 硬件支持简 述: 翻译一篇关于和 wayland 有关的技术文章, 其英文标题为Wayland Architecture .Wayland 架构若是想要更好的理解 Wayland 架构及其与 X (X11 or X Window System) 结构;一种很好的方法是将事件从输入设备就开始跟踪, 查看期间所有的屏幕上出现的变化。这就是我们现在对 X 的理解。 内核是从一个输入设备中获取一个事件,并通过 evdev 输入_wayland