科学网

 找回密码
  注册
科学网 标签 校正 相关日志

tag 标签: 校正

相关日志

Fluke 17B 校正手记
dingsir 2019-12-7 21:49
2014年的时候,我花了大价钱买了一块好点的万用表,左挑右选相中了Fluke 的17B,到手后到也感觉不错,自动量程,3位半,小电压可以显示到毫伏(半位上的),平时测测感觉非常够用了。 过了两年,电池用完了,我拆开发现是5号电池用完了,于是给它换了电池。结果电池的簧片意外的断了。这么少的更换操作,簧片竟然断了,我心中一凉,这块表的质量应该不咋的。于是把剩下的簧片用尖嘴钳再弯一次,保证电池能够接触正常,算是解决了。 几年过去,前段时间也搞了块更精准的万用表,对比一下,发现这块17B已经很不准确了,测量一节锂电池,一个是3.8226,一个是3.784,相差感觉比较大(实际误差百分比1.01%,还算好,绝对误差在38.6mv),顿时对其精度产生了高度的怀疑。于是网上搜索了下这个型号的信息,发现这块表是可以校正的。于是找了一些背景资料,发现一些意外 1)原来这块表早期的版本是比较好的焊接芯片的版本, 后来廉价化之后,简化成了集成在电路板上的芯片(万用表行业戏称为“牛屎芯”),我的这块正是牛屎芯这种版本,怪不得连簧片质量也这么差劲,拆开拍个图对比一下。 为了校正,找到了一张指导校正的图片 好在电路板的布局还是相同的,可调电位器的位置和编号还是相同的。趁今天有空,就动手校正一下。先是打开直流电源,调了一下3.2V的电压,将输出的鳄鱼夹接到拆开的万用表的正负极上(因为要拆掉背面的胶盖,电池拆开来了不能使用电池供电)。 1)先断开直流电源,保证无电状态下,测试R34与R64的阻抗。这一步花了点时间,由于对这种贴片式可调电阻不熟悉,找两端的引脚费了点时间。后来弄明白了,就是一边有两个脚的(单脚的那一边是焊接用的)。测量R34+R64的电阻(用精确的表)在9.805MΩ,按指导说明,用小型的梅花起子调整R64将整体电池调整到10.00MΩ。 2)开启直流电源,用精确的表测量一节电池(我用了一节新的南孚碱性电池,用高精度表测量值为1.5952V。再切换到用17B测量其电压,用梅花起子调整R18,使之显示调整到1.595V,校正就结束了。再来测量下一节锂电池的电压,高精度表显示4.2409V,17B显示4.24V,比较准确了。再看看更大量程的,用直流电源设定一个12V的输出,高精度表显示为11.9959V,17B显示了12.01V,相差值14.1mv,与测量值的百分比为0.12%,足够准确了。 满意收工。 附一张17B的电路图,有兴趣的研究下,反正我是看不懂。 Fluke 17B circuit.rar
个人分类: 其它兴趣|11839 次阅读|0 个评论
“地基”探测的空间
xuxfyuwp 2018-4-22 19:56
在数值预报模式中,卫星资料所占比重越来越大,已达90%以上。且随着时间推移,仍在增加。而地基观测或地对空的观测资料占比则持续减少,“天基”明显占了上风,“地基”节节退守。这是否意味着地基观测的重要性也在减弱或可能逐步退出呢?结论应不会这么简单。 在EGU-2018大会上,来自英国里丁大学的学者介绍了一项工作,从一侧面表明地对空的探测依然重要。报告的题目是:“用全球测温激光雷达探测网矫正ECMWF再分析资料(ERA-5)和季节预测在平流层上层的偏差”。老外也能忽悠,题目既长又大,但注意了一下,所谓全球测温激光雷达探测网好像也就是由布设在欧美区域的五部激光雷达组成。 欧洲中期天气预报中心ECMWF的再分析资料应是国际一流水准,不少国家在天气气候分析研究中都在使用。但这些资料对大气的描述是否准确,需要做哪些改进,还要不断通过与实测资料的对比进行验证。 通过5部激光雷达资料的对比,可以看出无论是再分析资料还是季节预报模式在平流层高层都确实存在明显偏差。再分析资料的偏差在1-5度,月季预报偏差可达10度。 这些资料不仅可以验证再分析资料或季节预报存在的偏差,对于日常的数值模式天气预报的改进也同样重要。英国之所以要做这样的验证,是因为平流层的突然增温会给欧洲带来冷空气的影响,这可不是件小事,关系到灾害性天气预报准确率。 对已完成的成果不断进行评估和改进,这本是科技工作者的责任,而从这项工作可以得到的另一个启示是在大力发展卫星(天基)观测时,不能忽视地面或地对空观测的作用。从里丁大学给出的报告看,几个LIDAR组成的观测网范围虽还很有限,但仍能起到对数值场的订正作用,而中国区域还没有这类测站,这是需要引起关注的。 另一方面,连平流层高层的观测尚需要地基资料验证,近地面层的许多要素就更离不开实测对比了。发展好地基观测,既是做好预报的基础,也可以通过相互校正进一步促进卫星探测的精度,可起到一举两得的效果。 \0 \0
个人分类: 杂记|7084 次阅读|0 个评论
《陋室铭》升级版
热度 1 duke01361 2014-9-12 10:26
原版 陋室铭 山不在高,有仙则名。水不在深,有龙则灵。斯是陋室,惟吾德馨。苔痕上阶绿,草色入帘青。谈笑有鸿儒,往来无白丁。可以调素琴,阅金经。无丝竹之乱耳,无案牍之劳形。南阳诸葛庐,西蜀子云亭,孔子云:何陋之有? 升级版 陋室铭 山任其高,不攀不登。水何其深,我素我行。斯是斗室,惟吾独尊。苔痕上阶绿,草色入帘青。他乡有鸿儒,网来无白丁。可以谈新知,阅人生。无世俗之乱耳,无名利之劳形。京城“诸葛”庐,保定“子云”亭,老子云:何陋之有?
个人分类: Life a Bit|1886 次阅读|2 个评论
[转载]球差校正扫描透射电子显微镜
rongrongwang 2014-2-11 03:23
相差校正功能 近年来,透射电子显微镜(TEM)的进步显著,最为突出的是球差校正技术的开发。电子显微镜中使用的磁场透镜,原理上因为只能采用凸透镜制作,并不能像光学显微镜一样通过组合凹透镜改善像差。因此,镜片的各种像差,尤其是3层球面像差(Cs)的影响,分辨率会有所限制。但是最近几年来,利用理论上一直被倡导的多极子实现了凹透镜装置的制作,并凭借电稳定性的提高和调节中必要实验的数据积累等而最终被广泛应用。这种球差校正构造结合了TEM的透射和成像功能,实现了即便加速电压值200KV也无法实现的0.1nm的分辨率,使单个原子等级上的位置锁定和元素识别成为可能。 Fig.1中,根据照射体系中有无像差校正结构,电子束会显示不同的路线图。由于偏离透镜中心区域通过的电子束会发生大幅度折射,所以通常电子束不会在样品表面的一点聚焦,电子的探测范围被放大。另一方面,在凸镜的上方导入凹镜(像差校正结构),将会完全消除折射角度不同的现象。因此电子束在样品表面局部区域内聚焦,可以形成极细的电子探针。 另一方面,依靠像差校正,可以利用偏离透镜中心的高角度电子束,从而可以使用大型聚光器(集束)的透镜光圈。通常,聚光器透镜光圈只能使用到10mrad左右,但是通过像差校正后,20~40mrad的大型聚光器透镜光圈也可以使用。如Fig.2中所示,与像差校正前相比像差校正后的电子探测电流值增大了十倍以上。 (a) (b) Fig.1 Schematic diagram of Cs-TEM    (a) without Cs-corrected, (b) with Cs-corrected (a) (b) (c) Fig.2  The Intensity profile of electron-beam (a)Before correction、(b) After correction、(c)Before and after correction 原文章地址:http://www.toray.cn/trc/kinougenri/keitai/kei_002.html
个人分类: 电子显微镜|3992 次阅读|0 个评论
中国电力《三菱FX系列PLC编程与应用》(控制水电生产)
hncszwj 2013-8-2 10:06
【内容提要】 倒数第7行:“2013年1月第一版 2013年1月北京第一次印刷”→“2013年 7 月第一版 2013年 7 月北京第一次印刷”; P003:倒数12行,“西门子”→“SIEMENS”; P004:第13行、倒数第6行,两处“西门子”→“SIEMENS”; 倒数第3行,“0.1ms/千步”→“ 0.065 ms/千步”; P005:第17行,“ 控制系统( FCS, Fie l d D us Control System 将 逐步取代DCS ” →“ 控制系统 FCS ( Fie l d B us Control System ) 将 逐步取代DCS ”; P020:第8行,“ 数是PLC内部定时器 …… ” → “ 常 数是PLC内部定时器 …… ”; 第15行,“ 通信接口和扩展接口等含 一般基本单元I/O数16~128点……” → “通信接口和扩展接口等 。 一般基本单元I/O数16~128点……”; P023:两个图不清晰,如下: 图 1-10 图 1-11 P024:三个图不清晰,如下: 图 1-12 图 1-13 图 1-14 P028:倒数11行,“系列名称” →“系列 序号 ”: P030 :表 2-4中,“尺寸”列内,“mm3” →“mm 3 ”或“mm”(全书表中尺寸单位,都作同一处理); P035:第13行,“ 17种基本单元 ( CPU单元或主单元 单元 )” → “ 17种基本单元 ( CPU单元或主单元 )”; P037:倒数8行,“、128 6 种情况” → “、128 六 种情况”; P045:不看表格及注释的小5号字,第1行(5号字),“FX3G系列产品有 、 ” → “FX3G系列产品有 : ”; P055:第5行,“条使用的顺序程序” → “ 调试 用的顺序程序”; P056:倒数12行,“ FX2N-2AD / 2D A 和主单元用电缆在主单元的右边进行连接 ” →“ FX2N-2AD / 2D A 用电缆连接在主单元的右边 ”; 第三章开始,每页最下面页面序号应和前面一致,都用三位码,如“ 60” → “ 0 60”; P062:第3行,“运算出错” → “运算出错 (M8068) ”; P073:图3-29 P076:第13行, “ 编程软件说明书上与 动合称为A触点 相应 ,动断 被 称作B触点 ” → “ 有时将 动合 触点 称为A触点,动断 触点 称 为 B触点 ”; P085:图3-47的注释,b)“执行PLC 协 写入” → “执行PLC写入”; P089:第12行,“ 例 如,” →“如 图3-57所示 ,”; P092:图3-61如下 P103:第17行,“D1 2 D12(32位)” → “D1 3 D12(32位)”; P107:图4-5中左上部,“X 303 ” →“X 030 ”; P122:图4-26中,“SF RD P”→“SF WR P”; P136:图4-49中,梯形图画错 P145:倒数第9行(不计图4-71),“正弦S in 、余弦C os ” → “正弦S IN 、余弦C OS ”; P151:倒数第5行(不计图),段前空格有四格,修改为两格; P215:第9行,“41 OO MW” → “41 00 MW”; P217:第2行,“ 可操作 S K ” →“ 可操作 选择开关 S A ”、“按钮 1AN ” →“ 按钮 SB1 ”; 第3行,“1Q” → “Q1”; 第8行,“1Q” → “Q1”; 第9行,“HD” → “HR”; 图6-18中,看左边:电动机“ M1 ” →“ 1M ”、 “ M2 ” →“ 2M ”(区别PLC内存地址); 图6-18中,看右边:第2排的“S B1 ” →“S A ”、第4排的 “S K ” →“S A ”; 倒数13行,“LD” → “Q1”; 倒数6行,“ M1 ” →“ 1M ”; 倒数5行, “ M2 ” →“ 2M ”; P220:第14行,“1Q” → “Q1”; 第15行,“ 1D ” →“ 1M(电动机) ”; 第16行, “2Q” → “Q2”、“ 2D ” →“ 2M(电动机) ” ; 第17行,“不论 1D (Y000)还是 2D (Y001)” →“ “不论 电动机1M (Y000)还是 2M (Y001)”; 第20行,“1 D 或2 D ” → “1 M 或2 M ”; 第22行,“不论是 1D (Y000)还是 2D (Y001)” →“ “不论是 1M (Y000)还是 2M (Y001)”; P230:第5行,“41 OO MW” → “41 00 MW”; P230:第15行,kN/m 2 前面不要空格; 第19行,“负”字前面不要空格; P231:第18、19行,两处“m 3 ” → “m 3 ”(3是指数); 第19行,两处“1/2”是指数; P240:第4行,“开度r 1 、r 2 ” → “开度r 1 、r 2 ”(1和2是下标); 倒数9行,“ ” → “ ” ; 倒数8 行,“ ” →“ ”、“n、n-1、 l ”→“n、n-1、 i ”; 倒数7行,“与”字前面不要空格。
2838 次阅读|0 个评论
[转载]谈谈BSSE校正与Gaussian对它的处理
hjlyyc 2013-7-18 04:44
FROM http://hi.baidu.com/sobereva/item/e1c1874c603580adde2a9fdc 谈谈BSSE校正与Gaussian对它的处理 ​文/Sobereva First release 2009-Aug-27 Last update: 2012-Dec-7 计算A、B分子间的弱相互作用能时,一般不能简单地通过E_interaction = E_AB - E(A) - E(B)来计算,因为E_AB能量相对于E(A) + E(B)的降低来自两方面,一方面是真实的A、B分子间的相互作用能,这是我们要求的;另一方面来自于A、B分子的基函数在复合物体系中重叠,相当于增大了复合物的基组而使E(AB)能量降低(严格来说前提是所用的理论方法是基于变分原理的),这个部分贡献如果也掺入E_interaction,则高估了相互作用能(即实际上结合能没有算出来的那么负),所以要去掉,它称为Basis Set Superposition Error(BSSE)。所以双分子的相互作用能应该表述为E_interaction = E_AB - E(A) - E(B) + E_BSSE。对于弱相互作用,E_BSSE所占E_interaction的比例往往不小,甚至超过它,如果不进行校正,可能正负号都不对。 基组越小,单体间相互作用越弱则E_BSSE越大(也有人认为对于很小的基组,由于基函数不容易延伸到相邻原子上,E_BSSE反倒不大)。E_BSSE会随基组趋于完备而逐渐减小至0,给基组加上弥散函数能有效减小E_BSSE。对于氢键复合物,由于相互作用不算很弱,所以用了带弥散的中上等基组后,不做BSSE校正无妨。而pi-pi相互作用,即便用了aug-cc-pVTZ级别的基组,BSSE仍然不很小,若是用aug-cc-pVDZ的话十分建议加上BSSE校正。 计算E_BSSE有多种方法,Gaussian用的是目前使用最广泛的Boys和Bernardi发展的counterpoise correction方法,应注意这种方法计算出来的只是实际E_BSSE的近似,并非完全严格、精确。设E_i为第i个分子在自身基组下的能量,E_i'为第i个分子在全部n个分子上的基函数都出现下的能量,则计算n个分子相互作用能中的E_BSSE = ∑ ( E_i - E_i' ),E_BSSE必为正值。注意计算E_i与E_i'时的分子几何结构必须与i处在复合物时的一致,Gaussian会自动这样处理。 要计算A、B两个分子的相互作用能,在Gaussian中使用counterpoise=2关键字(可简写为counter=2),会计算5个体系,输出的能量按照如下顺序: E_AB:A、B基组下AB复合物的能量 E_A,bAB:A、B基组下A的能量 E_B,bAB:A、B基组下B的能量 E_A:A基组下A的能量 E_B:B基组下B的能量 最后会输出corrected energy (E_corrected)和BSSE energy (E_BSSE)。E_BSSE是校正能,E_BSSE = (E_A - E_A,bAB) + (E_B - E_B,bAB) ;E_corrected就是消除了因单体基组重叠造成的能量降低后的AB复合物能量,E_corrected = E_AB + E_BSSE。 BSSE校正后的真实的相互作用能这样计算:E_interaction = E_corrected - (E_A` + E_B`),其中E_A`和E_B`分别是A和B在孤立状态下经过优化后的能量。注意,相互作用能是复合物(考虑了BSSE时)减去单体在孤立状态时的能量差,单体在孤立状态的结构与在复合物中的结构并不相同,尤其单体结构呈柔性、分子间相互作用比较强时差异会较为明显。在实际中为了方便、省时往往做这样的近似:E_A`=E_A,E_B`=E_B,即是说计算单体能量时使用处在复合物状态下的结构,这个能量在计算E_BSSE时就已经顺便计算了,所以省得再单独算,此时相互作用能公式也有另一种等价写法E_interaction = E_AB - E_A,bAB - E_B,bAB,可见E_A和E_B并没被用到,其实可以不用在counterpoise任务中计算。 计算过程中会输出类似这样的语句Counterpoise: doing DCBS calculation for fragment 1。这里就是说明接下来计算的是E_A,bAB(假设A分子为fragment 1),其中DCBS代表dimer centered basis set,说明以A、B分子为中心的基函数都出现,但是计算中并不纳入B的电子和原子核,这称为计算A的能量时添加了B的ghost轨道;如果是doing MCBS calculation for fragment 1,就是要计算E_A,MCBS代表monomer centered basis set,计算中只出现属于A分子的基函数。 若计算n个分子间的BSSE,则关键字为counterpoise=n,结果输出顺序与计算相互作用能的方法与双分子的情况是一样的。能量按如下顺序输出:E_complex,E_1',E_2'...E_n',E_1,E_2...E_n。E_BSSE = E_1 - E_1' + E_2 - E_2' + ... + E_n - E_n'。E_corrected = E_complex + E_BSSE。相互作用能公式即E_interaction = E_corrected - ( E_1` + E_2` + ... + E_n` ),其中`的含义与前面所述一致,也可以将E_i`近似为E_i。计算过程中也用DCBS和MCBS来说明接下来将要计算的是哪项,但此时DCBS中的D的含义就不是具体指Dimer了,而是多分子复合物。 在分子内相互作用能计算时一般也要考虑BSSE,比如一条长链分子,计算一字形和字母C形的能量差就不能忽略这个问题。但由于两个片段属于同一个分子,需要特殊处理,而不能直接用上述方法。比如可以将分子人为地切成两段,悬键用比如H来封闭,然后适当调整两个片段,即让切断的部位离得远一些(否则这部分也会对E_BSSE产生贡献),可能造成BSSE的部位相对位置保持不变,然后将两个片段当成两个分子来同上获得E_BSSE。 如果有特殊原因,需要手动进行上述BSSE计算的每步操作而不直接用Gaussian提供的Counterpoise,可通过设定Ghost原子来实现,事实上这也是老版本不支持Counterpoise任务的Gaussian获得E_BSSE所必需的。只要把某个原子名后面加上-Bq就说明它是Ghost原子,如Na-Bq。Ghost原子也就是指这个原子照常有基函数,但是没有原子核和电子。注意Bq原子若破坏了原有对称性,需加上nosymm,否则可能中途报错停止。这样靠每次手工修改输入文件分别得到E_AB、E_A,bAB、E_B,bAB、E_A、E_B也就能算得E_BSSE,显然没直接用Counterpoise关键词省事。不过也有好处,就是手动做counterpoise的每个步骤比直接用counterpoise关键词更省时,尤其是高度对称的体系甚至能省几倍。因为用counterpoise关键词时所有任务都关掉了对称性,而手动做时,对单体、复合物仍然可以利用对称性来节省时间。 实际上无论计算何种分子,BSSE总是存在的,但常被忽略,尤其是分子内BSSE。计算柔性分子构象能量差人们往往还会注意这一点,但计算普通小分子时很多人会认为没有BSSE,实际上此时也存在着BSSE。如计算乙醇和它的同分异构体二甲基醚,都存在分子内BSSE,并且数值不一样。可以考察其中的氧原子,在乙醇中它连着一个C一个H,而二甲基醚中连着两个C,其它原子带来的BSSE显然是不一样的,这会影响到两种异构体的相对能量差,尽管微小到被忽略也无妨。哪怕对于乙醇和乙烷的甲基碳,虽然都是连着三个H和一个C,但是键长是有略微差别的,BSSE仍存在极微小的差异。再比如计算乙烷交错式与重叠式的能量差,显然两种构象BSSE也是不同的。有时这不再是个小问题,例如2006年报道的用MP2等后HF方法结合一些pople基组,优化出的苯环结构是略微弯曲的,显然不对,经过分子内BSSE校正后避免了这个问题,有兴趣可以去看JCP,128,144108和JCTC,5,2574。 分子内BSSE校正的方法不如分子间校正的方法广泛、成熟,虽然按上述方法切割成分子片段(比如把苯环切成一堆C-H)是可行的,但是如何切割过于任意,操作也麻烦,不易普适化。有人提出了atomic counterpoise(ACP)方法,将体系的总BSSE由每个原子的BSSE加和得到,即E_BSSE=∑ E_A,bA-E_A,bS,其中E_A,bA是A原子孤立存在时的能量,E_A,bS是A原子在周围Ghost原子基组下的能量。显然不可能计算每个原子的BSSE时都要带着全部其它原子的轨道,这样太慢,实际计算E_A,bS时可以只让A原子附近7埃以内的其它Ghost原子的轨道出现,再靠外的原子的轨道由于伸展不到A所以不用考虑,这样效率是比较高的。 Counterpoise方法有几个问题值得注意: 1 由于Ghost原子和真实原子是不同的,所以会破坏体系对称性。故计算高对称性二聚体时E_A,bAB和E_B,bAB的耗时都可能大于计算复合物的耗时。 2 Counterpoise有时有过校正问题,例如计算E_A,bAB时,A感受到的是完全“空闲”的B的基函数,能被A充分利用;而在复合物中A感受到的B的基函数已经有一定占据了,不能被A充分利用。由于两种状况B的基函数的状态不同,用E_A - E_A,bAB作为复合物中A的校正显得校正过头了,故所得相互作用能会被低估,有时结果还不如不用Counterpoise,对于最需要E_BSSE的小基组,过校正问题反倒容易出现。故Counterpoise该不该用其实没有明确答案,这和体系有一定关系。过校正容易在氢键体系出现,所以建议用个好点的基组就行了,而不用Counterpoise。对于范德华复合物,尤其是pi-pi作用,建议用aug-cc-pVDZ及以上基组+Counterpoise。 3 在Counterpoise下进行优化时由于无法获能量的得解析导数,只能通过数值方式进行,因此优化很慢。 附:高斯中counterpoise输入文件的写法 每个原子说明最后需要有一个整数说明这个原子属于第几个片段。例如: # MP2/6-31G Counterpoise=2 Opt Counterpoise with Cartesian 0,1 -整体、片段1、片段2的电荷和自旋多重度都一样,只需写一次即可。 1 0.00 0.00 0.92 1 9 0.17 0.00 2.73 2 1 0.77 0.00 3.43 2 9 0.00 0.00 0.00 1 使用Z矩阵时写法有些特殊,下面是test562的输入文件: #p uhf/3-21g counter=2 test scf=tight opt freq Gaussian Test Job 562: Counterpoise with different charges and multiplicites 0,1,0,3,1,2 -总电荷,总自旋;片段1的电荷和自旋多重度;片段2的电荷和自旋多重度。 O,0.0,0.0,0.0,1 -第一个原子必须用笛卡尔坐标。 O,1,ROO,2 -片段2的原子。 X,1,1.,2,X3O H,1,RO1H,3,HOX3,2,90.,0,1 H,1,RO1H,3,HOX3,2,-90.,0,1 -片段1的原子。-90.与1之间必须写0来代表这是二面角 X,2,1.,1,52.5,3,180.,0 -定位用的虚原子,不写所属片段号。由于是二面角,仍得写0。 H,2,RO2H1,6,H7OX,1,180.,0,2 H,2,RO2H2,6,H8OX,1,0.,0,2 ROO=2.98308 RO1H=0.94839 X3O=120.2827 HOX3=52.90868 RO2H1=0.94686 RO2H2=0.95173 H7OX=52.98178 H8OX=51.9632 注:如果是g09,也可以在定义坐标的时候写成诸如H(Fragment=1) -0.046866 0. 0.586860来让这个原子处于第一个片段。 补充 :由于有网友看过此文后对BSSE的认识仍然有误,因此我将对他的回复也贴到了这里(修改了符号以与本文一致),希望读者看了之后能对BSSE问题根源了解得更明白一些。为简化讨论,这里假设不考虑单体在复合物状态下结构的改变。 对于计算相互作用能,E_interaction = E_AB - E_A - E_B这个式子没错,这是由结合能的定义直接得到的。而必须指出的是,在实际计算中,E_interaction = E_AB,bAB - E_A,bA - E_B,bB这个式子中(b代表基组), E_AB,bAB作为复合物能量是对的,E_A,bA和E_B,bB作为单体能量也是对的,但是,这么计算相互作用能,在基组完备性差的时候,是绝对错误的!因为计算复合物和计算单体时不在相同的基组级别下。例如,A单体中靠近B单体的原子的基组完备性在复合物状态下高(因为B的基函数侵入此原子空间),而在A单体单独计算时低,这就导致能量在求差值时基组误差无法被抵消,这就是BSSE的根源。只有在基组完备的情况下,或者计算复合物和计算单体时都在相同基组级别下才能直接通过计算出的复合物能量减去单体能量来计算结合能。即正确的写法是E_interaction = E_AB,bQ - E_A,bQ - E_B,bQ,Q可以代表平面波基函数,也可以代表一套原子中心的基函数,比如可以指AB的基函数,这正是为什么本文中给出了这样的式子:E_interaction = E_AB - E_A,bAB - E_B,bAB,即计算复合物和单体的能量时都在AB基组下(文中E_AB和E_AB,bAB在符号上是等价的)。
个人分类: 分子模拟|5 次阅读|0 个评论
[转载]环境小卫星高光谱数据FLAASH精确大气校正方法
zoelu 2013-7-4 22:17
高光谱遥感数据的特点是光谱分辨率高、 波段连续性强,能够获得地物在一定范围内连续的、 精细的光谱曲线,具有非常好的应用前景。第一个星载民用成像光谱仪是高光谱成像仪 Hyperion ,其平均光谱分辨率为 10 nm ,空间分辨率达到 30 米。 2008 年 9 月 6 日我国 HJ-1 顺利升空,其中 A 星搭载了我国自主研制的空间调制型干涉高光谱成像仪( HSI )。 HSI 对地成像幅宽为 50 km, 星下点像元地面分辨率为 100 m , 115 个波段,工作谱段: 459~ 956nm 。具有 30 度侧视能力和星上定标功能。 ENVI 扩展工具 一、 HDF5 读取补丁 从 HSI 数据中获取元数据信息。下载地址: http://bbs.esrichina-bj.cn/ESRI/viewthread.php?tid=37118extra=page=1 二、 HJ-1 数据预处理补丁 为 HSI 数据添加中心波长信息和波段宽度( FWHM ) . 下载地址: http://bbs.esrichina-bj.cn/ESRI/thread-75575-1-1.html 或者 http://vdisk.weibo.com/s/lO16q 拷贝 sav 文件到 ENVI 安装目录的 save_add 目录下。 数据 准备 数据可网上免费下载,下载地址: 中国资源卫星数据服务网: http://www.cresda.com 环保部环境星下载服务网: www.secmep.cn (1) 选择 ENVI-File-Open External File-HJ-1- HJ-1A /1B Tools 工具(下图)。 (2) 选择 HSI ,单击 Input Files 选择 HSI 文件( .xml ) (3) 选择输出路径,单击 Apply 执行。 (4) 在 ENVI 主模块中,选择 Basic Tools-Convert Data(BSQ,BIL,BIP) ,将刚才生成的文件转成 BIP 储存顺序的文件。 按照 FLAASH 工具的要求,已经将 HSI 数据转成 BIP 储存顺序、带有中心波长信息、 波段宽度信息的 ENVI 格式文件。 说明:如果不用 HJ-1 数据预处理补丁,可以在波段列表中手动输入中心波长和波段宽度信息。 图 1 HJ-1A /1B Tools FLAASH 参数设定 HSI 数据已经经过了定标,单位是 100* W × m -2 × sr -1 × m m -1 ,启动 FLAASH 工具,在输入辐射率数据时候,缩放系数填写: 1000 ,(符合 FLAASH 对辐射亮度单位的要求) 图 2 辐射亮度值缩放系数 一、 传感器基本信息设置 选择传感器类型,成像中心点经纬度,成像时间,高度信息(成像区域和传感器飞行高度)等都可以在 HDF5 文件中获取。 选择 ENVI-File-Open External File-Generic Formats-HDF5 ,打开 .h5 文件,如下是相应参数填写。 图 3 成像中心点经纬度信息 图 4 成像时间(选择开始扫描时间,开始和结束时间相差不大) 传感器信息是固定的:飞行高度—— 650km ,像元分辨率—— 100 米 二、 大气参数设置 l 大气模型( Atmospheric Model )根据经纬度和日期来定。 40.3 度和 8 月份,选择 MLS 。 表 1 六种标准的大气模型 Latitude (°N) Jan. March May July Sept. Nov. 80 SAW SAW SAW MLW MLW SAW 70 SAW SAW MLW MLW MLW SAW 60 MLW MLW MLW SAS SAS MLW 50 MLW MLW SAS SAS SAS SAS 40 SAS SAS SAS MLS MLS SAS 30 MLS MLS MLS T T MLS 20 T T T T T T 10 T T T T T T 0 T T T T T T -10 T T T T T T -20 T T T MLS MLS T -30 MLS MLS MLS MLS MLS MLS -40 SAS SAS SAS SAS SAS SAS -50 SAS SAS SAS MLW MLW SAS -60 MLW MLW MLW MLW MLW MLW -70 MLW MLW MLW MLW MLW MLW -80 MLW MLW MLW SAW MLW MLW l 水气反演设置 Water Retrieval , 设置为 yes 。 使用水气去除模型,数据必须具有 15nm 以上波谱分辨率,且至少覆盖以下波谱范围之一: 1050-1210 nm ( 对应 1135 nm) 870-1020 nm ( 对应 940 nm) 770-870 nm ( 对应 820 nm) HSI 选择 820 nm 。 l 气溶胶模型 根据影像所在地区的情况来设定,这里选择 Rural (乡村)。 提供两种气溶胶去除方法: 2-Band ( K-T )方法 ( 类似模糊减少法 ) ,如果没有找到适应的黑值(一般是阴影区或者水体),系统将采用能见度值来计算;所以即使选择了该选项也要给能见度。 HSI 数据由于缺少短波红外波段,这里选择 None ,采用能见度值参与气溶胶去除。 能见度值大约参考值参见下表 天气条件 能见度 晴朗 40 to 100 km 中等雾、阴霾 20 to 30 km 厚雾、阴霾 15 km 或者更少 三、 高光谱数据参数设置 1) 自动选择通道定义 ( 推荐 ) 图 5 高光谱设置 浏览结果 选择校正前后相同像素点的 Z profile ,如下图为相同点植被的波谱曲线,可以看到大气校正之后植被的波谱曲线趋向正确。 HSI 数据的应用研究科查阅:马灵玲 , 王新鸿 , 唐伶俐的“ HJ -1A 高光谱数据高效大气校正及应用潜力初探”。 常见错位及解决方法 1、 ACC error:autochan IDL error: ACC_AUTOCHAN:Error:Upper K-T band not found automatically HSI 数据缺少 K-T 气溶胶反演的所需的波段,因此不能使用 K-T 方法反演气溶胶。 2、 ACC error:autochan ACC_AUTOCHAN:No water absorption channels selected 这种错误一直没模拟出来,一般是由于 水汽反演、大气模型的选择和飞行时间的问题,如飞行时间使用的是北京时间等。 3、 ACC error:avrd IDL error:REBIN:Result dimensions must be integer factor of original dimensions 这个错位一般出现在你多次执行 FLAASH 工程文件时,没有重新设置输出文件而选择替代的时候出现。好的方法就是重新启动 FLAASH 和加载一次工程文件。 4、 校正结果全部为 0. 可能是因为 HSI 定标结果不是 100* W × m -2 × sr -1 × m m -1 而是 10* W × m -2 × sr -1 × m m -1 ,源自: http://bbs.esrichina-bj.cn/ESRI/thread-102714-1-4.html 其他错误可参考: http://blog.sina.com.cn/s/blog_764b1e9d0100pvrk.html
个人分类: 大气校正|0 个评论
[转载]​ENVI中基于统计学模型的大气校正方法详解
zoelu 2013-7-3 16:09
基于统计学模型的反射率反演的方法主要有平场域法( Flat Field , FF )、对数残差法( Log Residuals )、内部平均法( Internal Average Relative Reflectance , IARR )、经验线性法 (Empirical Line) 。集中在 Basic Tools-Preprocessing- Calibration Utilities 菜单下。 1. 平场域定标 (Flat Field Calibration) Flat Field 定标工具通过选择图像中一块具有高反射率、光谱变化平坦的区域,利用这个区域的平均光谱值来模拟飞行时的大气条件下的太阳光谱。将每个像元的 DN 值除以选择区域的平均光谱值得到相对反射率,以此来消除大气的影响。 在使用这个工具前,需要利用 ENVI 提供的感兴趣区绘制工具( ROI Tool )在被定标图像上选择感兴趣区作为平场域( Flat Field ),感兴趣区可选择沙漠、大块水泥地、沙地等区域。操作过程如下: (1) 在主菜单中,选择 Basic Tools-Preprocessing- Calibration Utilities- Flat Field 。 在打开的 Calibration Input File 对话框中,选择输入文件,单击“ OK ” 。 (2) 在打开的 Flat Field Calibration Parameters 面板中 (图 13.4 ),在标有“ Select ROI for Calibration ”一栏中,选择感兴趣区(只能选择一个),作为平场域定标的平均波谱区。 (3) 选择输出路径及文件名,单击“ OK ”执行定标处理。 图 13.4 Flat Field Calibration Parameters 面板 2. 对数残差( Log Residuals ) 对数残差定标工具将数据除以波段几何均值,后再除以像元几何均值,可以消除光照、大气传输、仪器系统误差、地形影响和星体反照率对数据辐射的影响。定标结果的值在 1 附近。操作过程如下: (1) 在主菜单中,选择 Basic Tools - Preprocessing - Calibration Utilities - Log Residuals ,在 Log Residuals Calibration Input File 对话框中,选择输入文件,单击“ OK ”。 (2) 在打开的 Log Residuals Calibration Parameters 面板中, 选择输出路径及文件名,单击 OK 执行定标处理 。 3. 内部平均( IAR )反射率定标 IAR (Internal Average Relative) Reflectance 定标工具假定整幅图像的平均光谱基本代表了大气影响下的太阳光谱信息。把图像 DN 值与整幅图像的平均辐射光谱值相除,得到的结果为相对反射率。该工具特别适用于没有植被的干旱区域。 操作过程如下: (1) 在主菜单中,选择 Basic Tools-Preprocessing-Calibration Utilities- IAR Reflectance 。 在打开的 Calibration Input File 对话框中,选择输入文件,单击“ OK ”。 (2) 在打开的 IARR Calibration Parameters 面板中,选择输出到“ File ”或“ Memory ”。点击“ OK ”执行定标处理。 4. 经验线性定标 (Empirical Line Calibration) Empirical Line 定标 方法是假设图像 DN 值与反射率之间存在线性关系: 反射率 = 增益 * DN 值 + 偏移 利用两个已知点的地面反射光谱值,再计算图像上对应像元点的平均 DN 值,然后利用线性回归求出增益和偏移值,建立 DN 值与反射率之间的相互关系式,进行反射率的定标。消除了太阳辐亮度和大气程辐射。 ENVI 的 Empirical Line 定标工具要求至少需要一个 已知区域 的地面反射光谱值( Field Spectra )作为参照波谱,以及图像上对应像元点的波谱曲线( Data Spectra )。它们可以来自波谱剖面或波谱曲线、波谱库、感兴趣区、统计文件和 ASCII 文件。输入的波谱将自动被重采样,以与选择的数据波长相匹配。也可以用已经存在的系数对数据集进行定标。 1) 计算系数并定标( Compute Factors and Calibrate ) Empirical Line 定标工具 时,一般可以在图像上选择一个暗区和一个亮区作为已知区域(假定这些区域中的参照波谱是可以获得的)。使用越多的已知波谱也可以提高定标精度,至少需要一组已知区域的波谱。操作过程如下: (1) 选择 Basic Tools-Preprocessing- Calibration Utilities- Empirical Line- Compute Factors and Calibrate 。在打开的 Empirical Line Input File 对话框中,选择输入文件。单击“ OK ”, (2) 在打开的 Empirical Line Spectra 面板中(图 13.5 ),需要选择地面反射光谱值( Field Spectra )以及图像上对应像元点的波谱值( Data Spectra )。 图 13.5 Empirical Line Spectra 面板 l 选择数据(图像)波谱 在 Empirical Line Spectra 对话框中,单击“ Import Spectra ”按钮,打开 Data Spectral Collection 对话框(图 13.6 )。 在 Data Spectral Collection 对话框中,选择 Import-from ROI/EVF from input file ,选择定义好的感兴趣区文件,点击“ Apply ”,波谱名被输入到 Empirical Line Spectra 面板中。点击“ Cancel ”,关闭 Data Spectra Collection 对话框。 图 13.6 Data Spectral Collection 对话框 l 选择参照波谱 在 Empirical Line Spectra 对话框中,单击“ Import Spectra ”按钮,打开 Feild Spectra Collection 对话框,这个对话框与 Data Spectra Collection 对话框类似。 在 Feild Spectra Collection 对话框中,选择 Import-from ASD binary file ,选择对应图像波谱区域用 ASD 波谱仪测量波谱文件。点击“ Apply ”,波谱名被输入到 Empirical Line Spectra 面板中。点击“ Cancel ”,关闭 Feild Spectra Collection 对话框。 (3) 回到 Empirical Line Spectra 面板中,在顶部的列表内点击波谱名选择数据波谱。在底部的列表中,点击相应的参照波谱名,单击击“ Enter Pair ”按钮使两个波谱相关联,相关联的波谱将被列在“ Selected Pairs ”文本框中。(说明:如果错误的选择了关联波谱,在“ Selected Pairs ”文本框中单击关联波谱可移除) (4) 重复( 3 )( 4 )步骤可选择关联波谱。 (5) 单击击“ OK ”。打开 Empirical Line Calibration Parameters 对话框。 (6) 在 Empirical Line Calibration Parameters 对话框中,选择输出定标结果文件路径及文件名,将定标系数保存在 ASCII 文件中,在“ Output Calibration Filename ”文本框中键入第二个文件名,定标系数文件的默认扩展名是 .cff 。 (7) 单击 OK 执行定标过程。 2) 使用现有系数定标 (Calibrating Using Existing Factors) Calibrate Using Existing Factors 工具可以使用另一个定标过程中存储的纠正系数来运行经验行定标功能。 (1) 在主菜单中,选择 Basic Tools- Preprocessing- Calibration Utilities- Empirical Line-Calibrate Using Existing Factors 。选择输入文件,点击“ OK ”。 (2) 在打开的 Enter Calibration Factors Filename 对话框中,选择一个之前定标过程中创建的定标系数文件( .cff )。点击“ OK ”。 在打开的 Empirical Line Calibration Parameters 对话框中,选择输出路径及文件名。单击“ OK ”执行定标过程。 摘自《 ENVI 遥感图像处理方法》科学出版社
个人分类: 大气校正|1 次阅读|0 个评论
[转载]ENVI FLAASH大气校正常见错误及解决方法
zoelu 2013-7-3 16:07
本文汇总了 ENVI FLAASH 大气校正模块中常见的错误,并给出解决方法,分为两部分:运行错误和结果错误。前面是错误提示及说明,后面是错误解释及解决方法。 FLAASH 对输入数据类型有以下几个要求: 1、 波段范围:卫星图像: 400 - 2500nm ,航空图像: 860nm-1135nm 。如果要执行水汽反演,光谱分辨率 =15nm ,且至少包含以下波段范围中的一个: l 1050-1210 nm l 770-870 nm l 870-1020 nm 2、 像元值类型:经过定标后的辐射亮度(辐射率)数据,单位是:( μ W ) / ( cm2*nm*sr )。 3、 数据类型:浮点型( Floating Point )、 32 位无符号整型( Long Integer )、 16 位无符号和有符号整型( Integer 、 Unsigned Int) ,但是最终会在导入数据时通过 Scale Factor 转成浮点型的辐射亮度( μ W ) / ( cm2*nm*sr )。 4、 文件类型: ENVI 标准栅格格式文件, BIP 或者 BIL 储存结构。 5、 中心波长:数据头文件中(或者单独的一个文本文件)包含中心波长( wavelenth )值,如果是高光谱还必须有波段宽度( FWHM ),这两个参数都可以通过编辑头文件信息输入( Edit Header )。 运行错误 1. Unable to write to this file.File or directory is invalid or unavailable 。 img name="image_operate_55591299488609298" src="http://s15.sinaimg.cn/middle/764b1e9dt9de0c1db7f1e690" real_src="http://s15.sinaimg.cn/middle/764b1e9dt9de0c1db7f1e690" width="368" height="221" alt="ENVI FLAASH大气校正常见错误及解决方法" title="ENVI FLAASH大气校正常见错误及解决方法" action-data="http%3A%2F%2Fs15.sinaimg.cn%2Fmiddle%2F764b1e9dt9de0c1db7f1e%26690" action-type="show-slide" style="margin:0px;padding:0px;border:0px;list-style:none;text-align:center;display:block;" / 没有设置输出反射率文件名。 解决方法是单击 Output Reflectance File 按钮,选择反射率数据输出目录及文件名,或者直接手动输入。 2. ACC Error : convert7 IDL Error : End of input record encountered on file unit:0. img name="image_operate_12681299488623728" src="http://s13.sinaimg.cn/middle/764b1e9dt9de0c37da18c690" real_src="http://s13.sinaimg.cn/middle/764b1e9dt9de0c37da18c690" width="489" height="374" alt="ENVI FLAASH大气校正常见错误及解决方法" title="ENVI FLAASH大气校正常见错误及解决方法" action-data="http%3A%2F%2Fs13.sinaimg.cn%2Fmiddle%2F764b1e9dt9de0c37da18c%26690" action-type="show-slide" style="margin:0px;padding:0px;border:0px;list-style:none;text-align:center;display:block;" / 平均海拔高程太大。 注意:填写影像所在区域的平均海拔高程的单位是 km : Ground Elevation ( Km )。 3. ACC error : avrd: IDL error : Unable to allocate memory:to make array Not enough space ACC_AVRD img name="image_operate_45381299488650854" src="http://s7.sinaimg.cn/middle/764b1e9dt9de0c5508ff6690" real_src="http://s7.sinaimg.cn/middle/764b1e9dt9de0c5508ff6690" width="489" height="374" alt="ENVI FLAASH大气校正常见错误及解决方法" title="ENVI FLAASH大气校正常见错误及解决方法" action-data="http%3A%2F%2Fs7.sinaimg.cn%2Fmiddle%2F764b1e9dt9de0c5508ff6%26690" action-type="show-slide" style="margin:0px;padding:0px;border:0px;list-style:none;text-align:center;display:block;" / 为了能处理大数据, ENVI 采用分块计算的方式,这个提示是分块( Tile )太大了。 在高级设置里面( Advanced Settings ), tile size :它默认是 File- preferences -miscellaneous:cache 的大小,这个值正常是 1~4Mb ( 0 背景很少的情况下) ; 如果 0 背景较多,这个值还是需要设置大一些,比如 100~200m 。 4. ACC error:avrd:No nonblank pixels found IDL error: OPENR: Error opening file. Unit: 100, img name="image_operate_12611299488679380" src="http://s16.sinaimg.cn/middle/764b1e9dt9de0c66e8aef690" real_src="http://s16.sinaimg.cn/middle/764b1e9dt9de0c66e8aef690" width="489" height="391" alt="ENVI FLAASH大气校正常见错误及解决方法" title="ENVI FLAASH大气校正常见错误及解决方法" action-data="http%3A%2F%2Fs16.sinaimg.cn%2Fmiddle%2F764b1e9dt9de0c66e8aef%26690" action-type="show-slide" style="margin:0px;padding:0px;border:0px;list-style:none;text-align:center;display:block;" / 为了能处理大数据, ENVI 采用分块计算的方式,当 Tile 设置太小,而且有背景值( 0 ),就会出现一个 Tile 中全部为 0 的情况,提示这个错误信息。 在高级设置里面( Advanced Settings ), tile size :设置稍微大一些,如 100~200M 等。 5. ACC error:lsmooth2: IDL error: ACC_LSMOOTH2:Cannot continue with smoothing calculation img name="image_operate_14111299488717518" src="http://s10.sinaimg.cn/middle/764b1e9dt9de0c7949449690" real_src="http://s10.sinaimg.cn/middle/764b1e9dt9de0c7949449690" width="489" height="374" alt="ENVI FLAASH大气校正常见错误及解决方法" title="ENVI FLAASH大气校正常见错误及解决方法" action-data="http%3A%2F%2Fs10.sinaimg.cn%2Fmiddle%2F764b1e9dt9de0c7949449%26690" action-type="show-slide" style="margin:0px;padding:0px;border:0px;list-style:none;text-align:center;display:block;" / ENVI 的 FLAASH 提供领域纠正功能,但是 MODIS 、 AVHRR 等图像分辨率比较低,领域效应区分不出来。 解决方法是在在高级设置里面( Advanced Settings ),将领域纠正( Use Adjacency Correction )设置为 No 。 6. ACC error:modrd5:Nonfinite numbers in coefficient array coef?xml:namespace prefix = o ns = urn:schemas-microsoft-com:office:office / IDL error: CDRIVER4V3R2:Wait for MODTRAN4 calculation to finish… img name="image_operate_23261299488798487" src="http://s15.sinaimg.cn/middle/764b1e9dt9de0c941eb5e690" real_src="http://s15.sinaimg.cn/middle/764b1e9dt9de0c941eb5e690" alt="ENVI FLAASH大气校正常见错误及解决方法" title="ENVI FLAASH大气校正常见错误及解决方法" action-data="http%3A%2F%2Fs15.sinaimg.cn%2Fmiddle%2F764b1e9dt9de0c941eb5e%26690" action-type="show-slide" style="margin:0px;padding:0px;border:0px;list-style:none;text-align:center;display:block;" / 提示传入 MODTRAN 模型参数有误,常常是由于太阳高度角太小或者太大引起的。提示这个错误之前会出现以下提示框。 img name="image_operate_90521299488812220" src="http://s1.sinaimg.cn/middle/764b1e9dt9de0ced22850690" real_src="http://s1.sinaimg.cn/middle/764b1e9dt9de0ced22850690" width="478" height="238" alt="ENVI FLAASH大气校正常见错误及解决方法" title="ENVI FLAASH大气校正常见错误及解决方法" action-data="http%3A%2F%2Fs1.sinaimg.cn%2Fmiddle%2F764b1e9dt9de0ced22850%26690" action-type="show-slide" style="margin:0px;padding:0px;border:0px;list-style:none;text-align:center;display:block;" / 解决方法是确认填写的影像中心位置经纬度信息(西经为负数、南纬为负数)、影像成像时间(格林威治时间)是否正确。 结果错误 1. 结果中某一个波段或者多个波段全部为 0 或者负值 主要产生的原因是输入的辐射亮度数据值偏小。 可能有以下几种情况: (1) 在传感器定标的时候选择的是表观反射率( Reflectance )而不是辐射率数据( Radiance ); (2) 没有做传感器定标,即没有将 DN 值转换为辐射率数据; (3) 选择了错误的波谱响应函数; (4) 用 BandMath 做了辐射亮度的单位换算,在 FLAASH 中导入辐射亮度数据时, Scale Factor 选择的不是 1.0000 。 2. 结果中极大值、极小值非常多,也就是 0~10000 之外的值。当选择 RGB 假彩色显示的时候,出现花花绿绿的情况。 img name="image_operate_27081299488837934" src="http://s12.sinaimg.cn/middle/764b1e9dt9de0d077266b690" real_src="http://s12.sinaimg.cn/middle/764b1e9dt9de0d077266b690" alt="ENVI FLAASH大气校正常见错误及解决方法" title="ENVI FLAASH大气校正常见错误及解决方法" style="margin:0px;padding:0px;border:0px;list-style:none;text-align:center;display:block;" / 辐射定标得到的辐射率数据单位与 FLAASH 要求的单位不一样。 可能有以下几种情况: (1) 用 BandMath 做了辐射亮度的单位换算,在 FLAASH 中导入辐射亮度数据时, Scale Factor 选择的不是 1.0000 。 (2) 没有做单位换算。 3. 结果中部分像元为负值 这个属于正常现象。 FLAASH 是采用 MODTRAN 辐射传输模型模拟成像中的大气过程,而且很多大气属性都是通过图像来估算,加上大气组成的非均一性,即使 MODTRAN4 模型精度很高,也不能完全表达大气辐射传输的真实状况。当影像上有强吸收或者高反射地物时候,就会出现部分像元为负值(如深水、高密度)或者大于 10000 。 解决方法可以手动修改,如用周围的像元的平均值代替,可以使用 ENVI 下的 The DEM Editing Tool 工具。 4. 结果图像以 RGB 显示比原图像视觉效果要差,如模糊。 这个需要了解下遥感软件 RGB 显示机制。遥感软件为了让遥感图像显示更加“艳丽”,方便解译,默认会对图像进行拉伸显示,一般是 2% 的线性拉伸。也就是我们常常会看到图像值有两个,一个是 Scrn 值,也就是拉伸之后的显示值,一个是 Data 值,也就是原始的图像 DN 值。另外一个方面,目前我们的 RGB 加色法显示都是基于 8bit 显示,也就是 0~255 。 FLAASH 大气校正之后的结果是 16bit 的整型,而且存在一些极小、极大值,这些对直方图整体形状有一定的影响,影响拉伸效果。解决方法是利用 ENVI 下的 Interactive stretching 工具,选择有效值范围进行拉伸。 还有一种情况是校正图像有很多背景值,比如经过几何校正的整景 TM 影像。背景及图像边缘处的像元在大气校正之后变成负值或者 0 值,由于这样的像元数量多,对图像的整体拉伸影响很大。解决方法就是将这些像元掩膜掉,如在 ROI Tool 中,利用 Option-Band Threshold to ROI 建立一个 ROI 进行图像裁剪。 5. 结果图像为什么像元值大多是大于 1 ENVI FLAASH 考虑到数据储存和后续处理,将大气校正得到的反射率结果乘以 10000 变成 16bit 整型。 如果想让反射率结果在 0~1 范围,可用 BandMath ,表达式为 b1/10000.0 。
个人分类: 大气校正|1 次阅读|0 个评论
ArcGIS的地理校正
JiuBaiYi 2012-10-30 20:47
别忘了先给要空间校正的栅格影像定义空间坐标系统。
个人分类: 科研|2791 次阅读|0 个评论
杜月笙正传
duke01361 2012-10-23 23:54
杜月笙正传 徐铸成 著 浙江人民出版社 1982年6月 第1版,1983年5月第2 次印刷 我只想为作者校正一个错误,本来这本书只是第2次印刷,但作者却给出了“再版后记”这似乎值得商榷 杜月笙简介 杜镛 表字 月笙, 籍贯:上海 职务:上海市参议会议员 上海地方协会会长 中华民国红十字会总会副会长 上海慈善团体联合会董事长 上海市社会救济事业协会理事长 中华民国机制纺织工业同业联合会董事长 上海市商会监事 上海市银行商业同业公会董事长 中国银行董事 复旦中学校董 《申报》董事长 以下略
个人分类: My Ideas|3134 次阅读|0 个评论
天地图截图之空间校正(傻瓜式)
热度 1 majiaping 2012-9-8 17:34
天地图截图之空间校正(傻瓜式)
注:所用软件有arcgis10.0、envi4.7和地图一把抓
6103 次阅读|3 个评论
多重假设检验中的p值校正
bioseq 2012-9-4 09:48
多重假设检验中的p值校正
在生物学特别是基因组学的研究工作中,经常会遇到多重假设检验(multiple testing)的问题;此时,得到的原始p值需要进行校正后才能使用,那么哪种校正方法更加适合自己的研究工作呢?p-values, false discovery rates(FDR) 和 q-values有什么不同?它们分别代表什么意义?对于统计科班的同学来说,这不过是小菜一碟;但对于纯生物出身的同学来说,别说去看公式了,光是听听就觉得头大!不过幸运的是,有牛人(William S Noble)了解我们的苦衷,于是一篇nature biotechnology的文章诞生了——《How does multiple testing correction work?》。这片文章不长,只有3页,用不了多长时间就可以看完。更加令人高兴的是,全篇没有一个让人头大的公式;了解基本的统计学知识、特别是p值的相关概念之后,阅读这片文章就不会有太大的困难了。作者以一个生物学例子贯穿全篇,这个例子对于大多数生物专业的同学来说都非常容易理解——在人的21号染色体上寻找CTCF(一个高度保守的锌指DNA结合蛋白)的潜在结合位点。作者先介绍了零假设(null hypothesis),进而引出了p-value的概念。之后,解释了为什么原始p值不能够直接使用,从而过渡到p值校正的话题。在这一部分,作者层层深入,以简洁明了的语言介绍、解释了Bonferroni adjustment、false discovery rate (FDR)、q-value和local FDR的概念、由来、意义等基本但非常重要的知识。最后作者给出了实际应用时的指导建议,并以点睛之笔概括总结了全文中的要点。如果你的工作涉及p值的校正、FDR、q值等概念,这篇文章绝对胜任引你入门的角色(但绝不仅限于此!)。 文章链接: http://www.seq.cn/forum.php?mod=viewthreadtid=3504 1 2 3 When prioritizing hits from a high-throughput experiment, it is important to correct for random events that falsely appear significant. How is this done and what methods should be used? Imagine that you have just invested a substantial amount of time and money in a shotgun proteomics experiment designed to identify proteins involved in a particular biological process. The experiment successfully identifies most of the proteins that you already know to be involved in the process and implicates a few more.
12937 次阅读|0 个评论
方案
yunxiangzhao 2012-8-31 16:02
中国不同群体杜洛克遗传参数评估研究及应用(提高生猪饲养效率) 目的:结合实用性和科学性,研究不同的遗传参数评估方法的准确性以及在不同群体不同性状的适用性。探索跨场或者跨区域遗传评估的有效方法。获得不同群体的遗传评估的遗传参数。确定不同群体的遗传评估的BLUP模型。获得不同群体的遗传评估用性状的校正公式。 具体思路: 1, 资源选择:至少5个群体、不同地域分布、不同品系,每个群体基础群在300头以上。三年以上的测定性能记录数据。优先选择有自动测定系统数据和B超测定数据的群体。 2, 基础资料的整理,挖掘当前中国杜洛克选育的水平和现状;分析引种以来种群遗传结构变化,评估育种措施、方案和效果。(国外引种效率比较以及影响因素) 不同育种组织的育种效果(联合育种,公司化育种,输入式育种) 3, 比较不同遗传评模型的相对效率,确定猪遗传评估的最佳模型。 4, 对初生重、成活率 、料肉比 、肌内脂肪含量等性状进行遗传参数估计。 母猪初情期的遗传参数估计,发情间隔遗传分析, 公猪不同日龄和体重终测结果的效果比较,公猪睾丸遗传参数,与精液品质关系 5,  对来自不同群体、不同世代的个体进行育种值估计和比较, 种猪遗传评估质量控制、、遗传联系度量与建立 6, 遗传和选择进展 7, 实现遗传力 realized heritability for shell length in the Manila clam Ruditapes philippinarum
1 次阅读|0 个评论
[转载]辐射校正
gosci 2012-8-3 14:55
可以先看本链接,再看下面文字 http://glaciershine.blog.163.com/blog/static/2963259201010301081697/ 一、辐射定标 1. 由于ENVI 4.4 中有专门进行辐射定标的模块,因此实际的操作十分简单。将原始TM 影像打开以后,选择 Basic Tools–Preprocessing–Calibration Utilities–Landsat TM 2. 进入下一步参数选择:根据传感器类型选择Landsat 4,5 或者7。从遥感影像的头文件中获取Data Acquisition 的时间,Sun elevation。如果你是用File–Open External File–Landsat–Fast 的方法打开header.dat(头文件) 的话,sun elevation 就已经填好了。这里Calibration Type 注意选择为Radiance。输出文件,定标就完成了。 二、大气校正 简单一点的大气校正可以采用ENVI的FLAASH模块,以下就是FLAASH操作的步骤: 1. FLAASH 模块的进入方法是Spectral–FLAASH,或者是Basic Tools–Preprocessing–Calibration Utilities–FLAASH。 2. FLAASH 模块的操作界面分为三块:最上部设定输入输出文件;中间设定传感器的参数;下部设定大气参数。 3. 首先设定输入输出文件。FLAASH 模块要求输入辐亮度图像,输出反射率图像。之前我们进行了辐射定标,得到辐亮度图像,在这里要把BSQ 格式的图像转换为BIL 或者BIP 格式的图像,然后再Input Radiance Image 中选择转换格式后的图像。(Basic Tools–Convert Data(BSQ,BIL,BIP))。这里注意,当输入图像后,程序会让你选择Scale Factor,即原始辐亮度单位与ENVI 默认辐亮度单位之间的比例。ENVI 默认的辐亮度单位是μW/cm2 •sr•nm,而之前我们做辐射定标时单位是W/m2 •sr•μm,二者之间转换的比例是10,因此在下图中选择Single scale factor,填写10.000。 4. 此外,如果TM 影像的头文件中没有波段的信息,在这里也要求你提供一个.txt 文件以包含此信息。那么,准备好一个.txt 文件,其中含有一列TM 每个波段中心波长的信息。 5. 在Output Reflectance File 和Output Directory for FLAASH files 里面设定输出文件的文件名和位置。 6. 设定传感器参数。首先是Scene Center Location,即遥感图像中心的坐标,以及Flight Date, Flight Time GMT,这三者都可以在TM 的头文件中找到,填入即可。 7. 在Sensor Type 菜单中选择Landsat TM5。此时Sensor altitude 自动填上为705km。而Pixel Size 填为30m。 8. 根据遥感影像研究区实际情况,填写Ground Elevation,比如华北平原可以写为0.05km。 9. 最关键的为大气参数部分: a) Atmospheric Model( 大气模式): 共有Sub-Arctic Winter (SAW) ,Mid-Latitude Winter (MLW),U.S. Standard (US) ,Sub-Arctic Summer(SAS), Mid-Latitude Summer (MLS) 和Tropical (T) 。根据经纬度和时间可以选定研究区的大气模式,见ENVI Help。 b) Aerosol Model(气溶胶模式):有Rural, Urban, Maritime 和Tropospheric四种选择。根据实际情况选择即可。关于此四种模式的解释见ENVI Help。 c) 当我们选择TM 时,可选的参数还有Aerosol Retrieval 和Initial Visibility。这两个参数对最后的结果又相当重要的影像,因此最好能调查到当地的Initial Visibility。此外,AERONET 在全世界各地有测定AOD(Atmospheric Optical Depth)的站点,可以查询AOD 以后转换为消光系数,通过消光系数估算能见度,此步骤比较繁琐,在此不予详述。如果采用Aerosol Retrieval 中的K-T算法计算Visibility,且能够计算出结果的话,则采用K-T 算法的能见度,否 则采用Initial Visibility 所指定的能见度。 d) 关于Aerosol Retrieval。如果选择了下拉菜单中的K-T method,那么需要在Multispectral Settings 中设定参数,在Assign Default Values Based on Retrieval Conditions 中选择Over-land Retrieval Standard (660:2100nm)即可。根据不同的研究区可以设定不同的模式。其他设定可以不改变。Apply即可。 FLAASH模块的大气校正 2008-06-24 22:21 FLAASH模块的大气校正 1.1 FLAASH模块简介 FLAASH是由世界一流的光学成像研究所-波谱科学研究所(Spectral Sciences)在美国空气动力实验室支持下开发的大气校正模块。波谱科学研究所在1989年大气辐射传输模型开发初期就广泛从事MODTRAN的研究工作,已成为大气辐射传输模型开发过程中不可缺少的一员。FLAASH适用于高光谱遥感数据(如HyMap,AVIRIS,HYIDCE,HYPERION,Probe-1,CASI和AISA)和多光谱遥感数据(如陆地资源卫星(landset),SPOT,IRS和ASTER)的大气校正。当遥感数据中包含合适的波段时,用FLAASH还可以反演水气、气溶胶等参数。 ENVI中大气校正模型FLAASH,是高光谱辐射能量影像反射率反演的首选大气校正模型。FLAASH能够精确补偿大气影响,其适用的波长范围包括可见光至近红外及短波红外,最大波长范围为3μm。其他的大气校正模型是计算方法基于查找表(Look-up Table)、利用插值方法计算,而FLAASH是直接移植了modtran4中的辐射传输计算方法。用户可以选取代表研究区的大气模型和气溶胶类型,并且对每景影像,Modtran都有独特的解决方案。 1.2 ASTER数据预处理 ASTER L1B数据是记录是DN(Digital Number)值,而基于FLAASH大气校正过程中,需要的是辐射能量值。因此,需要对ASTER L1B数据辐射定标,即把无量纲的DN值转换成有量纲的分辐辐射亮度值的过程(式1), Radiance=gain*DN+offset             (式1) 其中,gain是增益,offset是偏差。经辐射定标后,得到天顶辐射能量值,其量纲为W/(m2.sr.um)。 ASTER数据多以HDF格式储存,利用ENVI软件中Baisc Tools-Preprocessing-Data-Specific Utilities-View HDF Global Attribute功能,读取相应ASTER HDF文件中的增益、偏差、成像时间和中心点坐标信息,增益和偏差信息见表1。 根据ASTER数据特点[6],需要对短波红外波段数据重采样,以使其与可见光-近红外波段影像的像元数相同。然后,把可见光-近红外波段与短波红外段按波段顺序合并成新的文件。利用ENVI软件,编辑新文件的头文件。主要输入波长值见文献6、波长单位、增益、偏差和中心波长半极值宽度(FWHM,Full Width Half Maximum)。中心波长半极值宽度可以近似用每个波段的带宽代替,如果波长单位为微米,应该转成纳米单位。然后,利用ENVI软件中Baisc Tools-Preprocessing-General-Purpose Utilities-Apply Gain and Offset功能下,生成辐射能量数据。由于ENVI软件保存的数据格式通常为波段顺序格式(BSQ),而FLAASH大气校正模块使用的是波段逐行交叉顺序(BIL)或波段逐像元交叉顺序(BIP)格式的文件,因此,ASTER数据预处理的最后一步是把BSQ格式的文件转换成BIL或BIP格式的文件. FLAASH模块参数设置 1.3.1 尺度转换因子的计算 FLAASH模块中,在输入辐射能量数据时,同时要求输入尺度转换因子。尺度因子有两种输入方式①当各波段尺度转换因子不同时,选择每一种的输入方式,即事先把尺度因子输入到记事本文件中,然后,从记事本文件中直接读取;②当尺度转换因子相同时,选择第二种输入方式图2。由于模块中要求辐射能量的量纲是μW/(cm2·nm·sr),而经辐射定标ASTER数据的量纲为W/(m2·um·sr),所以后者还需通过换算关系式1μW/(cm2·nm·sr)=10 W/(m2·um·sr)进行量纲转换。因此,利用FLAASH模块校正ASTER数据时,其尺度转换因子为10。对其它类型遥感数据大气校正时,可以参照上述方法计算相应的尺度转换因子。判断尺度因子设置正确与否可的方法是依据图像的数据统计特征,即当统计数据没有负值和大于1×放大系数的数值,则可以认为尺度转换因子设置正确。 1.3.2 FLAASH其它参数的设置 (1)图像中心点坐标 可以从相应的HDF文件中找到,也可以从屏幕上直接读取影像的中心坐标,对反演结果影响不大。当影像位于西半球时,经度为负值; (2)传感器类型 当选择传感器类型时,模块会选择相应的类型的传感器波段响应函数,同时系统一般会自动设置传感器的高度和图像的空间分辨率; (3)海拔高度 海拔高度为研究区的平均海拔; (4)数据获取日期和卫星过境时间 卫星过境时间为格林尼治时间,可以从相应的HDF文件中找到; (5)大气模型 模块提供热带、中纬度夏季、中纬度冬季、极地夏季、极地冬季和美国标准大气模型,研究者根据数据获取时间选择相应的大气模型; (6)水气反演 大多数多光谱数据不推荐反演水汽含量; (7)气溶胶模型 可供选择的气溶胶模型有无气溶胶、城市气溶胶、乡村气溶胶、海洋气溶和对流层气溶胶模型。当能见度大于40Km时,气溶胶类型选择对反演没有太多影响,一般情况下利用ASTER数据不做气胶反演; 在高级设置中,①Modtran 分辨率(Modtran resolution):一般设置成5cm-1;②反射率输出的时尺度系数,默认尺度系数是10000,可以使用默认的尺度系数。若使用默认的尺度系数,大气校正后得到反射率图像的数值域为:0-10000。其余参数使用默认值。 出处华夏土地网: http://bbs.hxland.com/forum.php?mod=viewthreadtid=7403722 鲍平勇,河海大学,硕士论文p18
个人分类: ENVI|1002 次阅读|0 个评论
[转载]地形校正相关
sandy513 2012-6-20 18:54
基于ENVI的面向坡向的地形校正方法实现 地形校正的目的   地形校正的目的主要是补偿由于不规则的地形起伏而造成的地物亮度的变化。由于这种变化会导致相似或同种植被的反射率不一致而影响遥感影像的分类精度,因此,精确的地形校正不仅能提高影像分类的精度,而且还是遥感应用的前提。 实验方法: 1、数据准备 试验中使用的影像数据是云南省的Landsat5 TM 影像, 影像获取时间为2006年01月25日, 影像中心位置位于东经99.59E,北纬28.36N,太阳高度角和方位角分别为36.15和146.46。文中选取了一块2014 pixels ×934 pixels的试验区域,该地区的地物类型主要是山林。 2、数据预处理 2.1 数据定标处 图2.1Landsat TM元数据以及DEM、NUM数据列表 如上图所示,landsat TM数据总共包含了7个波段的数据,其中6波段为热红外波段,其他6个波段为可见光波段。在软件中逐一打开7个波段,再进行波段组合给用户带来了很大的不便。利用envi主菜单 open image file L5131041_04120060125_MTL.txt打开文件,该文件数据可自动进行波段组合以及分类,显示结果如下图2.2: 图2.2 波段自动组合分类 波段数据读入以后进行数据的辐射定标处理,在envi下有一个专门针对Landsat数据的定标处理模块,具体的操作步骤如下: ENVI主菜单 Basic Tools Preprocessing Calibration Utilities Landsat Calibration。 图2.3 Landsat定标模块自动读取参数 2.2 数据裁剪 数据裁剪的方式有很多,在本次试验中主要采用数据的规则采用——ROI rectangle裁剪。打开影像数据,在其image窗口中点击:Overlay Region of Interest,弹出ROI TOOL对话框: 图2.4 ROI裁剪对话框 感兴趣区划好之后进行影像数据的裁剪:ENVI主菜单 Basic Tools Subset Data via ROIs。 对于DEM数据而言,该数据的投影方式以及空间分辨率可能会与上述影像数据有差,为了能将二者的数据进行很好的吻合,在进行ROI裁剪的时候需要对其ROI进行一个Map的转换,具体转换方式如下: 图2.5 ROI在不同影像间的地图转换 在此之后选择DEM数据,并对其进行裁剪,得到本实验所需要数据。 2.3 坡度、坡向数据获取 在进行地形校正的过程中,坡度、坡向对于其反应地表真实的影响很大。在获得DEM数据以后,通过ENVI主菜单下 Topographic Topographic Modeling进行坡度、坡向数据的提取。 图2.6 坡度、坡向提取 在地形校正过程中影响最大的坡向为南坡,需要通过BandMath工具进行南坡数据(157.5,202.5)的提取以及对应的坡度数据提取,为后面的地形校正提供数据。具体的实现方法、步骤如下: ENVI主菜单 Basic Tools Band Math (b1 ge 157.5 and b1 le 202.5)*b1+(b1 lt 157.5 and b1 gt 202.5)*0,Add to List ,点击OK,配置b1对应的波段数据Aspect文件,即可完成南坡信息提取。 图2.7 南坡信息提取条件 关于南坡影像信息对应的坡度影像,TM影像区域信息的提取方法依然采用Band Math进行提取,提取公式如上图2.7中的(b1 ne 0)*b2+(b1 eq 0)*0,其中b1代表南坡影像信息,b2代表TM影像以及坡度影像信息。重复操作两次即可提取相对于的南坡区域的TM影像以及坡度影像。 3、地形校正——坡度匹配技术 C校正模型的基本思想是:对于任意波段影像的像素DN值和其对应的太阳入射角余弦值都遵循线性关系。理想情况下,当太阳入射角为零或小于零时,表明该点缺乏太阳光照,则该点的DN值应该为零,该拟合直线应通过原点。然而,实际情况是,由于大气散射和地表相邻点反射光折射的缘故,使像素DN值和太阳入射角α成一定的余弦关系。本文采用的模型是在二阶校正模型的基础上改进建立的,需经过二个阶段(二次校正)才能得到校正结果。 按照公式(1)计算太阳有效入射角余弦值 cos i= cosz*cosS+sinz*sinS*cos(ψx-ψn) (1) 式中, z 为太阳天顶角(即,90-太阳高度角), ψx 为太阳方位角, S 为坡度, ψn 为坡向。 在ENVI主菜单 Basic Tools Band Math 中输入: cos((90-36.1530740)*0.01745)*cos(b2*0.01745)+sin((90-36.1530740) *0.01745)*sin(b2*0.01745)*cos((146.4600750-b4) *0.01745) 其中,b2 代表坡度;b4代表坡向。 第一阶段(第一次校正),按公式(2)进行。 L H = L T + ( L Tmax -L Tmin ) * (( -X ) / ) (2) 式中,LH为第一次校正结果,LT为原始影像,LTmax为原始影像最大值,LT民为原始影像最小值, 为拉伸为0-255的阳坡(南坡)太阳有效入射角余弦值的平均值,X为拉伸为0-255的太阳有效入射角余弦值。 在ENVI主菜单 Basic Tools Band Math 中输入: b1+(max(b1)-min(b1))*((mean(b2)-b3)/mean(b2)) 其中,b1 代表Landsat 5TM原始数据;b2代表拉伸为0-255的阳坡(南坡)太阳有效入射角余弦值的平均值;b3代表拉伸为0-255的太阳有效入射角余弦值。 第二阶段(第二次校正),根据公式(3)求算模型修正系数。 Cλ= ( S′λ-Nλ ) / ( N′λ-Nλ ) (3) 式中:Cλ为图像各波段的模型修正系数,S′λ为第一次校正结果阳坡(南坡)平均值,N′λ为第一次校正结果阴坡(北坡)平均值,Nλ为原始图像阴坡(北坡)平均值。 在ENVI主菜单 Basic Tools Band Math 中输入: (mean(b1)-mean(b2))/(mean(b3)-mean(b2)) 其中,b1 代表第一次校正结果阳坡(南坡)平均值;b2代表原始图像阴坡(北坡)平均值;b3代表第一次校正结果阴坡(北坡)平均值。 把计算出的各波段模型修正系数带入公式(2)建立公式(4)进行第二次校正。 L H = L T + ( L Tmax -L Tmin ) * (( -X ) / ) * Cλ (4) 在ENVI主菜单 Basic Tools Band Math 中输入: b1+(max(b1)-min(b1))*((225.43-b3)/225.43)* (mean(b4)-mean(b5))/(mean(b6)-mean(b5)) 其中,其中,b1 代表Landsat 5TM原始数据;b2代表拉伸为0-255的阳坡(南坡)太阳有效入射角余弦值的平均值;b3代表拉伸为0-255的太阳有效入射角余弦值。b4 代表第一次校正结果阳坡(南坡)平均值;b5代表原始图像阴坡(北坡)平均值,b6第一次校正结果阴坡(北坡)平均值。 4、结果展示 图 4.1 原始影像、第一次校正结果以及坡向校正结果对比图 图 4.2 原始影像、第一次校正结果以及坡向校正结果对比图 根据上述两幅图的对比可以看出,在经过坡向地形校正以后其校正影像显示结果很显著。需要注意的是:利用该种方法进行地形校正的时候一定要保证dem数据的准确性,不然在后面进行坡度坡向信息提取的时候容易出现与实际情况不符的状况,最终影像到纠正结果。 结果分析 采用该种方法很好的进行了地形校正,是不是最终的结果能将山体阴影下的植被给予很好的显示,需要对其校正前后获取的NDVI值进行对比分析,从而确定该种方法的精度。 图5.1 校正前后的NDVI值对比 图5.2 山体阴影处经过校正后的NDVI值的变化 上述两幅图是从山体阴影上分析校正前后NDVI值的变化情况,从这个DN值的变化上可以看出,经过校正后的影像获取的NDVI值更为接近实际地表的植被覆盖情况。对于裸露的地表,其NDVI值的变化是不是很大,下面从裸土着手分析,坡向地形校正对于不存在阴影的地区的校正结果如何。 图 5.3裸露地表NDVI值的变化 从上图5.3的结果对比分析可以看出,裸露地表的NDVI值的变化在0.03范围,变化不是很大,这与开始设想的坡向地形校正方法主要是针对山体阴影区域的校正,对于平坦区域的作用影像还是比较小的。 综合上述的结果对比分析以及影像显示效果分析得出,该方法能很好的帮助我们进行山体阴影区域的植被生长情况的恢复,为以后的研究工作提供了更为便捷的校正方法。 坡向工作原理制定 用于识别出从每个像元到其相邻像元方向上值的变化率最大的下坡方向。坡向可以被视为坡度方向。输出栅格中各像元的值可指示出各像元位置处表面的朝向的罗盘方向。将按照顺时针方向进行测量,角度范围介于 0(正北)到 360(仍是正北)之间,即完整的圆。不具有下坡方向的平坦区域将赋值为 -1。 坡向数据集中每个像元的值都可指示出该像元的坡度朝向。 从概念上讲, 坡向 工具将根据要处理的像元或中心像元周围一个 3 x 3 的像元邻域的 z 值拟合出一个平面。该平面的朝向就是要处理的像元的坡向。 下图显示的是输入高程数据集和输出坡向栅格。 1. 坡向算法 移动的 3 x 3 窗口会访问输入栅格中的每个像元,而每次位于窗口中心的像元的坡向值将通过一种将纳入八个相邻像元值的算法进行计算。这些像元使用字母 a 至 i 进行标识,其中 e 表示当前正在计算坡向的像元。 像元 e 在 x 方向上的变化率将通过以下算法进行计算: = (( c + 2 f + i ) - ( a + 2 d + g )) / 8 像元 e 在 y 方向上的变化率将通过以下算法进行计算: = (( g + 2 h + i ) - ( a + 2 b + c )) / 8 代入像元 e 在 x 方向和 y 方向上的变化率,坡向将通过以下算法进行计算: a spect = 57.29578 * atan2 ( , - ) 然后,坡向值将根据以下规则转换为罗盘方向值(0 到 360 度): 2. 坡向计算示例 示例中,将计算移动窗口内中心像元的坡向值。 中心像元 e 在 x 方向上的变化率为: =(( c + 2 f + i )-( a + 2 d + g ))/8=((85 + 170 + 84))-(101 + 202 + 101))/8=-8.125 像元 e 在 y 方向上的变化率为: = (( g + 2 h + i )-( a + 2 b + c ) / 8=((101 + 182 + 84)-(101 + 184 + 85))/ 8 = -0.375 坡向计算如下: a spect = 57.29578 * atan2 ( , - )= 57.29578 * atan2 (-0.375, 8.125) = -2.64 由于计算得出的值小于零,则根据最终规则得出: c ell = 90.0 - aspect = 90 - (-2.64)= 90 + 2.64 = 92.64 中心像元 e 的值 92.64 表明它的坡向为朝东。 3. 参考文献 Burrough, P. A., and McDonell, R. A., 1998. Principles of Geographical Information Systems (Oxford University Press, New York), 190 pp.
个人分类: 地形校正相关|1 次阅读|0 个评论
[转载]几何校正与几何配准的区别与联系
lishumin 2012-3-8 22:58
我在其他论坛上看到 “ 其实几何校正和影像配准原理是一样的,几何校正是借助一组地面控制点,对一幅图像进行地理坐标的校正;把影像纳入一个投影坐标系中,有坐标信息地理参考;影像配准是用一影像对另一幅图像的校准,一式两幅图像的同名像元配准 ” ,如果写的对的话,我觉得几何配准其实也是一种几何校正。 我的感觉配准是相对不同影像之间的处理,几何校正是对数据自己的处理。为了进一步明确两者的区别,在 Baidu , Google 上进行了搜索, 校正和影像配准原理是一样的,几何校正是借助一组地面控制点,对一幅图像进行地理坐标的校正;把影像纳入一个投影坐标系中,有坐标信息地理参考;影像配准是用一影像对另一幅图像的校准,一式两幅图像的同名像元配准 几何校正是 指消除或改正遥感影像几何误差的过程。遥感影像的几何畸变,大体分为两类:①内部畸变。由传感器性能差异引起,主要有:比例尺畸变( a ),可通过比例尺系 数计算校正;歪斜畸变( b ),可经一次方程式变换加以改正;中心移动畸变( c ),可经平行移动改正;扫描非线性畸变( d ),必须获得每条扫描线校正数据才 能改正;辐射状畸变( e ),经 2 次方程式变换即可校正;正交扭曲畸变( f ),经 3 次以上方程式变换才可加以改正;②外部畸变。由运载工具姿态变化和目标物引起。包括:由运载工具姿态变化(偏航、俯仰、滚动)引起的畸变,如因倾斜引起的投影畸变( g ),可用投影变换加以校正;因高度变化引起的比例尺不一致 ( h ),可用比例尺系数加以改正;由目标物引起的畸变,如地形起伏引起的畸变( i ),需要逐点校正;若因地球曲率引起的畸变( j ),则需经 2 次以上高次方 程式变换才能加以改正。多光谱、多时相影像配准和遥感影像制图,必须经过上述几何校正。因人们已习惯于用正射投影地图,故多数遥感影像的几何校正以正射投 影为基准进行。 配准是指同一区域内以不同成像手段所获得的不同图像图形的地理坐标的匹配。包括几何纠正、投影变换与统一比例尺三方面的处理。在多时相、多信息的复合综合分析时常需进行各种配准处理。是产生一个空间校准集合或匹配某一区域图像的过程。一般步骤为:①选取同名地物控制点;②输入计算机,实现控制点的相应配准。即作几何纠正、投影变换和比例尺配准处理。配准方法有:①相互配准。以多图像的一个分量作为参考图像,其它图像与其配准;②绝对配准。即定义一个控制格网,使所有图像与其配准。此外,在多光谱影像进行彩色合成时,必须进行不同波段影像的配准,以保证相同景物的有关像元能一一对应。 辐射校正( radiometric correction )是指对由于外界因素,数据获取和传输系统产生的系统的、随机的辐射失真或畸变进行的校正,消除或改正因辐射误差而引起影像畸变的过程。 遥感影像产生辐射误差(即灰度失真)的因素主要有:①大气对电磁波辐射的散射和吸收;②太阳高度与传感器观察角的变化;③地形起伏引起的辐射强度变化;④传感器探测系统性能差异,如光学系统或不同探测器在灵敏度、光谱响应和透光性能上的差异;⑤影像处理,如摄影处理等。影像灰度失真与影像空间频率有关。空间频率愈高,即目标愈小时,辐射误差愈大。辐射校正实际上是影像恢复(或称复原)的一个内容。校正方式有两类:①传感器辐射校正。通常采用内部校准 光源和校准楔,如陆地卫星多光谱扫描仪的辐射校正;②影像辐射畸变校正。常采用物理或数学(校正曲线或各种算法)方法,如空间滤波、平滑化,校正各种灰度 失真及疵点、灰点、条纹、信号缺失等分布在整个影像上的离散形式的辐射误差。其中大气影响的校正还可通过实测反射辐射通量和影像密度,并对数据进行回归分 析来进行校正。 配准与几何校正的区别与联系.doc
个人分类: 专业转载|0 个评论
[转载]TM卫星数据的预处理-大气校正
gosci 2012-2-8 17:10
一、主要内容: 1、ENVI软件常规功能的操作 2、图片重新采样:使不同波段、不同分辩率卫星资料合成可供计算机分析的同一文件内不同波段、同一分辩率的数据文件(ENVI里面的Layer Stacking功能) 3、图片切割:将卫星数据与地图有机配合起来,将卫星数据按照地图所要求的地理空间切割出来,让计算机只分析我们感兴趣的地点(ENVI里面的Masking功能) 4、辐射校正:将卫星传感器捕获的瞬间视场角观测值转化为传感器上空观测的反射辐射值(ENVI里面的Data Special-utilities功能) 5、大气校正:将计算的点的辐射值进行大气校正,消除大气层吸收、散射的影响,获取地面对照地物的反射率(ENVI里面的FLAASH功能) 6、大气校正前后彩色图片的比较 7、不同地物波谱特征的绘制 二、加入分析文件: ①、选择 File Open Image File. ②、当出现 Enter Data Filename 对话框,点击文件名,再点击“OK” 或“Open” 以打开选择的文件。 主图像 Image 窗口:( 400*400 ) 100 %显示(全分辨率显示) scroll 的方框,可交互式分析、查询信息。 滚动 Scroll 窗口: 全局,重采样 ( 降低分辨率 ) 显示一幅图像。 缩放 Zoom 窗口:( 200*200 )显示 image 的方框。 三、图片重新采样 将打开的六个文件的数据放入统一的文件中。 方式:用鼠标点击Basic Tools---Layers Stacking,会弹出下面的窗口: A 、点 Import File 按钮将 6 个文件全选,然后按确定,这样 6 个文件全部选入 selected Files for Layer Stacking 列表中。 B 、这里的顺序颠倒了,是 7 号波段在前面, 1 号波段在后面,这对我们后面的分析不利,为此用鼠标点击 Recorder Files 按钮,对出现的窗口中用鼠标拖动文件上下位置,使波段按 1-7 号排列。 排列好后在“Enter output Filesname”中输入文件名(这里取名为p124r39_5t19920728_1_7),然后按确定键,系统会自动将六个文件重新取样,生成一个文件。 为方便起见,将以前的六个文件关闭,步骤如下: 在Available Band List 窗口中,选中文件p124r39_5t19920728_nn1,用鼠标点击右键,在弹出的窗口中选Close Selected File;再依次关掉2,3,4,5,7文件。最后Available Band List 窗口便将所有六个文件集中为一个包含多光谱的文件。 四、头文件的修改 : 由于进行重新采样,文件里的相关信息都消失了,为了便于计算机分析,我们要将相关的信息(卫星类型、波段)进行修改,步骤如下: A、打开头文件修改窗口:方式一:在Available Band List 窗口中,选中文件p124r39_5t19920728_1_7,用鼠标点击右键,在弹出的窗口中选Edit Header.方式二:打开菜单 File---Edit ENVI Header,在窗口中选取文件p124r39_5t19920728_1_7。 B、修改波段名称: 点击Edit Attributes按钮,选取Band Names 。在修改框中分别写入Band1、Band2、……Band5、Band7。 C、修改各波段波长: TM卫星各波段波长分布如下表: 波段号 波长范围(微米) 中间值(微米) 波段对应颜色 1 0.45-0.52 0.485 蓝(B) 2 0.52-0.60 0.560 绿(G) 3 0.63-0.69 0.660 红(R) 4 0.76-0.90 0.830 近红外(NIR) 5 1.55-1.75 1.650 中红外(MIR) 7 2.08-2.35 2.215 中红外(MIR) 点击Edit Attributes按钮,选取WaveLengths弹出如下窗口,用鼠标选中待改的波段号,输入波段范围的中间值,最后选择单位,MicroMeters,按OK键完成。 五、显示彩色图片: 可以根据B、G、R三种颜色合成彩色图(也可以任选其它三个波段合成伪彩色图) 方式一:在Available Band List 窗口中,选中文件p124r39_5t19920728_1_7,用鼠标点击右键,在弹出的窗口中选Load True Color. 方式二:在Available Band List 窗口中点击RGB color,分别选取3、2、1号波段,点display按钮。 六、与电子地图结合,将分析区域提取出来。 A、点击主菜单Window-----Available Vector List,弹出图。在当前窗口选取File----open Vector File弹出窗口,下拉文件类型(*.Shp),选取要求分析的县名,便可,系统会弹出一个类型转换对话框将SHP文件转成evf类型,你按OK系统便自动转换,参数为计算机缺省,这样一个矢量文件便调入ENVI软件中。 B、Masking生成:MASK是一个模板,可以运用这个模板将要求的数据取出来,下面将建一个县行政边界的模板,表示县内的数据要取出来分析。 在主菜单中点击Basic Tools-----Masking---Build MASK,在弹出的对话框中Display1操作。 ①在Options---Imports EVFs ②选取矢量文件后按OK ③选择保存在内存,按Apply便可。 C、MASK应用:在主菜单中点击Basic Tools----Masking---Apply MASK,弹出窗口。 ①点击 Spatial SubSet( 空间子集 ) ,表示不将全图抄本,只计算里面一小部分,在弹出的窗口中点击 ROI/EVF ,选取矢量图,按 OK。 ②在 Select Mask Band 中选中刚才做好的 MASK ( memory1 ) ③这样我们便将所有参数设置好了,按 OK ,取一个名 p124r39_5t19920728_1_7_Jz ,按 OK ,便好了。 这样我们便将只有这个县的分析图片做出来了,同样为了分析方便 在 Available Band List 窗口中将 p124r39_5t19920728_1_7 、 memory1 关掉,这样 Available Band List 窗口。 点文件的右键,在弹出的窗口中选中 display True Color。 七 、反射辐射值的计算: 上图显示的值是卫星传感器的感应值,并非辐射值,在SCROLL或者IMAGE窗口中点鼠标右键,在弹出的窗口中点选Current Location/Value,你会发现窗口,当你的鼠标在图片上移动时,鼠标点处的位置(经、纬度)和观测值就会实时显示出来。其DATA值都是整数,是传感器的感应值,而非辐射的绝对值,我们可以利用ENVI中的辐射较正功能进行转换,具体步骤如下: 在主菜单中选中 Basic Tools----Preprocessing----Calibration Utilities----Landsat TM,在弹出的窗口选中p124r39_5t19920728_1_7_JingZhou后出现对数设置窗口,按如下设置: LandSat Satellite 选5 时间选:1992年7月28日 太阳高度角选:51.25 校正类型选:辐射校正 保存文件:p124r39_5t19920728_1_7_Jz_R(可以随便取名)按OK键便可,同样为了分析方便在Available Band List 窗口中将p124r39_5t19920728_1_7_JZ关掉,这样Available Band List 口。 点文件的右键,在弹出的窗口中选中display True Color,并在图片显示窗口中用鼠标点右键,在弹出的菜单中选择Z Profile(spectrum),系统会出现一个光谱窗口,实时显示某一地物传感器观测到的光谱图,这时观测值为实时,单位是:mW/(cm2*sr*mm)。 八 、地物反射率计算: 上面是传感器得到的反射辐射值,它有可能是来自地物,也有可能来自大气的散射辐射(这部分是无用的),为止我们要进行大气校正,通过传感器得到的反射辐射值来计算地物的反射率,即FLASSH功能。操作步骤如下: A、由于FLAASH功能转换时要求卫星数据以BIP或者BIL方式保存,所以应首先将p124r39_5t19920728_1_7_JZ_R转成BIP格式的文件p124r39_5t19920728_1_7_JZ_RB 在主菜单中选中 Basic Tools----Convert Data (BSQ、BIP、BIL),选中文件p124r39_5t19920728_1_7_JZ_R后出现窗口: 参数选取 格式:BIL 文件名:p124r39_5t19920728_1_7_JZ_RB 按OK,便生成了以BIL格式存放的卫星数据 B、FLASSH功能:在主菜单中选中 Basic Tools----Preprocessing----Calibration Utilities----Flaash 按Input Radiance Image按钮输入用于计算的源图片:p124r39_5t19920728_1_7_JZ_RB,在弹出下面的参数选择窗口中选下面的值: 输出文件选择:p124r39_5t19920728_1_7_Finale 工作目录:选择p124r39_5t19920728_1_7_Finale所在的目录 FLAASH前缀选择TEMP_ 图片中心位置选择:纬度30度9分25微分 经度:112度1分45微分 卫星类型:多光谱的LandSat TM5 日期时间:1992年7月28日02时36分 大气模式:中纬度夏天 微粒模式选无 最后按apply便可。 在Available Band List 窗口点p124r39_5t19920728_1_7_Finale文件右键,在弹出窗口中选Load True Color. 在图片显示窗口中用鼠标点右键,在弹出的菜单中选择Z Profile(spectrum),再在图片显示窗口中用鼠标点右键,在弹出的菜单中选择Current Location/Value选取一个你确认的地物(城市、水体、植物)点击一下,转入光谱窗口:点菜单EDIT---Data Values,将数据记录下来: 选取几个特定地物后,在实验纸上以图的形式绘制出来. from:http://lisyangtzeu.blog.sohu.com/151554656.html
个人分类: ENVI|1135 次阅读|0 个评论
[转载]BFGS算法C程序
ywmucn 2011-11-13 00:24
/* * BFGS.c * * Created on: 2009-6-4 * Author: Administrator */ //×××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××× //××××××××××××××××××××××××××××BFGS算法C程序××××××××××××××××××××××××××××××× //×××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××× //×××××××××××本程序适用于2设计变量的函数优化问题,对于不同的设计变量××××××××× //***********个数可以改变维数,对于DFP算法只需修改校正矩阵即可。对于×××××××××× //×××××××××××BFGS算法的C++程序与之类似。 ×××××××× //*××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××× //×××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××× //************************************************************************ #includestdio.h #includestdlib.h #includemath.h #includeconio.h #define tt 0.01 //-------------------一维搜索初始步长--------------------------- #define ff 1.0e-6 //-------------------差分法求梯度时的步长----------------------- #define ac 1.0e-6 //-------------------终止迭代梯度模收敛精度--------------------- #define ad 1.0e-6 //-------------------一维搜索收敛精度--------------------------- #define ax 1.0e-6 //-------------------终止迭代点距收敛精度----------------------- #define n 2 //-------------------设计变量的维数----------------------------- double ia; /*FILE *fp;*/ //=================================输入目标函数表达式======================================== double fny(double *x) { double x1=x ,x2=x ; double f; f=x1*x1+2*x2*x2-4*x1-2*x1*x2; return f; } //========================================================================================== //==================================迭代更新最优点=========================================== double * iterate(double *x,double a,double *s) { double *x1; int i; x1=(double *)malloc(n*sizeof(double)); for(i=0;in;i++) x1 =x +a*s ; return x1; } //========================================================================================== //===================================计算更新后节点函数值==================================== double func(double *x,double a,double *s) { double *x1; double f; x1=iterate(x,a,s); f=fny(x1); return f; } //========================================================================================== //===================================确定初始单谷区间的退进法================================ void finding(double a ,double f ,double *xk,double *s) { double t=tt; int i; double a1,f1; a =0;f =func(xk,a ,s); for(i=0;;i++) //循环次数i,受break执行的影响 { a =a +t; f =func(xk,a ,s); if(f f ) break; if(fabs(f -f )=ad) { t=-t; a =a ;f =f ; } else { if(ia==1) return; //break t=t/2;ia=1; } } for(i=0;;i++) { a =a +t; f =func(xk,a ,s); if(f f ) break; t=2*t; a =a ;f =f ; a =a ;f =f ; } if(a a ) { a1=a ; f1=f ; a =a ; f =f ; a =a1; f =f1; } return; } //====================================================================================== //=======================================一维搜索======================================== double lagrange(double *xk,double *ft,double *s) { int i; double a ,f ; double b,c,d,aa; finding(a,f,xk,s); for(i=0;;i++) { if(ia==1) { aa=a ; *ft=f ; break; } d=(pow(a ,2)-pow(a ,2))*(a -a )-(pow(a ,2)-pow(a ,2))*(a -a ); if(fabs(d)==0) break; c=((f -f )*(a -a )-(f -f )*(a -a ))/d; if(fabs(c)==0) break; b=((f -f )-c*(pow(a ,2)-pow(a ,2)))/(a -a ); aa=-b/(2*c); *ft=func(xk,aa,s); if(fabs(aa-a )=ad) {if(*ftf ) aa=a ;break;} if(aaa ) { if(*ftf ) {a =aa;f =*ft;} else if(*ftf ) {a =a ;a =aa;f =f ;f =*ft;} else if(*ft==f ) { a =aa;a =a ; f =*ft;f =f ; a =(a +a )/2; f =func(xk,a ,s); } } else { if(*ftf ) {a =aa;f =*ft;} else if(*ftf ) {a =a ;a =aa;f =f ;f =*ft;} else if(*ft==f ) {a =aa;a =a ; f =*ft;f =f ; a =(a +a )/2; f =func(xk,a ,s); } } } if(*ftf ) {*ft=f ;aa=a ;} return aa; } //======================================================================================= //====================================计算向量梯度======================================== double *gradient(double *xk) { double *g,f1,f2,q; int i; g=(double*)malloc(n*sizeof(double)); f1=fny(xk); for(i=0;in;i++) {q=ff; xk =xk +q; f2=fny(xk); g =(f2-f1)/q; xk =xk -q; } return g; } //======================================================================================== //=======================================BFGS算法主程序======================================================= double * bfgs(double *xk) { double u ,v ,h ,dx ,dg ,s ; double aa,ib; double *ft,*xk1,*g1,*g2,*xx,*x0=xk; double fi,fxk1; int i,j,k,step=1; ft=(double *)malloc(sizeof(double)); xk1=(double *)malloc(n*sizeof(double)); //=====================================数据初始化============================================== for(i=0;in;i++) //--------------尺度矩阵h 、搜索方向s 初始化-------- { s =0; for(j=0;jn;j++) { h =0; if(j==i) h =1; } } //--------------------------------------------------------- g1=gradient(xk); //----------------向量梯度初始化g1 ---------------------- fi=fny(xk); x0=xk; printf("%d %f %f %f \n",step,x0 ,x0 ,fi); //============================================================================================ //=========================================迭代求解最优点====================================== for(k=0;kn;k++) { ib=0; if(ia==1) { xx=xk; break; } ib=0; for(i=0;in;i++) s =0; //----------------求搜索方向s --------------------- for(i=0;in;i++) for(j=0;jn;j++) s += -h *g1 ; //--------------------------------------------------- aa=lagrange(xk,ft,s); //----------------计算步长a-------------------------- xk1=iterate(xk,aa,s); //----------------迭代计算X ----------------------- g2=gradient(xk1); //----------------计算向量梯度g ------------------- fxk1=fny(xk1); step+=1; printf("%d %f %f %f \n",step,xk1 ,xk1 ,fxk1); if(sqrt((g2 -g1 )*(g2 -g1 )+(g2 -g1 )*(g2 -g1 ))=ac||sqrt((xk1 -xk )*(xk1 -xk )+(xk1 -xk )*(xk1 -xk ))=ax) //---------------------判断迭代是否结束---------------------- {ib=ib+1;} if(ib==0) { xx=xk1; break; } //---------------------------------------- fi=*ft; if(k==0)//------------------------------求解校正矩阵------------------------------- { int j; double a1=0,a2=0; for(i=0;in;i++) { dg =g2 -g1 ; dx =xk1 -xk ; } for(i=0;in;i++) { int j; u =0;v =0; for(j=0;jn;j++) { u =u +dg *h ; v =v +dg *h ; } } for(j=0;jn;j++) { a1+=dx *dg ; a2+=v *dg ; } if(fabs(a1)!=0) { a2=1+a2/a1; for(i=0;in;i++) for(j=0;jn;j++) h +=(a2*dx *dx -v *dx -dx *u )/a1; } //-------------------------------------------------------------------------------- xk=xk1; g1=g2; //------------------------------------更新求解向量和向量梯度----------------------------------- k=-1; } } //=============================================================================== if(*ftfi) { *ft=fi; xx=xk;} xk=x0; return xx; } //===================================================================================== main () { int k; double *xx,f; double xk ={1,1}; xx=bfgs(xk); f=fny(xx); printf("\n\nThe Optimal Design Result Is:\n"); for(k=0;kn;k++) {printf("\n\tx *=%f",k+1,xx );} printf("\n\tf*=%f",f); getch(); }
个人分类: 编程|3218 次阅读|0 个评论
[转载]环境星HSI数据FLAASH大气校正
lishumin 2011-10-12 11:09
[转载]环境星HSI数据FLAASH大气校正
来自: http://bbs.esrichina-bj.cn/ESRI/thread-102714-1-1.html 看到有不少同学在用环境星影像做分析,在大气校正这块碰到不少问题,我自己闲暇也做了一些工作,这里写个帖子与大家分享下。 关于环境星数据的介绍大家可以点击这里: http://www.cresda.com/n16/n1130/n1582/8384.html 由于 CCD 影像的 FLAASH 大气校正与 TM 影像的大同小异,我就不做介绍了,版主写过一篇关于 CCD 帖子: http://bbs.esrichina-bj.cn/ESRI/thread-77857-1-1.html 这主要说下关于 HSI 数据 的。 HSI 为高光谱数据 ,在资源卫星应用中心下载的是 2 级产品, 数据的格式为 HDF5 ,由于大部分同学的 ENVI 没有安装 HDF 插件,无法进行直接操作,可以下载一个转换工具, 从 HDF5 转到 TIF ,或者在资源卫星中心那下载软件( http://www.cresda.com/n16/n1115/n1522/n2164/84990.html )。 2 级数据已经是辐亮度 ,所以 不需要再做辐射定标 了,从 DN 值转到辐亮度这个过程就可以省去 。 使用 FLAASH 大气校正之前,数据必须转为 BIL 或 BIP 格式 。环境卫星的传感器的类型在 FLASSH 中并没有设置,所以要 手动在影像的头文件中添加中心波长和 FWHM 。选择 Edit Header… 点击 Edit Attributes 之后选择 FWHM 之后会出现下面的对话框,点击 ImportASCII… , 出现一个选择文件的对话框,选择 XXXXXXX_BandInfo.txt 这个文件,之后出现右面的对话框,在 Wavelength Column 后面填 1 ,在 FWHM Column 后面填写 3 。要求的是中心波长,我直接选择的是起始波长,不知道影响大不 \(^o^)/~ 。。毕竟只是做个试验,如果大家要求比较精确的话,可以求下中心波长。 前期的处理是做完了,下面就进行大气校正,要求严格的话,其实应该先做下条纹去除。毕竟环境星的条纹还是挺多的。我这就不做这项工作了。 \(^o^)/~ 资源卫星中心在生成辐亮度产品时乘上了 100 ,单位是( W ) / ( m2*nm*sr ), FLASSH 要求的单位为( μW ) / ( cm2*nm*sr ),根据 1μW/ ( cm2*nm*sr ) =10 W/(m2*um*sr) , 所以这个因子应该是 1000 !!!!! 但是在我处理的过程中,我发现情况并不是这样的,有两景影像都为 2009-12-06 获得的,其中一幅影像名称为: HJ1A-HSI-2-66-A1-20091206-L20000213407 , 不知道是什么原因因子是 10 ,而不是 100 !!所以在填写这个因子的时候应该是 100 ,如果还写 1000 的话,会发现结果全部为 0 。 有一个比较简单的方法判断就是做下统计,如果发现 DN 值没有超过 1000 的,就可以基本可以确定因子是 10 而不是 100 了!我想资源卫星不会还是手动处理这么大批量的数据吧?? 之前我在用 IRS 数据的时候也发现了一个问题, IRS 数据的第四波段说明是 300M 分辨率,不知为何有的数据竟然是 150M 分辨率。 影像的中心经纬度,扫描时间 可以从随带的 XML 文件中得到 ,但是不知扫描时间给出的是否为格林尼治时间,都是在凌晨 3 点左右,如果加上 8 的话,正好是北京时间 11 点多,估计应该是吧。 卫星的传感器海拔为 650KM 。研究区的平均海拔高度可以根据 DEM 数据获得,这里直接写 0.1 。 大气模式和气溶胶模式根据情况来进行选择。 HSI 的波段只达到 950 ,无法进行气溶胶反演,所以 气溶胶反演选择 NONE 。在水汽反演中,可以选择 820 波长的进行。在高光谱的设置中,可以选择 AutomaticSelection ,让其根据波谱信息自动进行云和水汽的反演。 高级设置就选择默认的吧。之后选择 APPLY 就可以了。 !环境星HSI数据FLAASH大气校正.doc
个人分类: 专业转载|0 个评论
[转载]FLAASH大气校正常见错误及解决方法
lishumin 2011-9-22 09:55
[转载]FLAASH大气校正常见错误及解决方法
http://bbs.esrichina-bj.cn/ESRI/thread-83384-1-6.html 本文汇总了 ENVI FLAASH 大气校正模块中常见的错误,并给出解决方法,分为两部分:运行错误和结果错误。前面是错误提示及说明,后面是错误解释及解决方法。 FLAASH 对输入数据类型有以下几个要求: 1、 波段范围:卫星图像: 400 - 2500nm ,航空图像: 860nm-1135nm 。如果要执行水汽反演,光谱分辨率 =15nm ,且至少包含以下波段范围中的一个: l 1050-1210 nm l 770-870 nm l 870-1020 nm 2、 像元值类型:经过定标后的辐射亮度(辐射率)数据,单位是:( μ W ) / ( cm2*nm*sr )。 3、 数据类型:浮点型( Floating Point )、 32 位无符号整型( Long Integer )、 16 位无符号和有符号整型( Integer 、 Unsigned Int) ,但是最终会在导入数据时通过 Scale Factor 转成浮点型的辐射亮度( μ W ) / ( cm2*nm*sr )。 4、 文件类型: ENVI 标准栅格格式文件, BIP 或者 BIL 储存结构。 5、 中心波长:数据头文件中(或者单独的一个文本文件)包含中心波长( wavelenth )值,如果是高光谱还必须有波段宽度( FWHM ),这两个参数都可以通过编辑头文件信息输入( Edit Header )。 运行错误 1. Unable to write to this file.File or directory is invalid or unavailable 。 没有设置输出反射率文件名。 解决方法是单击 Output Reflectance File 按钮,选择反射率数据输出目录及文件名,或者直接手动输入。 2. ACC Error : convert7 IDL Error : End of input record encountered on file unit:0. 平均海拔高程太大。 注意:填写影像所在区域的平均海拔高程的单位是 km : Ground Elevation ( Km )。 3. ACC error : avrd: IDL error : Unable to allocate memory:to make array Not enough space ACC_AVRD 为了能处理大数据, ENVI 采用分块计算的方式,这个提示是分块( Tile )太大了。 在高级设置里面( Advanced Settings ), tile size :它默认是 File- preferences -miscellaneous:cache 的大小,这个值正常是 1~4Mb ( 0 背景很少的情况下) ; 如果 0 背景较多,这个值还是需要设置大一些,比如 100~200m 。 4. ACC error:avrd:No nonblank pixels found IDL error: OPENR: Error opening file. Unit: 100, 为了能处理大数据, ENVI 采用分块计算的方式,当 Tile 设置太小,而且有背景值( 0 ),就会出现一个 Tile 中全部为 0 的情况,提示这个错误信息。 在高级设置里面( Advanced Settings ), tile size :设置稍微大一些,如 100~200M 等。 5. ACC error:lsmooth2: IDL error: ACC_LSMOOTH2:Cannot continue with smoothing calculation ENVI 的 FLAASH 提供领域纠正功能,但是 MODIS 、 AVHRR 等图像分辨率比较低,领域效应区分不出来。 解决方法是在在高级设置里面( Advanced Settings ),将领域纠正( Use Adjacency Correction )设置为 No 。 6. ACC error:modrd5:Nonfinite numbers in coefficient array coef IDL error: CDRIVER4V3R2:Wait for MODTRAN4 calculation to finish… 提示传入 MODTRAN 模型参数有误,常常是由于太阳高度角太小或者太大引起的。提示这个错误之前会出现以下提示框。 解决方法是确认填写的影像中心位置经纬度信息(西经为负数、南纬为负数)、影像成像时间(格林威治时间)是否正确。 结果错误 1. 结果中某一个波段或者多个波段全部为 0 或者负值 主要产生的原因是输入的辐射亮度数据值偏小。 可能有以下几种情况: (1) 在传感器定标的时候选择的是表观反射率( Reflectance )而不是辐射率数据( Radiance ); (2) 没有做传感器定标,即没有将 DN 值转换为辐射率数据; (3) 选择了错误的波谱响应函数; (4) 用 BandMath 做了辐射亮度的单位换算,在 FLAASH 中导入辐射亮度数据时, Scale Factor 选择的不是 1.0000 。 2. 结果中极大值、极小值非常多,也就是 0~10000 之外的值。当选择 RGB 假彩色显示的时候,出现花花绿绿的情况。 辐射定标得到的辐射率数据单位与 FLAASH 要求的单位不一样。 可能有以下几种情况: (1) 用 BandMath 做了辐射亮度的单位换算,在 FLAASH 中导入辐射亮度数据时, Scale Factor 选择的不是 1.0000 。 (2) 没有做单位换算。 3. 结果中部分像元为负值 这个属于正常现象。 FLAASH 是采用 MODTRAN 辐射传输模型模拟成像中的大气过程,而且很多大气属性都是通过图像来估算,加上大气组成的非均一性,即使 MODTRAN4 模型精度很高,也不能完全表达大气辐射传输的真实状况。当影像上有强吸收或者高反射地物时候,就会出现部分像元为负值(如深水、高密度)或者大于 10000 。 解决方法可以手动修改,如用周围的像元的平均值代替,可以使用 ENVI 下的 The DEM Editing Tool 工具。 4. 结果图像以 RGB 显示比原图像视觉效果要差,如模糊。 这个需要了解下遥感软件 RGB 显示机制。遥感软件为了让遥感图像显示更加“艳丽”,方便解译,默认会对图像进行拉伸显示,一般是 2% 的线性拉伸。也就是我们常常会看到图像值有两个,一个是 Scrn 值,也就是拉伸之后的显示值,一个是 Data 值,也就是原始的图像 DN 值。另外一个方面,目前我们的 RGB 加色法显示都是基于 8bit 显示,也就是 0~255 。 FLAASH 大气校正之后的结果是 16bit 的整型,而且存在一些极小、极大值,这些对直方图整体形状有一定的影响,影响拉伸效果。解决方法是利用 ENVI 下的 Interactive stretching 工具,选择有效值范围进行拉伸。 还有一种情况是校正图像有很多背景值,比如经过几何校正的整景 TM 影像。背景及图像边缘处的像元在大气校正之后变成负值或者 0 值,由于这样的像元数量多,对图像的整体拉伸影响很大。解决方法就是将这些像元掩膜掉,如在 ROI Tool 中,利用 Option-Band Threshold to ROI 建立一个 ROI 进行图像裁剪。 5. 结果图像为什么像元值大多是大于 1 ENVI FLAASH 考虑到数据储存和后续处理,将大气校正得到的反射率结果乘以 10000 变成 16bit 整型。 如果想让反射率结果在 0~1 范围,可用 BandMath ,表达式为 b1/10000.0 。
个人分类: 专业转载|0 个评论
[转载]环境小卫星CCD影像的精确大气校正
lishumin 2011-9-22 09:50
http://blog.sina.com.cn/s/blog_764b1e9d0100qsnl.html 国产卫星技术越来越成熟,使用范围也越来越广。包括环境小卫星和资源卫星都拥有多光谱传感器,常用于定量遥感。本专题介绍利用 ENVI FLAASH 工具完成环境小卫星(资源卫星类似)的大气校正 . 1 数据定标 环境小卫星的 CCD 相机,利用绝对定标系数将 DN 值图像转换为辐亮度图像的公式为 L=DN/a + L 0 式中 L 为辐亮度, a 为绝对定标系数增益, L 0 为偏移量,转换后辐亮度单位为 W ⋅ m − 2 ⋅ sr − 1 ⋅ μ m − 1 。 定标系数直接可以从元数据文件( .xml )中找到如下列: absCalibType(gain1,Fielddata,L=DN/g+L0,W*m^(-2)*sr^(-1)*um^(-1))B1:g 0.5782, L0 3.4608, B2:g 0.5087, L0 5.8769, B3:g 0.6825, L0 8.0069, B4:g 0.6468, L0 8.8583/absCalibType 也可以参照表 1 所示。 表 1 HJ1A/B 星 CCD 绝对辐射定标系数 卫星 传感器 增益 参数 定标系数 Band1 Band2 Band3 Band4 HJ1A CCD1 1 a ( DN/W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 0.5763 0.5410 0.6824 0.7209 L 0 (W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 9.3183 9.1758 7.5072 4.1484 2 a ( DN/W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 0.9160 0.9228 1.1277 1.0753 L 0 (W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 7.3250 6.0737 3.6123 1.9028 CCD2 1 a ( DN/W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 0.6360 0.5910 0.8142 0.8768 L 0 (W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 7.5575 7.0944 4.1319 1.2232 2 a ( DN/W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 0.9997 1.0016 1.3777 1.3043 L 0 (W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 4.6344 4.0982 3.7360 0.7385 HJ1B CCD1 1 a ( DN/W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 0.5329 0.52895 0.68495 0.72245 L 0 (W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 1.6146 4.0052 6.2193 2.8302 2 a ( DN/W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 0.8685 0.9367 1.2433 1.3002 L 0 (W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 3.0089 4.4487 3.2144 2.5609 CCD2 1 a ( DN/W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 0.5782 0.5087 0.6825 0.6468 L 0 (W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 3.4608 5.8769 8.0069 8.8583 2 a ( DN/W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 0.9076 0.8502 1.1635 0.9800 L 0 (W ⋅ m − 2 ⋅ sr − 1 ⋅ μm − 1 ) 2.2219 4.0683 5.2537 6.3497 注:增益 1 的定标系数是通过实验室定标系数得到的,增益 2 状态的定标系数为场地替代定标获取。 (一) 有了绝对定标参数和定标公式,选择 Basic Tool- Band math 工具很容易进行传感器定标。 (二) 由于是单个波段文件的定标,选择 Basic Tool-Layer Stacking 将定标后的单波段文件组合成一个多波段文件。 (三) 选择 Basic Tool-Convert Data ( BSQ,BIL,BIP ),将定标后的文件转换为 BIL 格式文件。 (四) 在波段列表中,右键打开 Edit header-Edit Attributes-Wavelength ,将每个波段的中心波长输入—— b1(475nm) , b2(560nm) , b3(660nm) , b4(830nm) 目前为止,影像数据准备工作完成了。 2 波谱响应函数 波段响应函数,英文名为 spectral response function ( SRF )或者叫 Relative Spectral Response (RSR) ,与宽波段传感器出现。我们知道每一个波段都有一个波段范围,比如 HJ-b1(475nm) 波段为 520-430 纳米,实际上传感器的感光元件在这个波段范围内的每一点所感应的强度都是不一样。在成像中,原则上讲应该根据波段响应函数来进行加权平均,但由于处理起来比较麻烦,而且一般的精度要求不太高,所以大多数图像都是直接取了波段范围内的中点值来运算。波段响应函数是描述一定波长范围内(超出波段范围)的量子效应,当需要精确计算像元响应时候,比如大气校正反演真实地表反射率,就需要使用波谱响应函数。 任何传感器在设计时都会给出严格的波段响应函数,环境小卫星也不例外。如下图中为波谱响应函数数字表达的一部分,第一列表示波长,后面四列分别表示 4 个波段对应波长的波谱响应值。 图 1 波谱响应值 在 ENVI 中,使用波谱曲线来描述波谱响应函数,也就是以波长作为 x 轴,波谱响应值作为 y 轴,存储格式为 ENVI 波谱库文件 (.sli) 。下面介绍 ENVI 中的波谱响应函数的制作。 (一) 选择 Window-Start New Plot Window , ENVI Plot Window 窗口中,选择 File-Input Data-ASCII, 如图 2 所示,自动将第一列作为 X 轴,后面 4 列作为 Y 轴。单击 OK 。 (二) 如图 3 所示,生成了 4 条曲线。选择 Edit-Data Parameters ,更改每一条曲线的名称: b1,b2,b3,b4 ,便于区分。 (三) 选择 File-Save Plot As-Spectral Library ,将波谱曲线保存为波谱库文件。 图 2 导入 ASCII 文件 图 3 波谱响应函数的曲线表达 3 FLAASH 大气校正 启动 FLAASH 模块,将准备好的数据输入,注意定标后的单位是 W ⋅ m −2 ⋅ sr −1 ⋅ μm −1 ,与 FLAASH 要求的单位( μ W ) / ( cm 2 *nm*sr )相差 10 倍关系,因此在 Radiance Scale Factors 中输入 10 缩放系数。 其他特殊参数设置主要如下: l 一般参数设置 传感器高度( Sensor Altitude ): 650km 像元大小( pixel Size ): 30m l 多光谱数据参数设置 Fileter Function File :选择前面生成的波谱响应文件。 其它参数参考 FLAASH 大气校正的参数说明。 4 总结 从上面可以看到,对于未知的多光谱传感器的 FLAASH 大气校正,关键是获取波谱响应函数,而波谱响应函数一般都是公开的。类似的方法可用于资源卫星的 FLAASH 大气校正。
个人分类: 专业转载|0 个评论
[转载]LANDSAT5-TM 的ENVI FLAASH模块大气校正
gosci 2011-7-17 22:59
LANDSAT5-TM的ENVIFLAASH模块大气校正【转载】 相关文献1. http://dufan20086.blog.163.com/blog/static/6616452320110229957123/ 一、辐射定标 1. 由于ENVI 4.4 中有专门进行辐射定标的模块,因此实际的操作十分简单。将原始TM 影像打开以后,选择 Basic Tools–Preprocessing–Calibration Utilities–Landsat TM 2. 进入下一步参数选择:根据传感器类型选择Landsat 4,5 或者7。从遥感影像的头文件中获取Data Acquisition 的时间,Sun elevation。如果你是用File–Open External File–Landsat–Fast 的方法打开header.dat 的话,sun elevation 就已经填好了。这里Calibration Type 注意选择为Radiance。输出文件,定标就完成了。 二、大气校正 简单一点的大气校正可以采用ENVI的FLAASH模块,以下就是FLAASH操作的步骤: 1. FLAASH 模块的进入方法是Spectral–FLAASH,或者是Basic Tools–Preprocessing–Calibration Utilities–FLAASH。 2. FLAASH 模块的操作界面分为三块:最上部设定输入输出文件;中间设定传感器的参数;下部设定大气参数。 3. 首先设定输入输出文件。FLAASH 模块要求输入辐亮度图像,输出反射率图像。之前我们进行了辐射定标,得到辐亮度图像,在这里要把BSQ 格式的图像转换为BIL 或者BIP 格式的图像,然后再Input Radiance Image 中选择转换格式后的图像。(Basic Tools–Convert Data(BSQ,BIL,BIP))。这里注意,当输入图像后,程序会让你选择Scale Factor,即原始辐亮度单位与ENVI 默认辐亮度单位之间的比例。ENVI 默认的辐亮度单位是μW/cm2 •sr•nm,而之前我们做辐射定标时单位是W/m2 •sr•μm,二者之间转换的比例是10,因此在下图中选择Single scale factor,填写10.000。 4. 此外,如果TM 影像的头文件中没有波段的信息,在这里也要求你提供一个.txt 文件以包含此信息。那么,准备好一个.txt 文件,其中含有一列TM 每个波段中心波长的信息。 5. 在Output Reflectance File 和Output Directory for FLAASH files 里面设定输出文件的文件名和位置。 6. 设定传感器参数。首先是Scene Center Location,即遥感图像中心的坐标,以及Flight Date, Flight Time GMT,这三者都可以在TM 的头文件中找到,填入即可。 7. 在Sensor Type 菜单中选择Landsat TM5。此时Sensor altitude 自动填上为705km。而Pixel Size 填为30m。 8. 根据遥感影像研究区实际情况,填写Ground Elevation,比如华北平原可以写为0.05km。 9. 最关键的为大气参数部分: a) Atmospheric Model( 大气模式): 共有Sub-Arctic Winter (SAW) ,Mid-Latitude Winter (MLW),U.S. Standard (US) ,Sub-Arctic Summer(SAS), Mid-Latitude Summer (MLS) 和Tropical (T) 。根据经纬度和时间可以选定研究区的大气模式,见ENVI Help。 b) Aerosol Model(气溶胶模式):有Rural, Urban, Maritime 和Tropospheric四种选择。根据实际情况选择即可。关于此四种模式的解释见ENVI Help。 c) 当我们选择TM 时,可选的参数还有Aerosol Retrieval 和Initial Visibility。这两个参数对最后的结果又相当重要的影像,因此最好能调查到当地的Initial Visibility。此外,AERONET 在全世界各地有测定AOD(Atmospheric Optical Depth)的站点,可以查询AOD 以后转换为消光系数,通过消光系数估算能见度,此步骤比较繁琐,在此不予详述。如果采用Aerosol Retrieval 中的K-T算法计算Visibility,且能够计算出结果的话,则采用K-T 算法的能见度,否 则采用Initial Visibility 所指定的能见度。 d) 关于Aerosol Retrieval。如果选择了下拉菜单中的K-T method,那么需要在Multispectral Settings 中设定参数,在Assign Default Values Based on Retrieval Conditions 中选择Over-land Retrieval Standard (660:2100nm)即可。根据不同的研究区可以设定不同的模式。其他设定可以不改变。Apply即可。 http://blog.sina.com.cn/s/blog_4d3ef47e0100cdxy.html
个人分类: ENVI|1284 次阅读|0 个评论
[转载]ENVI中TM的辐射定标和大气校正
KeXinli 2011-4-26 08:50
一、辐射定标 1. 由于ENVI 4.4 中有专门进行辐射定标的模块,因此实际的操作十分简单。将原始TM 影像打开以后,选择 Basic Tools–Preprocessing–Calibration Utilities–Landsat TM 2. 进入下一步参数选择:根据传感器类型选择Landsat 4,5 或者7。从遥感影像的头文件中获取Data Acquisition 的时间,Sun elevation。如果你是用File–Open External File–Landsat–Fast 的方法打开header.dat 的话,sun elevation 就已经填好了。这里Calibration Type 注意选择为Radiance。输出文件,定标就完成了。 二、大气校正 简单一点的大气校正可以采用ENVI的FLAASH模块,以下就是FLAASH操作的步骤: 1. FLAASH 模块的进入方法是Spectral–FLAASH,或者是Basic Tools–Preprocessing–Calibration Utilities–FLAASH。 2. FLAASH 模块的操作界面分为三块:最上部设定输入输出文件;中间设定传感器的参数;下部设定大气参数。 3. 首先设定输入输出文件。FLAASH 模块要求输入辐亮度图像,输出反射率图像。之前我们进行了辐射定标,得到辐亮度图像,在这里要把BSQ 格式的图像转换为BIL 或者BIP 格式的图像,然后再Input Radiance Image 中选择转换格式后的图像。(Basic Tools–Convert Data(BSQ,BIL,BIP))。这里注意,当输入图像后,程序会让你选择Scale Factor,即原始辐亮度单位与ENVI 默认辐亮度单位之间的比例。ENVI 默认的辐亮度单位是μW/cm2 ?sr?nm,而之前我们做辐射定标时单位是W/m2 ?sr?μm,二者之间转换的比例是10,因此在下图中选择Single scale factor,填写10.000。 4. 此外,如果TM 影像的头文件中没有波段的信息,在这里也要求你提供一个.txt 文件以包含此信息。那么,准备好一个.txt 文件,其中含有一列TM 每个波段中心波长的信息。 5. 在Output Reflectance File 和Output Directory for FLAASH files 里面设定输出文件的文件名和位置。 6. 设定传感器参数。首先是Scene Center Location,即遥感图像中心的坐标,以及Flight Date, Flight Time GMT,这三者都可以在TM 的头文件中找到,填入即可。 7. 在Sensor Type 菜单中选择Landsat TM5。此时Sensor altitude 自动填上为705km。而Pixel Size 填为30m。 8. 根据遥感影像研究区实际情况,填写Ground Elevation,比如华北平原可以写为0.05km。 9. 最关键的为大气参数部分: a) Atmospheric Model( 大气模式): 共有Sub-Arctic Winter (SAW) ,Mid-Latitude Winter (MLW),U.S. Standard (US) ,Sub-Arctic Summer(SAS), Mid-Latitude Summer (MLS) 和Tropical (T) 。根据经纬度和时间可以选定研究区的大气模式,见ENVI Help。 b) Aerosol Model(气溶胶模式):有Rural, Urban, Maritime 和Tropospheric四种选择。根据实际情况选择即可。关于此四种模式的解释见ENVI Help。 c) 当我们选择TM 时,可选的参数还有Aerosol Retrieval 和Initial Visibility。这两个参数对最后的结果又相当重要的影像,因此最好能调查到当地的Initial Visibility。此外,AERONET 在全世界各地有测定AOD(Atmospheric Optical Depth)的站点,可以查询AOD 以后转换为消光系数,通过消光系数估算能见度,此步骤比较繁琐,在此不予详述。如果采用Aerosol Retrieval 中的K-T算法计算Visibility,且能够计算出结果的话,则采用K-T 算法的能见度,否 则采用Initial Visibility 所指定的能见度。 d) 关于Aerosol Retrieval。如果选择了下拉菜单中的K-T method,那么需要在Multispectral Settings 中设定参数,在Assign Default Values Based on Retrieval Conditions 中选择Over-land Retrieval Standard (660:2100nm)即可。根据不同的研究区可以设定不同的模式。其他设定可以不改变。Apply即可。
个人分类: 他山之石|6470 次阅读|0 个评论
谈谈大地电磁测深中一些基本概念(二):静位移的校正
热度 2 陈小斌 2010-8-17 01:52
自去年写了介绍大地电磁视电阻率静位移现象的博文以后,至今已经一年多了,一直不愿意写(二)。事情太多太杂。最近好像看到说中国科学家只有三分之一的时间做科研,好多人还惊讶不已。事实是,对于那些大牛们而言,有三分之一的时间已经是很不错的了! 言归正传。本来我想接着介绍大地电磁局部畸变和张量分解的问题,但最近一年多来接触到的情况,让我想继续介绍一下有关静位移校正的问题。正如我前面提到过,静位移问题是大地电磁中最难克服的问题。王家映教授去年、今年就这个问题做了几次报告,也持这个观点。王老师归纳了当今几乎所有用于静位移校正的方法,并指出都有其缺陷。但王老师最后将希望寄托在张量分解技术上面,这是个问题。我后续的文章(如果还写的话)将详细介绍张量分解与静位移之间的关系。这里需要强调:正是因为静位移因子无法确定,才导致张量分解技术得到的区域阻抗张量无法直接使用,大大限制了张量分解技术的实用价值。张量分解技术要想成为一项强有力的技术手段,还有待静位移校正技术的解决和支持。 当前,在工程电磁法勘探中,采用空间滤波、首枝重合等技术手段对观测的大地电磁视电阻率数据进行批量式静位移校正,似已成为常用的资料处理手段,这是非常危险的。确实,经过批量式的滤波处理后,视电阻率剖面图看起来更加舒服、更能符合我们的视觉需要,但这不一定是地球的真实响应。不分青红皂白地滤波处理很可能会去除有用的区域构造响应,而产生假的构造。 为什么这么说? 因为静位移的确认至少满足两个条件:第一、视电阻率发生了整体性偏移;第二、两个方向的相位(Pxy、Pyx)不受影响,是重合的。当前校正方法中,无论是采用首枝重合法,还是空间滤波法,只要是批量式处理,都只利用了视电阻率整体性偏移这一个条件,而没有用到相位不变的条件。那么,是否存在这样的区域构造模型,也可使得视电阻率发生了很明显的偏移?如果存在,则对这类资料进行滤波处理,就会破坏我们需要的区域构造的响应。 这样的模型是否存在呢?答案是肯定的,不仅存在,而且很常见。 切穿到地表的直立断层,由于其两侧电阻率的突变,就会产生这样的响应。图1所示为三个切穿到地表的直立断层某一频率的剖面响应,可见TM(蓝点)模式的视电阻率曲线发生了很大跳跃性变化,这是断层的响应,不是小异常体的静位移效应。 图2为该二维理论模型中某断层附近一测点的视电阻率和相位随频率变化的曲线。可见视电阻率和相位首支都不重合,都产生了明显的偏移现象,但同时,相位曲线也产生了分离现象,清楚地表明这不是静位移效应。 如果采用空间滤波的方法将图1、图2所示的视电阻率曲线进行光滑,则反演结果中不可能再是三条直立断层了,因为其所用的反演数据已经不再是三条直立断层的响应数据了。 以上模型使得视电阻率和相位都产生了偏移或分离现象,从而可以比较容易区分静位移和区域断层构造响应。有没有这样的区域构造响应,其高频处视电阻率表现为偏移现象,而相位重合的?我只能说不知道。不过就我目前所做过的理论模型、所阅读的文献看,还没有遇到这种情况。实测资料中比较难于处理的是数据受到了噪声的影响,有时候很难判断相位曲线究竟是分离的还是合拢的。 无论如何,对于观测数据进行批量式静位移校正,一定要慎之又慎! 图1 单频率二维模型的正演响应沿剖面的变化。二维模型中,从左到右共有三条直立断层切穿地表。蓝色为TM极化模式数据,红色为TE极化模式数据;上为视电阻率,下为相位曲线。 图2 切穿地表的断层附近,首支分开的视电阻率和相位曲线。视电阻率的静位移发生在低频段10S以下部分,高频部分的曲线分离并不是真的静位移引起的。
个人分类: 专业探讨|10828 次阅读|6 个评论

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-6-16 17:38

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部