[Hector学习笔记]GNSS时间序列处理软件Hector使用备忘（批处理脚本）

由Machiel Bos 和Rui Fernandes共同研发的Hector软件可以利用时间相关噪声来估计时间序列线性趋势，是与优秀的GNSS时间序列处理软件之一。
软件安装教程：Hector软件安装[附加软件环境变量设置]
Hector软件用户手册中包含了教程，示例文件等等，并且有非常详细的说明。在教程中，大多采用“设置驱动文件+在终端中逐行输入代码”进行计算的方式进行解算，这适用于新手，或者是需要处理的站点较少的情况。对于超过5个站点的时间序列的分析，如果还用逐个修改驱动文件和终端逐个输入的方式，未免太浪费时间。在学习该软件的过程中，逐步总结出批量处理站点时间序列的脚本，以备日后查询和学习。根据教程中的处理流程，分别从：剔除粗差、估计趋势项和估计功率谱密度三个方面进行总结，以上三部分具体的参数设置在用户手册中均有说明，此处不再赘述。

剔除粗差
输入文件格式：[site].neu
输出文件格式：[site].outliers.[direction].mom（此处文件名称“outliers”可自定义，[direction]用来表示方向，mom是固定的）
说明：由于Hector处理测站各个方向的时间序列是分开的，所以每处理一个方向所生成的文件名对方向的标注是非常有必要的；另外上述输出格式中文件名称“outliers”部分可以自己修改成任何其他名称，只是为了后面更方便的辨认；最后，重要的是Hector后面的步骤的输入文件均是mom格式，所以最后一定要以mom格式结尾。
下面是脚本：

#!/bin/bash
#remove outliers use Hector
#列出所有测站(由于测站较多，省略了大部分)
sites=(AHAQ AHBB BJGB BJSH BJYQ... )
for ((i=0;i<${#sites[*]};i++))
dosite=${sites[$i]}n_raw=$site.neuc_raw=`grep "DataFile" removeoutliers.ctl | awk '{print $2}'`
#   终端显示此时正在处理的文件，以便修改
#   echo "$c_raw $n_raw"
#   将驱动控制文件中文件名修改成要处理的文件名sed -i "s/$c_raw/$n_raw/g" removeoutliers.ctl
#   列出所有方向directions=(North East Up)for ((j=0;j<${#directions[*]};j++))dodirection=${directions[$j]}c_dir=`grep "component" removeoutliers.ctl | awk '{print $2}'`
#       将驱动控制文件中坐标分量改成要处理的分量，并通过循环实现三个坐标分量的替换sed -i "s/$c_dir/$direction/g" removeoutliers.ctlout_file=$site.outliers.$direction.mom
#       修改控制文件中的输出文件，此时输出文件需要区分开站的各个方向c_out_file=`grep "OutputFile" removeoutliers.ctl | awk '{print $2}'`echo "$c_out_file"sed -i "s/$c_out_file/$out_file/g" removeoutliers.ctlremoveoutliers removeoutliers.ctldone
done

【注】脚本中测站名的制作方法将单独介绍。如果要处理的测站有好几百个，列出测站名不能采用手动输入的方式，会有更好的办法进行制作。
更新：制作方法在这：GNSS时间序列分析笔记Ubuntu下制作站点列表
上述脚本可以放在与要处理的文件同一个文件夹下，当然也可以放在其他文件夹下，在removeoutliers.ctl文件中有对文件位置做详细的设置，比如我的是这样：（此时脚本与含有待处理文件的deformation_make在同一文件夹下，这样的逻辑关系可以自行体会）

经过处理后生成的文件与脚本在同一个文件夹下，生成的文件：

是Hector要求的mom格式。接下来进行趋势项估计。
2. 趋势项估计
输入文件格式：[site].outliers.[direction].mom
输出文件格式：[site].[direction].[noise].trend.mom
说明：在估计趋势项和周期项时，可以加入各种噪声以及噪声组合，此时输出文件需要对不同噪声进行区分。且输出问价后缀仍需是mom格式，以便后期的功率谱估计。
趋势项估计的脚本如下：

#!/bin/bash
#trend estimation by Hector
#列出所有测站(由于测站较多，省略了大部分)
sites=(AHAQ AHBB BJGB BJSH BJYQ... )
for ((i=0;i<${#sites[*]};i++))
dosite=${sites[$i]}directions=(North East Up)for ((j=0;j<${#directions[*]};j++))dodirection=${directions[$j]}m_names=(wh wh_fn wh_rw wh_pl wh_fn_rw)models=("White" "White FlickerGGM" "White RandomWalkGGM" "White Powerlaw" "White FlickerGGM RandomWalkGGM")for((k=0;k<${#m_names[*]};k++))dom_name=${m_names[$k]}model=${models[$k]}out_file=$site.$direction.$m_name.trend.momc_out_file=`grep "OutputFile" estimatetrend1.ctl | awk '{print $2}'`sed -i "s/$c_out_file/$out_file/g" estimatetrend.ctlecho "Outputting $out_file"c_model=`grep "NoiseModels" estimatetrend1.ctl | awk '{print $2}'`sed -i "s/$c_model/$model/g" estimatetrend.ctlre_file=$site.$direction.$m_nameestimatetrend estimatetrend1.ctl > $re_filedonedone
done

以上为完整的带坐标分量和噪声组合循环的脚本，我在解算是仅用到了测站的up方向，且仅考虑白噪声影响。如有类似这样的需要，可以将相应的循环部分取消，令其直接等于某个分量、噪声模型即可。此时解算得到的文件有mom文件和.wh文件（只考虑白噪声）。文件分别长这样：

这里的.wh文件是estimatetrend的log文件，里面记载了解算的结果，包括估计噪声的最大似然估计值、周年与半周年项等等信息，此时如果我们需要提取其中的信息，仍然可以采用脚本的方法。下面先以一个log文件为例：


************************************estimatetrend, version 1.6
************************************
0) White
generator type: mt19937
seed = 0
first value = 3176401122
Data format: MJD, Observations, Model
Filename              : ./AHAQ_GPS.outliers.Up.mom
Number of observations: 5569
Percentage of gaps    : 0.8799----------------AmmarGrag
----------------
No Polynomial degree set, using offset + linear trend
No extra periodic signal is included.Number of CPU's used (threads) = 1performing Ordinary Least-SquaresWhite:
fraction  = 1.00000
sigma     = 2.88355 mm
No noise parameters to showLikelihood value
--------------------
min log(L)=-13678.350
k         =4 + 0 + 1 = 5
AIC       =27366.701
BIC       =27399.781STD of the driving noise: 2.884
bias : -0.085 +/- 0.039 mm (at 2009/11/13, 12:0:0.000)
trend: -0.343 +/- 0.009 mm/year
cos yearly : -2.982 +/- 0.055 mm
sin yearly : 1.843 +/- 0.055 mm
Amp yearly : 3.506 +/- 0.055 mm
Pha yearly : 148.276 degrees
--> AHAQ_GPS.Up.wh.trend.mom
Total computing time: 8.00000 sec

此时利用下面的脚本完成振幅相位的提取（我在解算时仅估计了周年项，若想估计半周年项，需要在ctl文件中进行相关设置）

#!/bin/bash
#grep amplitude and phase from log file sites=(AHAQ AHBB BJGB BJSH BJYQ CHUN CQCS CQWZ DLHA DXIN FJPT FJWY FJXP GDSG GDST GDZH GDZJ GSAX GSDH GSDX GSGL GSGT GSJN GSJT GSJY GSLX GSLZ GSMA GSML GSMQ GSMX GSPL GSQS GSTS GSWD GUAN GXBH GXBS GXGL GXHC GXNN GXWZ GZFG GZGY GZSC HAHB HAJY HAQS HBES HBJM HBXF HBZG HECC HECD HECX HELQ HELY HETS HEYY HEZJ HIHK HISY HLAR HLFY HLHG HLMH HLWD HNLY HNMY HRBN JIXN JLCB JLCL JLYJ JSLS JSLY JSNT JSYC JXHK JXJA KMIN LALB LALX LHAS LNDD LNHL LNJZ LNSY LNYK LUZH NMAG NMAL NMAY NMAZ NMBT NMDW NMEJ NMEL NMER NMTK NMWH NMWJ NMWL NMWT NMZL NXHY NXYC NXZW QHBM QHDL QHGC QHGE QHLH QHMD QHME QHMQ QHMY QHQL QHTR QHTT QHYS QION SCBZ SCDF SCGZ SCJL SCJU SCLH SCLT SCMB SCML SCMN SCMX SCNC SCNN SCPZ SCSM SCSN SCSP SCTQ SCXC SCXD SCXJ SCYX SCYY SDCY SDJX SDLY SDQD SDRC SDYT SDZB SHA2 SNAK SNHY SNMX SNTB SNXY SNYA SUIY SXCZ SXDT SXGX SXKL SXLF SXLQ SXTY SXXX SXYC TAIN TJBD TJBH TJWQ WUSH XIAG XIAM XJAL XJBC XJBE XJBL XJBY XJDS XJFY XJHT XJJJ XJKC XJKE XJKL XJML XJQH XJQM XJRQ XJSH XJSS XJTC XJTZ XJWL XJWQ XJWU XJXY XJYC XJYN XJYT XJZS XNIN XZAR XZBG XZCD XZCY XZDX XZGE XZGZ XZNM XZRK XZRT XZYD XZZB XZZF YANC YNCX YNDC YNGM YNHZ YNJD YNJP YNLA YNLC YNLJ YNMH YNMJ YNML YNMZ YNRL YNSD YNSM YNTC YNTH YNWS YNXP YNYA YNYL YNYM YNYS YNZD ZHNZ ZJJD ZJWZ ZJZS )
#for ((i=0;i<${#sites[*]};i++))
for ((i=0;i<2;i++))
dosite=${sites[$i]}_GPSfile_name=$site.Up.whecho "new doing $file_name..."cos=`grep "cos yearly" $file_name | awk '{print $3 $4 $5}'`sin=`grep "sin yearly" $file_name | awk '{print $3 $4 $5}'`amp=`grep "Amp yearly" $file_name | awk '{print $3 $4 $5}'`pha=`grep "Pha yearly" $file_name | awk '{print $4}'`printf "%-1s %-1s %-1s %-1s\n" $file_name $cos $sin $amp $pha>> amp.dat
done

提取出来的结果文件如下：

可以对以上文件在excel中进行下一步分析。
3. 功率谱估计
按照之前的方法先用脚本运行estimatespectrum，脚本如下：（与上面脚本所用到的测站不一样）

#!/bin/bashsites=(AQSS AQWJ AQYX CZLA CZMG CZTC FYFN FYJS HBSX HSQM LAHQ LAHS MASM SZDS SZLB WHWH XCJD XCLX XCNG)
for ((i=0;i<${#sites[*]};i++))
dosite=${sites[$i]}directions=(North East Up)dirs=(N E U)for ((j=0;j<${#directions[*]};j++))dodirection=${directions[$j]}dir=${dirs[$j]}
#       m_names=(wh wh_fn wh_rw wh_pl wh_fn_rw)
#       for((k=0;k<${#m_names[*]};k++))
#       do
#           m_name=${m_names[$k]}file_name=$site.$direction.wh_fn.trend.momout_file=$site.$dir.wh_fn.spectrumlog_file=$site.$direction.wh_fn.logc_file=`grep "DataFile" estimatespectrum.ctl | awk '{print $2}'`sed -i "s/$c_file/$file_name/g" estimatespectrum.ctlc_out_file=`grep "OutputFile" estimatespectrum.ctl | awk '{print $2}'`sed -i "s/$c_out_file/$out_file/g" estimatespectrum.ctlestimatespectrum estimatespectrum.ctl > $log_file
#       donedone
done

生成的文件有：

[Hector学习笔记]GNSS时间序列处理软件Hector使用备忘（批处理脚本）相关推荐

软件调试学习笔记（五）—— 软件断点内存断点
软件调试学习笔记(五)-- 软件断点&内存断点调试的本质软件断点软件断点的执行流程分析INT 3执行流程实验:处理软件断点内存断点内存断点的执行流程实验:处理内存断点调试的本 ...
汽车电子学习笔记—AutoSAR之基础软件层（BSW）
汽车电子学习笔记-AutoSAR之基础软件层(BSW) - 1.概述如之前autosar概述笔记中说明,BSW按照层级结构可以分为服务层.ECU抽象层.硬件抽象层(MCAL)和复杂驱动层(CDD). ...
学习笔记之编程达到一个高的境界就是自制脚本语言（图）
学习笔记之编程达到一个高的境界就是自制脚本语言(图) 编程达到一个高的境界就是自制脚本语言,通过这可以精通编程里面的高深的技术,如编译原理.语言处理器.编译器与解释器,这些都是代表一个程序员实力的技术 ...
苹果手机什么日程安排提醒软件可以事件备忘和随时提醒？
苹果手机上什么日程安排提醒软件可以事件备忘和随时提醒? iPhone手机上有备忘录用来事件备忘,有提醒事项可以日程提醒,iPhone手机上想用一款集备忘与提醒于一体的桌面日程安排提醒软件,可以添加云便 ...
影像组学视频学习笔记(41)-如何使用软件提取组学特征、Li‘s have a solution and plan.
作者:北欧森林链接:https://www.jianshu.com/p/72186eb3e395 来源:简书,已获授权转载本笔记来源于B站Up主: 有Li 的影像组学系列教学视频本节(41)主要 ...
影像组学视频学习笔记(18)-使用MRIcroGL软件格式转换、勾画ROI、Li‘s have a solution and plan.
本笔记来源于B站Up主: 有Li 的影像组学系列教学视频本节(18)主要介绍: 使用MRIcroGL软件进行格式转换.勾画ROI并保存为mask文件 MRIcroGL 是一款免费.开源的轻量级软件: ...
freeRtos学习笔记（6）软件定时器
freeRtos学习笔记 freeRtos软件定时器软件定时器需要注意事项软件定时器的精度基于时钟节拍,例如系统时钟节拍为10ms, 软件定时器定时时间必须是10ms的整数倍,因此软件定时器一般用 ...
学习笔记 51单片机通用软件延时方法
对于STC51单片机来说,延时函数,想必都不陌生.而用的最多的延时基本都是通过软件方法实现的,但由于STC51不同系列的芯片所采用的指令集不同,各指令执行所用机器周期不同.例如STC12Cx的一个振荡 ...
学习笔记整理：Photoshop软件应用-图层混合与样式
以下内容为个人的学习笔记整理,如有错误,请指出,谢谢~ 一.混合模式的种类(最少要两个图层以上) (1)溶解:用于图层中的图像出现透明像素时,根据图像中透明像素的量显示出颗粒化效果. (2)正片叠底: ...

[Hector学习笔记]GNSS时间序列处理软件Hector使用备忘（批处理脚本）

[Hector学习笔记]GNSS时间序列处理软件Hector使用备忘（批处理脚本）相关推荐

最新文章

热门文章