CN105631133B - 基于频域平滑算法的chirp扫频信号处理方法 - Google Patents
基于频域平滑算法的chirp扫频信号处理方法 Download PDFInfo
- Publication number
- CN105631133B CN105631133B CN201511009339.9A CN201511009339A CN105631133B CN 105631133 B CN105631133 B CN 105631133B CN 201511009339 A CN201511009339 A CN 201511009339A CN 105631133 B CN105631133 B CN 105631133B
- Authority
- CN
- China
- Prior art keywords
- frequency
- signal
- chirp
- swept
- frequency domain
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
- G06F17/141—Discrete Fourier transforms
- G06F17/142—Fast Fourier transforms, e.g. using a Cooley-Tukey type algorithm
Landscapes
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Theoretical Computer Science (AREA)
- Computational Mathematics (AREA)
- Data Mining & Analysis (AREA)
- General Engineering & Computer Science (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Mathematical Analysis (AREA)
- Computer Hardware Design (AREA)
- Discrete Mathematics (AREA)
- Evolutionary Computation (AREA)
- Algebra (AREA)
- Geometry (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Complex Calculations (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
基于频域平滑算法的chirp扫频信号处理方法,涉及chirp扫频信号处理技术。是为了适应对chirp扫频信号处理的需求。本发明中,首先对原始输入信号进行快速傅里叶变换;并对频域信号进行截短;然后通过加窗取平均和三次样条插值后输出。本发明适用于chirp扫频信号处理。
Description
技术领域
本发明涉及chirp扫频信号处理技术。
背景技术
对于一个线性时不变的因果系统,初始条件为零时,其零输入响应为零。在时域,如果用e(t)表示激励,r(t)表示响应,h(t)表示单位冲激响应,则激励和响应的关系可写为:
r(t)=e(t)*h(t) (1)
由时域卷积性质,可得到其频域的傅里叶变换:
这是频响函数的原始定义。
传统的扫频方法采用正弦信号做激励,该激励的幅值不变,频率随时间的变化以某一固定步长△f增加。假设从t0到t0+△t之间信号源发出起始频率为f0的正弦,则在t0+△t到t0+2△t之间信号源发出的频率为f0+△f的正弦波,依次类推。
假设扫描的时间为T=N△t,则扫描的频带范围为f0~f0+(N-1)△f。通过设定起始频率f0、扫描时间T和频率变化步长△f,就可以在指定的频率范围内测定所研究电网络的频率特性。
根据傅里叶变换的公式,如果使用一个带宽信号作为激励,那么仅通过一次扫描就能得到信号在带宽范围内的频率特性,而无需向传统的扫频方法那样要逐点分布测量,从而加快了测量的速度。
采用无幅度调制的线性调频脉冲作为激励信号(在雷达技术领域,该信号又称chirp信号),其表达式为:
在测量中实际用到的是chirp信号的实部:
cos(βt2+ω0t) (4)
由图1所示的波形可见,该信号是基于预先的频率不断随时间增长的信号。
当chirp信号进入到实际系统中,由于实际器件中的噪声等因素的影响,会对频率响应产生非常大的影响。由图2和3所示的波形所示。
由图像,可知信号的大概频谱特征,但是由于未知因素的影响,还无法从原始数据中得到可以直接使用的平滑信息。因此,对chirp信号(切普信号)进行扫频是非常重要的。
发明内容
本发明是为了适应对chirp扫频信号处理的需求,从而提供一种基于频域平滑算法的chirp扫频信号处理方法。
基于频域平滑算法的chirp扫频信号处理方法,它由以下步骤实现:
步骤一、对原始输入信号进行快速傅里叶变换,具体为:
对于离散的输入信号e(n)和输出信号r(n),首先测量出该两列信号的长度,然后进行快速傅里叶变换,变换到频域后的信号分别表示为E(n)和R(n);
步骤二、对频域信号进行截短,具体为:
对于两列离散信号E(n)和R(n),根据公式:
获得预设的起始频率f1和预设的截止频率f2对应的数据点数n1和n2;
式中:N是输入数据的点数,fs是采样频率;
若n1或n2不是整数时,将n1或n2向下取整;
步骤三、对于数据H(ni),ni∈(n1,n2),根据公式:
对截短后的信号进行平均处理;
式中,nj∈(n1,n2);
步骤四、对步骤三得到的平均处理结果进行插值处理,具体为:
在(n1,n2)之间,取出C个值作为插值对像,C为正整数;对该C个值进行三次样条插值后,作为最终处理结果输出。
本发明充分适应了对chirp扫频信号处理的需求。
附图说明
图1是背景技术中的chirp信号仿真示意图;
图2是背景技术中的chirp信号幅频特性仿真示意图;
图3是背景技术中的chirp信号相频特性仿真示意图;
图4是本发明具体实施方式中信号的幅频特性仿真示意图;
图5是本发明具体实施方式中信号的相频特性仿真示意图;
图6是本发明的信号处理流程示意图;
图7是本发明具体实施方式中性能分析获得的0-18Hz幅频特性仿真示意图;其中曲线71为原始chirp信号;曲线72为经本发明处理过的信号;
图8是本发明具体实施方式中性能分析获得的0-18Hz相频特性仿真示意图;其中曲线81为原始chirp信号;曲线82为经本发明处理过的信号;
图9是本发明具体实施方式中性能分析获得的0.1-16Hz幅频特性仿真示意图;其中曲线91为原始chirp信号;曲线92为经本发明处理过的信号;
图10是本发明具体实施方式中性能分析获得的0.1-16Hz相频特性仿真示意图;其中曲线101为原始chirp信号;曲线102为经本发明处理过的信号;
图11是本发明具体实施方式中性能分析获得的0.1-26Hz幅频特性仿真示意图;其中曲线111为原始chirp信号;曲线112为经本发明处理过的信号;
图12是本发明具体实施方式中性能分析获得的0.1-26Hz相频特性仿真示意图;其中曲线121为原始chirp信号;曲线122为经本发明处理过的信号;
具体实施方式
具体实施方式一、结合图6说明本具体实施方式,基于频域平滑算法的chirp扫频信号处理方法,它由以下步骤实现:
步骤一、对原始输入信号进行快速傅里叶变换,具体为:
对于离散的输入信号e(n)和输出信号r(n),首先测量出该两列信号的长度,然后进行快速傅里叶变换,变换到频域后的信号分别表示为E(n)和R(n);
步骤二、对频域信号进行截短,具体为:
对于两列离散信号E(n)和R(n),根据公式:
获得预设的起始频率f1和预设的截止频率f2对应的数据点数n1和n2;
式中:N是输入数据的点数,fs是采样频率;
若n1或n2不是整数时,将n1或n2向下取整;
步骤三、对于数据H(ni),ni∈(n1,n2),根据公式:
对截短后的信号进行平均处理;
式中,nj∈(n1,n2);
步骤四、对步骤三得到的平均处理结果进行插值处理,具体为:
在(n1,n2)之间,取出10个值作为插值对像;对该10个值进行三次样条插值后,作为最终处理结果输出。
原理:对输入数据进行快速傅里叶变换FFT,根据输入数据点数N和采样频率fs,可以计算采样分辨率:
根据预设的起始频率f1和预设的截止频率f2,可以得到分别对应的数据点数:
由此,对f1到f2之间的信息进行频域平滑就等价于对n1到n2之间的数据点信息进行频域平滑。
定义频响为H(ω),对应于FFT变换后的频响为H(n)。以下为对频响H(n1)到H(n2)之间的数据点进行处理。
对于ni∈(n1,n2),对H(ni)数据的处理方式如下:
如图4和图5所示的波形可见,加窗取均值的做法能够很好的于频响的趋势相吻合,但是无法实现很好的平滑效果,究其原因是抽样点过密造成的。
对于这种抽样点过于密集所导致的不平滑现象,通常采用的手段是用插值的方法对已处理数据进行再处理。
现给出样条函数的定义:对于给定区间[a,b]的一个划分:a=x0<x1<…<xn=b,如果函数S(x)同时满足:
1)、在每个区间[xi,xi+1],i=0,1,…,n-1上S(x)是m次多项式;
2)、S(x)在给定区间[a,b]上具有m-1阶连续导数;
则称S(x)为关于上述划分的m次样条函数。
由此,可以推导出关于三次样条插值的定义:
已知函数y=f(x)在给定区间[a,b]上的n+1个节点a=x0<x1<…<xn=b上的值yj=f(xj),j=0,1,…,n,求插值函数S(x)使其同时满足:
1)、S(xj)=yj,j=0,1,…,n;
2)、在每个区间[xj,xj+1]上S(x)是三次项式,记为Sj(x),j=0,1,…,n-1;
3)、S(x)在给定区间[a,b]上二阶连续可微;
则S(x)称为f(x)的三次样条插值函数,它通过上述的给定点,为二阶连续可导的分段三次多项式函数。
根据上述三次样条函数插值的叙述,采用对处理区间采用10点插值。
即对向下取整,得到的数值为插值间隔,从而去除差值过程中的已知点(xj,yj)。通过插值,得到一条平滑的曲线。
以下以具体的仿真试验验证本发明的效果:
结合图7至图12,采用三组数据对设计好的系统进行测试,其中曲线72、82、92、102、112、122表示chirp信号的原始数据,而曲线71、81、91、101、111、121表示处理过的数据。
从图7至图12可以看到,整体趋势是可以准确描述的,但是在图9和10中,可以显著的发现其中的不准确性,这主要是插值方式的选取所造成的。
经过测试发现,由于数据之间存在误差,可以采用自适应的方法去自我调节固定的窗宽度5和插值已知点数10;
样条插值这种方法可以用更高级的,适应性更强的方法代替,比如:采用非监督学习的方法,可以有效的通过自我学习得到理想的结果。
对数据本身可以采用一些方法进行预处理,比如:对一些幅度变化较大的部分进行幅值的截短,再通过这些相对平滑的幅值画出包络形状。通过两条包络的信息可以有效规避单纯用样条插值得到的不好结果。
Claims (4)
1.基于频域平滑算法的chirp扫频信号处理方法,其特征是:它由以下步骤实现:
步骤一、对原始输入信号进行快速傅里叶变换,具体为:
对于离散的输入信号e(n)和输出信号r(n),首先测量出该两列信号的长度,然后进行快速傅里叶变换,变换到频域后的信号分别表示为E(n)和R(n);
步骤二、对频域信号进行截短,具体为:
对于两列离散信号E(n)和R(n),根据公式:
获得预设的起始频率f1和预设的截止频率f2对应的数据点数n1和n2;
式中:N是输入数据的点数,fs是采样频率;
若n1或n2不是整数时,将n1或n2向下取整;
步骤三、对于数据H(ni),ni∈(n1,n2),根据公式:
对截短后的信号进行平均处理;
式中,nj∈(n1,n2);
步骤四、对步骤三得到的平均处理结果进行插值处理,具体为:
在(n1,n2)之间,取出C个值作为插值对像,C为正整数;对该C个值进行m次样条插值后,作为最终处理结果输出。
2.根据权利要求1所述的基于频域平滑算法的chirp扫频信号处理方法,其特征在于C=10。
3.根据权利要求1所述的基于频域平滑算法的chirp扫频信号处理方法,其特征在于步骤四中所述的m次样条插值的含义为:
对于给定区间[a,b]的一个划分:a=x0<x1<…<xn=b,如果函数S(x)同时满足:
1)、在每个区间[xi,xi+1],i=0,1,…,n-1上S(x)是m次多项式;
2)、S(x)在给定区间[a,b]上具有m-1阶连续导数;
则称S(x)为关于上述划分的m次样条函数。
4.根据权利要求3所述的基于频域平滑算法的chirp扫频信号处理方法,其特征在于步骤四中所述的m次样条插值,当m=3时的含义为:
已知函数y=f(x)在给定区间[a,b]上的n+1个节点a=x0<x1<…<xn=b上的值yj=f(xj),j=0,1,…,n,求插值函数S(x)使其同时满足:
1)、S(xj)=yj,j=0,1,…,n;
2)、在每个区间[xj,xj+1]上S(x)是三次项式,记为Sj(x),j=0,1,…,n-1;
3)、S(x)在给定区间[a,b]上二阶连续可微;
则S(x)称为f(x)的三次样条插值函数。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201511009339.9A CN105631133B (zh) | 2015-12-28 | 2015-12-28 | 基于频域平滑算法的chirp扫频信号处理方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201511009339.9A CN105631133B (zh) | 2015-12-28 | 2015-12-28 | 基于频域平滑算法的chirp扫频信号处理方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN105631133A CN105631133A (zh) | 2016-06-01 |
CN105631133B true CN105631133B (zh) | 2019-09-13 |
Family
ID=56046062
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201511009339.9A Active CN105631133B (zh) | 2015-12-28 | 2015-12-28 | 基于频域平滑算法的chirp扫频信号处理方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN105631133B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106844292A (zh) * | 2017-01-12 | 2017-06-13 | 天津大学 | 基于快速傅里叶变换的室内空气数据异常值平滑方法 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101753497A (zh) * | 2009-11-27 | 2010-06-23 | 哈尔滨工业大学 | 基于环境感知的切普信号与余弦信号联合的信号调制和解调方法及其信号发射和接收方法 |
CN102724155A (zh) * | 2012-05-17 | 2012-10-10 | 哈尔滨工程大学 | 基于分数傅里叶变换的高频域能量集中度同步方法 |
CN102778674A (zh) * | 2012-05-25 | 2012-11-14 | 安徽理工大学 | 非均匀采样的Chirp脉冲时延估计方法 |
CN103888404A (zh) * | 2014-04-10 | 2014-06-25 | 厦门大学 | 一种基于频谱搬移的全频谱载波调制方法 |
-
2015
- 2015-12-28 CN CN201511009339.9A patent/CN105631133B/zh active Active
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101753497A (zh) * | 2009-11-27 | 2010-06-23 | 哈尔滨工业大学 | 基于环境感知的切普信号与余弦信号联合的信号调制和解调方法及其信号发射和接收方法 |
CN102724155A (zh) * | 2012-05-17 | 2012-10-10 | 哈尔滨工程大学 | 基于分数傅里叶变换的高频域能量集中度同步方法 |
CN102778674A (zh) * | 2012-05-25 | 2012-11-14 | 安徽理工大学 | 非均匀采样的Chirp脉冲时延估计方法 |
CN103888404A (zh) * | 2014-04-10 | 2014-06-25 | 厦门大学 | 一种基于频谱搬移的全频谱载波调制方法 |
Also Published As
Publication number | Publication date |
---|---|
CN105631133A (zh) | 2016-06-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105229431A (zh) | 具有改进的距离确定的料位测量 | |
CN109061302A (zh) | 一种基于EEMD和Hilbert变换的并网风力发电机组谐波测量系统 | |
US9960862B2 (en) | Method and device for detecting standing-wave ratio | |
KR101687658B1 (ko) | 처프-지 역변환 방법 및 시스템 | |
CN109359633B (zh) | 基于希尔伯特-黄变换和小波脊线的信号联合分类方法 | |
CN111912521A (zh) | 一种非平稳信号的频率检测方法和存储介质 | |
CN105631133B (zh) | 基于频域平滑算法的chirp扫频信号处理方法 | |
CN107991591A (zh) | 一种基于Kaiser窗FFT单峰插值修正的图像融合方法 | |
CN107209259B (zh) | 用于测距的方法和装置 | |
CN106471748A (zh) | 传输介质中的阻抗的估计 | |
CN104392141A (zh) | 电力系统低频振荡扰动源定位方法和装置 | |
CN109117816A (zh) | 基于六阶样条插值小波的信号奇异点检测方法 | |
CN115389877A (zh) | 电缆绝缘故障的定位方法、装置、终端及存储介质 | |
CN104748704A (zh) | 薄壁结构超声共振测厚频谱分析内插校正方法 | |
Sternharz et al. | Comparative performance assessment of methods for operational modal analysis during transient order excitation | |
CN103201639B (zh) | 使用连续扫描频率的系统频率响应测试 | |
GB2594760A (en) | Acoustic resonance fluid flow measurement device and method | |
CN115494303A (zh) | 一种emi接收机信号转化方法、装置及存储介质 | |
CN114422318B (zh) | 一种信号预失真的系统响应估计方法、装置及系统 | |
CN105372493B (zh) | 基于三条dft复数谱线的信号幅值和相位测量方法 | |
JP2006523848A (ja) | 複数の光学特性の単掃引測定 | |
Giaquinto et al. | Accuracy analysis in the estimation of ToF of TDR signals | |
Sottek et al. | High-resolution spectral analysis (HSA) vs. discrete fourier transform (DFT) | |
CN109490853A (zh) | 一种线性调频脉冲信号中心频率处谱线值确定方法 | |
Wolf et al. | Amplitude and frequency estimator for aperiodic multi-frequency noisy vibration signals of a tram gearbox |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |