博文

xvg曲线中心移动平均平滑脚本

已有 2951 次阅读 2017-5-11 08:21 |系统分类:科研笔记

2017年05月10日 18:25:47

在查看GROMACS分析程序生成的xvg文件时, 有时候曲线不够平滑, 上下波动很厉害. 为美观或其他目的有时候需要对曲线进行平滑. 常见的绘图软件一般都支持曲线平滑, 可选方法也很多. 其中最简单的是移动平均算法(moving average). 但简单的移动平均算法平滑出来的曲线存在滞后性, 对于有峰和谷的曲线不适合. 一个解决方法是采用基于中心点的移动平均算法. 简单地说, 每个点的平滑值是其自身及前后若干点数值的平均值. 下面这段脚本专门用于xvg文件的平滑, 供需要的人参考.

gmx_smxvg.bsh
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46	usage=">>>>>>>>>>>>>>>>gmx_smxvg<<<<<<<<<<<<<<<<>>>>>>>>>>>>>>>>JicunLi<<<<<<<<<<<<<<<<>>>>>>>>>>>>2017-05-1014:44:21<<<<<<<<<<<<>>Usage:gmx_smxvg{-n5}{-col2}<File.xvg...>>>-n:NumberofPointstoRunSmoothingAverage,ShoulebeODD>>-c:ColumnNumberoftheDatatoSmooth"[[$#-lt1]]&&{echo"$usage";exit;}pnt=5;col=2opt=($@);N=${#opt[@]}for((i=0;i<N;i++));doarg=${opt[$i]};j=$((i+1))[[$arg=~-n]]&&{pnt=${opt[$j]};opt[$i]="";opt[$j]="";}[[$arg=~-c]]&&{col=${opt[$j]};opt[$i]="";opt[$j]="";}done[[$((pnt%2))-eq0]]&&{echo'ERROR!!!NumberofPointsshouldbeODD!';exit;}forfilein${opt[@]};doawk-vpnt=$pnt-vcol=$col'\/^[#@]/{print;next}/^[^#@]/{n++;j=n%pntif(n<pnt){sma+=$colif(n<(pnt+1)/2){for(i=1;i<=col;i++)printf"%12.6f",$iprintf"%12.6f",sma/nfor(i=col+1;i<=NF;i++)printf"%12.6f",$iprint""}}else{sma+=$col-f[j]now=(n-(pnt-1)/2)%pntsplit(txt[now],dat)printf"%12.6f",dat[1]for(i=2;i<=col;i++)printf"%12.6f",dat[i]printf"%12.6f",sma/pntfor(i=col+1;i<=NF;i++)printf"%12.6f",dat[i]print""}f[j]=$col;txt[j]=$0}'$file>${file%.xvg}~sm.xvgdone

gmx_smxvg.bsh

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46

usage=">>>>>>>>>>>>>>>>gmx_smxvg<<<<<<<<<<<<<<<<>>>>>>>>>>>>>>>>JicunLi<<<<<<<<<<<<<<<<>>>>>>>>>>>>2017-05-1014:44:21<<<<<<<<<<<<>>Usage:gmx_smxvg{-n5}{-col2}<File.xvg...>>>-n:NumberofPointstoRunSmoothingAverage,ShoulebeODD>>-c:ColumnNumberoftheDatatoSmooth"[[$#-lt1]]&&{echo"$usage";exit;}pnt=5;col=2opt=($@);N=${#opt[@]}for((i=0;i<N;i++));doarg=${opt[$i]};j=$((i+1))[[$arg=~-n]]&&{pnt=${opt[$j]};opt[$i]="";opt[$j]="";}[[$arg=~-c]]&&{col=${opt[$j]};opt[$i]="";opt[$j]="";}done[[$((pnt%2))-eq0]]&&{echo'ERROR!!!NumberofPointsshouldbeODD!';exit;}forfilein${opt[@]};doawk-vpnt=$pnt-vcol=$col'\/^[#@]/{print;next}/^[^#@]/{n++;j=n%pntif(n<pnt){sma+=$colif(n<(pnt+1)/2){for(i=1;i<=col;i++)printf"%12.6f",$iprintf"%12.6f",sma/nfor(i=col+1;i<=NF;i++)printf"%12.6f",$iprint""}}else{sma+=$col-f[j]now=(n-(pnt-1)/2)%pntsplit(txt[now],dat)printf"%12.6f",dat[1]for(i=2;i<=col;i++)printf"%12.6f",dat[i]printf"%12.6f",sma/pntfor(i=col+1;i<=NF;i++)printf"%12.6f",dat[i]print""}f[j]=$col;txt[j]=$0}'$file>${file%.xvg}~sm.xvgdone