CN115587288A - 一种非2次幂特殊点数频谱计算方法及系统 - Google Patents

一种非2次幂特殊点数频谱计算方法及系统 Download PDF

Info

Publication number
CN115587288A
CN115587288A CN202211562557.5A CN202211562557A CN115587288A CN 115587288 A CN115587288 A CN 115587288A CN 202211562557 A CN202211562557 A CN 202211562557A CN 115587288 A CN115587288 A CN 115587288A
Authority
CN
China
Prior art keywords
spectrum
calculation
channel
calculating
unit
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.)
Granted
Application number
CN202211562557.5A
Other languages
English (en)
Other versions
CN115587288B (zh
Inventor
刘飞
段然
马晓耘
樊蕊蕊
张来宇
廖国福
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
National Astronomical Observatories of CAS
Original Assignee
National Astronomical Observatories of CAS
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by National Astronomical Observatories of CAS filed Critical National Astronomical Observatories of CAS
Priority to CN202211562557.5A priority Critical patent/CN115587288B/zh
Publication of CN115587288A publication Critical patent/CN115587288A/zh
Application granted granted Critical
Publication of CN115587288B publication Critical patent/CN115587288B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • 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/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • 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/15—Correlation function computation including computation of convolution operations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F7/00—Methods or arrangements for processing data by operating upon the order or content of the data handled
    • G06F7/38—Methods or arrangements for performing computations using exclusively denominational number representation, e.g. using binary, ternary, decimal representation
    • G06F7/48—Methods or arrangements for performing computations using exclusively denominational number representation, e.g. using binary, ternary, decimal representation using non-contact-making devices, e.g. tube, solid state device; using unspecified devices
    • G06F7/544—Methods or arrangements for performing computations using exclusively denominational number representation, e.g. using binary, ternary, decimal representation using non-contact-making devices, e.g. tube, solid state device; using unspecified devices for evaluating functions by calculation
    • G06F7/552—Powers or roots, e.g. Pythagorean sums
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02D—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN INFORMATION AND COMMUNICATION TECHNOLOGIES [ICT], I.E. INFORMATION AND COMMUNICATION TECHNOLOGIES AIMING AT THE REDUCTION OF THEIR OWN ENERGY USE
    • Y02D30/00—Reducing energy consumption in communication networks
    • Y02D30/70—Reducing energy consumption in communication networks in wireless communication networks

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • General Engineering & Computer Science (AREA)
  • Software Systems (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Computing Systems (AREA)
  • Probability & Statistics with Applications (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Operations Research (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Complex Calculations (AREA)

Abstract

本发明公开了一种非2次幂特殊点数频谱计算方法及系统，所述方法包括以下步骤：S1.系统参数计算单元接收观测参数，生成计算框架参数值及相位权重因子；S2.多相寄存器组单元对采样数据进行多相分解和载入缓存；S3.频谱计算单元根据计算框架参数值计算通道频谱和精细频谱；S4.通道选择单元根据计算框架频段通道参数值选择相应通道；S5.功率谱计算单元根据频谱计算单元的频谱输出计算相应的功率谱。本技术方案能够达到最优精细频谱分辨率；可对频谱计算单元进行灵活伸缩修改，在满足最优通道分辨率和精细频谱分辨率的前提下，达到“最小计算延迟”或“最少计算资源”约束目标。

Description

一种非2次幂特殊点数频谱计算方法及系统
技术领域
本发明涉及天文信号处理领域，尤其是涉及一种非2次幂特殊点数频谱计算方法及系统。
背景技术
探测到原初引力波，将会是首次探测到引力波后，“又一个首次”的重大发现。我国原初引力波探测计划将是北半球第一个地面宇宙微波背景辐射（Cosmic MicrowaveBackground，CMB）偏振观测实验，是我国“引力波探测计划”的地面实验部分，目标实现对北半球天区的 CMB 偏振信号的首次观测，给出北天区对原初引力波最为精确的测量，同时结合南半球观测数据实现 CMB 极化地面观测的全天覆盖。
原初引力波望远镜室温电子学系统是原初引力波望远镜所必须的关键环节和我国自主研制的核心技术。原初引力波探测所需的读出电子学系统需要处理极端微弱的 CMB信号（对噪声有极致的要求）、处理高速的数据流（最大程度的保存科学数据）、实时计算并应用复杂算法（最大限度优化每一个探测器性能并满足时间、频率精度要求）、保持高稳定性（长时间积分），并在此基础上读出数万个超导转变边沿探测器（Transition EdgeSensor, TES）来实现高灵敏度要求。
大规模超导相变边缘探测器（TES）和动态电感探测器（KID）室温读出，往往采用现场可编程逻辑阵列（Field-Programmable Gate Array，FPGA）与数模和模数转换器件（Analog-to-Digital Converter，ADC 和 Digital to Analog Converter，DAC）相结合的方式。TES 和动态电感探测器（Kinetic Inductive Detector，KID）具有很宽的频率适应范围和多领域应用场景，全世界各类项目已经累计应用超过10万像素。项目拟研发的室温电子学读出系统，是此类探测器系统的关键环节。项目规划的技术指标与世界水平相当，项目的技术积累不仅可以应用于我国原初引力波望远镜和相关探测器读出的应用场景，更具有很好的推广价值和前景。
引力波探测的电子学读出系统对天文超导探测器阵列输出的宽带观测信号进行实时高速采集和频谱通道化，选择探测器谐振频率所在频率通道，获取更高分辨率的精细频谱。与此对应的问题是：1) 尽管FPGA适于信号实时处理，但其内部IP核只能计算2的幂次方信号点数的频谱，不能达到观测所需的最优频率通道分辨率；2) 另一方面，直接使用内部IP核想要获得精细频谱分辨率，需要超大点数2的幂次方信号频谱计算，这就意味着将消耗大量计算资源。
发明内容
针对现有技术存在的问题，本发明的目的在于提出一种非2次幂特殊点数频谱计算方法，用于解决上述技术问题。此种特殊点数频谱计算方法的核心在于根据观测目标要求的采样率、通道频谱分辨率、选择频段、精细频谱分辨率、时间分辨率做总体目标规划，利用“多相分解-抽取频谱”计算方法达到“最优”通道频谱分辨率和精细频谱分辨率；在此前提下，根据观测时间分辨率和计算平台资源的约束，考虑对计算单元进行“最小计算延迟”或是“最少计算资源”的具体实现。
为实现上述目的，本发明提供一种非2次幂特殊点数频谱计算方法，所述方法包括以下步骤：
S1.系统参数计算单元接收观测参数，生成计算框架参数值及相位权重因子；
S2.多相寄存器组单元对采样数据进行多相分解和载入缓存；
S3. 频谱计算单元根据计算框架参数值计算通道频谱和精细频谱；
S4. 通道选择单元根据计算框架频段通道参数值选择相应通道；
S5. 功率谱计算单元根据频谱计算单元的频谱输出计算相应的功率谱。
进一步，步骤S1中，观测参数包括采样率、选择频段、通道分辨率、精细频谱分辨率、和时间分辨率，系统参数计算单元通过非线性整数规划约束算法计算得到计算框架所需的参数值；并计算多相寄存器组单元所需的相位权重因子。
进一步，步骤S1中包括以下子步骤：
步骤S11.非线性整数规划约束算法计算参数
Figure 673513DEST_PATH_IMAGE001
具体如下：
计算参数
Figure 308894DEST_PATH_IMAGE001
和采样率
Figure 137173DEST_PATH_IMAGE002
决定通道分辨率和精细频谱分辨率分别为
Figure 465386DEST_PATH_IMAGE003
和
Figure 480746DEST_PATH_IMAGE004
；观测目标给定通道分辨率
Figure 123080DEST_PATH_IMAGE005
，以及精细频谱分辨率
Figure 664920DEST_PATH_IMAGE006
，“最优”通道分辨率和精细频谱分辨率可按最小均方误差定义，其中
Figure 803515DEST_PATH_IMAGE007
代表模块化频谱计算的点数；
Figure 430806DEST_PATH_IMAGE008
为可实现的频谱计算模块数量，另外在均方误差最小时还需满足频谱分辨率不差于
Figure 876830DEST_PATH_IMAGE005
和
Figure 679701DEST_PATH_IMAGE006
；
因此可以定义非线性整数规划问题如下：
I)
Figure 349717DEST_PATH_IMAGE009
Figure 74091DEST_PATH_IMAGE010
II)
Figure 323806DEST_PATH_IMAGE011
Figure 840238DEST_PATH_IMAGE012
对于I)求解出最优值
Figure 822101DEST_PATH_IMAGE013
，
Figure 361667DEST_PATH_IMAGE014
；
注意到
Figure 211811DEST_PATH_IMAGE015
和
Figure 962510DEST_PATH_IMAGE016
对应的解空间数量最大为
Figure 177591DEST_PATH_IMAGE017
，因此可用穷举法搜索求最优解，以问题I)为例，具体过程包括如下步骤：
S111.构造解空间集合
Figure 266769DEST_PATH_IMAGE018
（
Figure 530392DEST_PATH_IMAGE019
为自然数集；
Figure 21416DEST_PATH_IMAGE020
为整数集）
S112.搜索解空间，找到满足约束的解空间子集
Figure 407398DEST_PATH_IMAGE021
，
Figure 593660DEST_PATH_IMAGE022
为满足I)中
Figure 785607DEST_PATH_IMAGE023
和
Figure 6504DEST_PATH_IMAGE024
的特定解
S113.对于子集
Figure 360124DEST_PATH_IMAGE025
，计算目标函数集合
Figure 361579DEST_PATH_IMAGE026
，
Figure 731118DEST_PATH_IMAGE027
为取特定值
Figure 931155DEST_PATH_IMAGE028
时,和I)相关的目标函数
S114.按照升序对
Figure 658939DEST_PATH_IMAGE029
进行排序，得到索引向量
Figure 85373DEST_PATH_IMAGE030
Figure 619122DEST_PATH_IMAGE031
即为所求最优解，
Figure 876928DEST_PATH_IMAGE032
表示
Figure 713297DEST_PATH_IMAGE033
在索引向量
Figure 486081DEST_PATH_IMAGE034
中取值时的
Figure 433309DEST_PATH_IMAGE035
和
Figure 545621DEST_PATH_IMAGE036
的特定解
S115.构造解空间集合
Figure 677525DEST_PATH_IMAGE037
S116.搜索解空间，找到满足约束的解空间子集
Figure 311507DEST_PATH_IMAGE038
，
Figure 187059DEST_PATH_IMAGE039
为满足II)中
Figure 153878DEST_PATH_IMAGE040
和
Figure 332049DEST_PATH_IMAGE041
的特定解
S117.对于子集
Figure 79425DEST_PATH_IMAGE038
，计算目标函数集合
Figure 368455DEST_PATH_IMAGE042
，
Figure 455360DEST_PATH_IMAGE043
为取特定值
Figure 194646DEST_PATH_IMAGE044
时的和II)相关的目标函数
S118.按照升序对目标函数集内的所有元素
Figure 39105DEST_PATH_IMAGE045
进行排序，得到索引向量
Figure 256460DEST_PATH_IMAGE046
Figure 197871DEST_PATH_IMAGE047
即为
Figure 216380DEST_PATH_IMAGE048
所求的最优解，
Figure 938349DEST_PATH_IMAGE049
表示
Figure 569181DEST_PATH_IMAGE050
在索引向量
Figure 427416DEST_PATH_IMAGE051
中取值时的
Figure 852712DEST_PATH_IMAGE052
和
Figure 999660DEST_PATH_IMAGE041
的特定解
S12.计算多相寄存器组单元所需的相位权重因子如下：
多相寄存器组单元通道相位权重因子计算
第一级多相分解频谱计算的通道相位因子定义为：
Figure 558817DEST_PATH_IMAGE053
，
类似可定义第一级多相分解频谱计算的通道相位因子：
Figure 146924DEST_PATH_IMAGE054
。
进一步，步骤S2中，多相寄存器组单元接收系统参数计算单元计算得到的P、Q参数值，形成相应寄存器组，缓存多相分解后的采样数据，作为频谱计算单元的输入；多相寄存器组单元包括第一级多相寄存器组和第二级多相寄存器组，其长度分别为
Figure 398914DEST_PATH_IMAGE055
和
Figure 932358DEST_PATH_IMAGE056
，
Figure 29627DEST_PATH_IMAGE057
满足“最优频谱分辨率”，输入信号采样率为
Figure 737820DEST_PATH_IMAGE058
；将输入信号进行多相分解后载入到第一级多相寄存器组，第一级多相寄存器组后的第一级频谱计算单元将实现频谱通道化计算和通道选择，通道分辨率为
Figure 895132DEST_PATH_IMAGE059
；第二级多相寄存器组及第二级频谱计算单元在通道选择后实现精细频谱计算，精细频谱分辨率为
Figure 751093DEST_PATH_IMAGE060
。
进一步，步骤S2中，多相寄存器组单元对采样数据进行多相分解计算频谱的方式包括以下步骤：
S21. 多相分解计算频谱。对
Figure 527419DEST_PATH_IMAGE061
点窗函数截断后的采样序列作频谱计算，将
Figure 480331DEST_PATH_IMAGE061
分解为
Figure 683911DEST_PATH_IMAGE062
，
Figure 27167DEST_PATH_IMAGE063
为每一多相分支上的快拍数(tap)，
Figure 466239DEST_PATH_IMAGE064
为多相分解后需计算的频谱通道数。对于N点多相分解后频谱计算来讲，输出的N个频率格(bin)和直接作M点频谱计算后输出的M个bin之间存在以下抽取关系：
对于M个bin，从首个bin开始，每
Figure 647559DEST_PATH_IMAGE065
个bin抽取1个，得到
Figure 146674DEST_PATH_IMAGE064
个抽取后的bin；从bin[1]开始抽取，然后每
Figure 649330DEST_PATH_IMAGE065
个bin抽取1个，得到另一组
Figure 892093DEST_PATH_IMAGE064
个抽取后的bin；依次类推，如果从
Figure 163805DEST_PATH_IMAGE066
分别开始抽取，最终可得到
Figure 37084DEST_PATH_IMAGE067
个抽取后的bin；
在用多相分解频谱计算时，需要引入额外的相位因子实现抽取初始位置的改变，从而使
Figure 151670DEST_PATH_IMAGE067
个bin和
Figure 339069DEST_PATH_IMAGE061
个bin完全对应；
Figure 589922DEST_PATH_IMAGE068
为采样序列，
Figure 804740DEST_PATH_IMAGE069
为
Figure 406623DEST_PATH_IMAGE061
点截断窗函数，计算
Figure 866554DEST_PATH_IMAGE070
的
Figure 440755DEST_PATH_IMAGE061
点频谱可作多相分解如下：
Figure 452573DEST_PATH_IMAGE071
Figure 682697DEST_PATH_IMAGE072
Figure 805374DEST_PATH_IMAGE073
Figure 906185DEST_PATH_IMAGE074
多相分解计算后的
Figure 354484DEST_PATH_IMAGE075
点频谱
Figure 868642DEST_PATH_IMAGE076
是
Figure 168911DEST_PATH_IMAGE061
点原始频谱
Figure 780021DEST_PATH_IMAGE077
的
Figure 9008DEST_PATH_IMAGE065
倍抽取，抽取初始位置为
Figure 807200DEST_PATH_IMAGE078
；
Figure 412625DEST_PATH_IMAGE079
和
Figure 815924DEST_PATH_IMAGE080
对应着
Figure 340447DEST_PATH_IMAGE081
和
Figure 501301DEST_PATH_IMAGE082
。
定义抽取索引数
Figure 35050DEST_PATH_IMAGE083
，
Figure 469355DEST_PATH_IMAGE084
对应着抽取初始位置为
Figure 633620DEST_PATH_IMAGE085
的多相分解计算频谱bin；
Figure 406404DEST_PATH_IMAGE086
相对于
Figure 822473DEST_PATH_IMAGE087
的计算，
Figure 872468DEST_PATH_IMAGE088
需要在第r个快拍(tap)上施加
Figure 942055DEST_PATH_IMAGE089
和第n个通道上施加
Figure 467715DEST_PATH_IMAGE090
；
Figure 280950DEST_PATH_IMAGE091
可以被吸收到窗函数的多相分支
Figure 949566DEST_PATH_IMAGE092
中，令
Figure 455634DEST_PATH_IMAGE093
第n个通道上的因子定义为：
Figure 203010DEST_PATH_IMAGE094
。
进一步，步骤S2中，还包括以下步骤：
S22. 修正窗函数值的计算；
多相分解后，窗函数被修正为：
Figure 819936DEST_PATH_IMAGE095
对于第一级多相分解频谱计算，其窗函数的多相分解形式为
Figure 844524DEST_PATH_IMAGE096
，修正后的窗函数多相分解形式为
Figure 318231DEST_PATH_IMAGE097
，
Figure 490586DEST_PATH_IMAGE098
为频率通道序号。
对于第二级多相分解频谱，其窗函数的多相分解形式为
Figure 583307DEST_PATH_IMAGE099
，修正后的窗函数多相分解式为：
Figure 790297DEST_PATH_IMAGE100
，
Figure 700485DEST_PATH_IMAGE101
为精细频谱通道序号。
进一步，步骤S3中，频谱计算单元接收来自系统参数计算单元的
Figure 360136DEST_PATH_IMAGE102
参数值和相位权重因子、来自多相寄存器组单元的多相分解采样数据，频谱计算单元中的FFT模块对相位加权后的多相信号进行频谱计算，频谱计算的结果被输入到功率谱计算单元中，最终得到通道功率和精细频谱功率。
进一步，根据约束目标的不同，如“最小计算延迟”和“最少计算资源”，频谱计算单元可做灵活伸缩修改，采用全并行频谱计算实现，或者时分复用计算实现。
另一方面，本发明提供一种非2次幂特殊点数频谱计算系统，所述系统用于实现根据本发明所述的特殊点数频谱计算方法，系统包括系统参数计算单元、多相寄存器组A、频谱计算单元A、多相寄存器组B、频谱计算单元B、通道选择单元和功率谱计算单元。
进一步，采样信号首先被输入到多相寄存器组A中，在信号通路上依次经过频谱计算单元A、通道选择单元、多相寄存器组B、频谱计算单元B，最终由功率谱计算单元输出通道功率和精细频谱功率。
本发明的有益效果在于：能够实现最优通道频率分辨率，其与特殊点数信号的直接DFT频谱分辨率的误差最小；能够灵活选择所需观测频段，对观测频段覆盖通道进行超精细频谱分辨，并达到最优精细频谱分辨率；可对频谱计算单元进行灵活伸缩修改，在满足最优通道分辨率和精细频谱分辨率的前提下，达到“最小计算延迟”或“最少计算资源”约束目标。
附图说明
图1示出了根据本发明实施例一种非2次幂特殊点数频谱计算方法及系统的架构示意图；
图2示出了根据本发明实施例一种非2次幂特殊点数频谱计算方法及系统的技术框图；
图3示出了根据本发明实施例实现“最小计算延迟”的技术框图；
图4示出了根据本发明实施例实现“最少计算资源”的技术框图。
具体实施方式
下面将结合附图，对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
在本发明的描述中，需要说明的是，术语“中心”、“上”、“下”、“左”、“右”、“竖直”、“水平”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限制。此外，术语“第一”、“第二”、“第三”仅用于描述目的，而不能理解为指示或暗示相对重要性。
在本发明的描述中，需要说明的是，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以根据具体情况理解上述术语在本发明中的具体含义。
以下结合图1-图4对本发明的具体实施方式进行详细说明。应当理解的是，此处所描述的具体实施方式仅用于说明和解释本发明，并不用于限制本发明。
为实现原初引力波探测计划，实现对北半球天区的 CMB 偏振信号的观测，给出北天区对原初引力波最为精确的测量，本发明提出一种特殊点数频谱计算方法，所述特殊点数为非2次幂点数，其核心在于根据观测目标要求的采样率、通道频谱分辨率、选择频段、精细频谱分辨率和时间分辨率做总体目标规划，利用“多相分解-抽取频谱”计算方法达到“最优”通道频谱分辨率和精细频谱分辨率；在此前提下，根据观测时间分辨率和计算平台资源的约束，考虑对频谱计算单元进行“最小计算延迟”或是“最少计算资源”的具体实现。
图1示出根据本发明优选实施例的非2次幂特殊点数频谱计算系统，整个计算系统构架包括的基本功能单元有：系统参数计算单元，多相寄存器组A，频谱计算单元A，多相寄存器组B，频谱计算单元B，通道选择单元，功率谱计算单元。
系统参数计算单元：接收观测参数，生成计算框架参数值。
多相寄存器组单元：对采样数据进行多相分解和载入缓存。
频谱计算单元：根据计算框架参数值计算通道频谱和精细频谱。
通道选择单元：根据计算框架频段通道参数值选择相应通道。
功率谱计算单元：根据频谱计算单元的频谱输出计算相应的功率谱。
系统参数计算单元其输入为观测参数，观测参数主要包括通道分辨率，选择频段，精细频谱分辨率，信号采样率；多相寄存器组单元包括多相寄存器A和多相寄存器B；频谱计算单元包括频谱计算单元A和频谱计算单元B。采样信号首先被输入到多相寄存器组A中，在信号通路上依次经过频谱计算单元A、通道选择单元、多相寄存器组B、频谱计算单元B，最终由功率谱计算单元输出通道功率和精细频谱功率。
如图2所示“多相分解-抽取频谱”的关键是利用参数
Figure 990969DEST_PATH_IMAGE102
形成两级多相寄存器组单元和频谱计算单元。其中多项寄存器组单元A对采样信号进行长度为
Figure 849203DEST_PATH_IMAGE103
的多相分解和加载，频谱计算单元A将多相分解信号进行相位加权后经
Figure 38614DEST_PATH_IMAGE104
点FFT模块计算抽取频谱，由于每个频谱计算单元A得到的是
Figure 185561DEST_PATH_IMAGE105
个频谱通道经过不同相位偏移的
Figure 744719DEST_PATH_IMAGE035
倍抽取，所有频谱计算单元A的频谱计算结果综合起来可得到完整的
Figure 598405DEST_PATH_IMAGE105
个频谱通道。对于多相寄存器组单元B和频谱计算单元B，依据同样的方法可得到
Figure 319237DEST_PATH_IMAGE106
个精细频谱通道。
“多相分解-抽取频谱”的频谱计算部分可利用2次幂小点数FFT IP核，从而兼顾了频谱分辨率的观测目标和频谱计算实现效率。
对射电天文宽带观测信号进行4G/s的高速采样，采样带宽为4~6 GHz。观测要求在采样带宽进行“粗”频谱通道化，通道分辨率
Figure 953480DEST_PATH_IMAGE107
，并选出总带宽的1/4实现精细频谱分辨，精细频谱分辨率达到通道带宽的1/40，观测模式为“快速捕捉模式”时，时间分辨率<50us。
依据不同计算实施平台的特点和实际观测需求，可将约束目标设定为“最小计算延迟”和“最少计算资源”。在实施例中，观测需求对时间分辨率有一定要求，经过实测，在FPGA中，FFT内核在1024点和8点计算时的延迟上限为12 us和1us，频谱计算单元和数据流水寄存器导致的计算延时在全并行情况下可控制在20us以内，从而满足观测时间分辨率要求。
此方法为天文观测场景中由“最优”观测分辨率所引入的特殊点数频谱计算问题提供了灵活高效的实现方案，适于在各类计算平台架构上(如CPU、FPGA、GPU等)有效实现。
本发明优选实施例的一种非2次幂特殊点数频谱计算方法，包括以下步骤：
S1.系统参数计算单元接收观测参数，生成计算框架参数值及相位权重因子；
S2.多相寄存器组单元对采样数据进行多相分解和载入缓存；
S3. 频谱计算单元：根据计算框架参数值计算通道频谱和精细频谱；
S4. 通道选择单元根据计算框架频段通道参数值选择相应通道；
S5. 功率谱计算单元根据频谱计算单元的频谱输出计算相应的功率谱。
步骤S1中，系统参数计算单元根据采样率、选择频段、通道频谱分辨率、精细频谱分辨率、时间分辨率等输入观测参数，通过非线性整数规划约束算法(最小均方误差准则)计算得到计算框架所需的
Figure 926116DEST_PATH_IMAGE108
参数值；并计算多相寄存器组单元所需的相位权重因子。其中，包括以下子步骤：
步骤S11.非线性整数规划约束算法计算参数
Figure 758942DEST_PATH_IMAGE057
，其中
Figure 791620DEST_PATH_IMAGE035
为第一级滤波器组多相寄存器组的数量，
Figure 709898DEST_PATH_IMAGE104
为第一级滤波器组频谱计算的点数；类似
Figure 814120DEST_PATH_IMAGE052
为第二级滤波器组多相寄存器组的数量，
Figure 875355DEST_PATH_IMAGE109
为第二级滤波器组频谱计算的点数。具体计算过程具体如下：
计算参数
Figure 203568DEST_PATH_IMAGE102
和采样率
Figure 546825DEST_PATH_IMAGE058
决定通道分辨率和精细频谱分辨率分别为
Figure 126842DEST_PATH_IMAGE059
和
Figure 934261DEST_PATH_IMAGE060
；观测目标给定通道分辨率
Figure 308741DEST_PATH_IMAGE110
，以及精细频谱分辨率
Figure 139294DEST_PATH_IMAGE111
，“最优”通道分辨率和精细频谱分辨率可按最小均方误差定义，其中
Figure 382056DEST_PATH_IMAGE112
代表模块化频谱计算的点数，在实际工程可取8~1024，对应
Figure 919348DEST_PATH_IMAGE113
；
Figure 589364DEST_PATH_IMAGE114
为可实现的频谱计算模块数量，考虑在FPGA中实现的资源约束，取上限为50(也可以为更大的值)。另外在均方误差最小时还需满足频谱分辨率不差于
Figure 907213DEST_PATH_IMAGE110
和
Figure 79566DEST_PATH_IMAGE111
。
因此可以定义非线性整数规划问题如下：
I)
Figure 595998DEST_PATH_IMAGE115
Figure 577861DEST_PATH_IMAGE116
II)
Figure 383006DEST_PATH_IMAGE117
Figure 967571DEST_PATH_IMAGE118
对于I)求解出最优值
Figure 479455DEST_PATH_IMAGE013
，
Figure 491273DEST_PATH_IMAGE014
。
注意到
Figure 518135DEST_PATH_IMAGE119
和
Figure 781757DEST_PATH_IMAGE120
对应的解空间数量最大为
Figure 272782DEST_PATH_IMAGE121
，因此可用穷举法搜索求最优解，以问题I)为例，具体过程包括如下步骤：
S111.构造解空间集合
Figure 829403DEST_PATH_IMAGE018
（
Figure 609140DEST_PATH_IMAGE019
为自然数集，
Figure 535507DEST_PATH_IMAGE020
为整数集)
S112.搜索解空间，找到满足约束的解空间子集
Figure 756404DEST_PATH_IMAGE021
，
Figure 375605DEST_PATH_IMAGE122
为满足I)中
Figure 377059DEST_PATH_IMAGE123
和
Figure 248063DEST_PATH_IMAGE024
的特定解
S113.对于子集
Figure 182521DEST_PATH_IMAGE025
，计算目标函数集合
Figure 175884DEST_PATH_IMAGE026
，
Figure 336738DEST_PATH_IMAGE124
为取特定值
Figure 604909DEST_PATH_IMAGE028
时的和I)相关的目标函数
S114.按照升序对
Figure 298933DEST_PATH_IMAGE029
进行排序，得到索引向量
Figure 259936DEST_PATH_IMAGE125
Figure 970403DEST_PATH_IMAGE031
即为所求
Figure 917630DEST_PATH_IMAGE126
的最优解，
Figure 92260DEST_PATH_IMAGE032
表示
Figure 99530DEST_PATH_IMAGE033
在索引向量
Figure 359610DEST_PATH_IMAGE034
中取值时的
Figure 844949DEST_PATH_IMAGE035
和
Figure 874085DEST_PATH_IMAGE036
的特定解。
S115.构造解空间集合
Figure 114573DEST_PATH_IMAGE127
S116.搜索解空间，找到满足约束的解空间子集
Figure 235851DEST_PATH_IMAGE038
，
Figure 649515DEST_PATH_IMAGE128
为满足II)中
Figure 674102DEST_PATH_IMAGE052
和
Figure 147809DEST_PATH_IMAGE129
的特定解
S117.对于子集
Figure 320165DEST_PATH_IMAGE130
，计算目标函数集合
Figure 147306DEST_PATH_IMAGE131
，
Figure 151034DEST_PATH_IMAGE132
为取特定值
Figure 405429DEST_PATH_IMAGE133
时的和II)相关的目标函数
S118.按照升序对
Figure 861818DEST_PATH_IMAGE045
进行排序，得到索引向量
Figure 262625DEST_PATH_IMAGE134
Figure 120859DEST_PATH_IMAGE047
即为所求
Figure 811735DEST_PATH_IMAGE048
的最优解
Figure 958682DEST_PATH_IMAGE135
表示
Figure 252260DEST_PATH_IMAGE050
在索引向量
Figure 105947DEST_PATH_IMAGE051
中取值时的
Figure 92358DEST_PATH_IMAGE052
和
Figure 398705DEST_PATH_IMAGE136
的特定解。
S12.计算多相寄存器组单元所需的相位权重因子如下：
多相寄存器组单元通道相位权重因子计算
第一级多相分解频谱计算的通道相位因子定义为：
Figure 495974DEST_PATH_IMAGE053
，
类似可定义第一级多相分解频谱计算的通道相位因子：
Figure 266484DEST_PATH_IMAGE137
步骤S2中，多相寄存器组单元接收系统参数计算单元计算得到的
Figure 797697DEST_PATH_IMAGE035
、
Figure 981554DEST_PATH_IMAGE052
参数值，形成相应寄存器组，缓存多相分解后的采样数据，作为频谱计算单元的输入。本发明设计实现两级多相寄存器组，分别为第一级多相寄存器组和第二级多相寄存器组（第一级多相寄存器组包含若干多相寄存器组A, 第二级多相寄存器组包含若干多相寄存器组B），其长度分别为
Figure 492301DEST_PATH_IMAGE055
和
Figure 382896DEST_PATH_IMAGE056
，
Figure 711110DEST_PATH_IMAGE138
满足“最优频谱分辨率”(最小均方误差准则)，输入信号采样率为
Figure 992049DEST_PATH_IMAGE139
；将输入信号进行多相分解后载入到第一级多相寄存器组，第一级多相寄存器组后的第一级频谱计算单元将实现频谱通道化计算和通道选择，通道分辨率为
Figure 696700DEST_PATH_IMAGE059
；第二级多相寄存器组及第二级频谱计算单元在通道选择后实现精细频谱计算，精细频谱分辨率为
Figure 441802DEST_PATH_IMAGE140
；此框架不仅能够实现非2次幂特殊点数频谱有效计算，同时还具有灵活挑选观测频段、并进行精细频谱分辨识别的能力。
具体地，步骤S2中多相寄存器组单元对采样数据进行多相分解计算频谱的方式如下：
S21. 多相分解计算频谱。对
Figure 816283DEST_PATH_IMAGE061
点窗函数截断后的采样序列作频谱计算，将
Figure 709153DEST_PATH_IMAGE061
分解为
Figure 889598DEST_PATH_IMAGE067
，其中
Figure 925425DEST_PATH_IMAGE065
为多相分解后每一分支上的快拍(tap)数，
Figure 595441DEST_PATH_IMAGE064
为多相分解后需要计算的频谱格数(bin)，对于
Figure 850973DEST_PATH_IMAGE064
点多相分解后频谱计算来讲，输出的
Figure 835109DEST_PATH_IMAGE064
个频率格(bin)和直接作
Figure 351541DEST_PATH_IMAGE061
点频谱计算后输出的
Figure 333404DEST_PATH_IMAGE061
个bin之间存在以下抽取关系：
对于
Figure 935286DEST_PATH_IMAGE061
个bin，从从首个bin(bin[0])开始，每
Figure 395218DEST_PATH_IMAGE065
个bin抽取1个，得到
Figure 969418DEST_PATH_IMAGE064
个抽取后的bin；
如果改变抽取的初始位置，从bin[1]开始抽取，然后每
Figure 981237DEST_PATH_IMAGE065
个bin抽取1个，得到另一组N个抽取后的bin。依次类推，如果从bin[2],…bin[R-1]分别开始抽取，最终可得到
Figure 444317DEST_PATH_IMAGE067
个抽取后的bin；
在用多相分解频谱计算时，需要引入额外的相位因子实现抽取初始位置的改变，从而使
Figure 35835DEST_PATH_IMAGE067
个bin和
Figure 261280DEST_PATH_IMAGE061
个bin完全对应；
假设
Figure 584945DEST_PATH_IMAGE141
为采样序列，
Figure 364682DEST_PATH_IMAGE142
为
Figure 291050DEST_PATH_IMAGE061
点截断窗函数，计算
Figure 511947DEST_PATH_IMAGE070
的
Figure 131147DEST_PATH_IMAGE061
点频谱
Figure 132601DEST_PATH_IMAGE143
可作多相分解如下：
Figure 472447DEST_PATH_IMAGE071
Figure 938063DEST_PATH_IMAGE144
Figure 931427DEST_PATH_IMAGE145
Figure 319378DEST_PATH_IMAGE074
这意味着多相分解后的N点计算频谱
Figure 853127DEST_PATH_IMAGE146
是
Figure 48616DEST_PATH_IMAGE061
点原始频谱
Figure 947302DEST_PATH_IMAGE077
的
Figure 720086DEST_PATH_IMAGE147
倍抽取，抽取初始位置为
Figure 932893DEST_PATH_IMAGE148
；
Figure 841943DEST_PATH_IMAGE149
和
Figure 177109DEST_PATH_IMAGE080
对应着
Figure 312556DEST_PATH_IMAGE087
和
Figure 188108DEST_PATH_IMAGE150
，因抽取跳过了
Figure 420506DEST_PATH_IMAGE151
这些点；
定义抽取索引数(decimation index)
Figure 97213DEST_PATH_IMAGE083
，
Figure 844589DEST_PATH_IMAGE152
对应着抽取初始位置为
Figure 461515DEST_PATH_IMAGE153
的多相分解计算频谱bin；
Figure 220524DEST_PATH_IMAGE086
相对于
Figure 428651DEST_PATH_IMAGE154
的计算，
Figure 538689DEST_PATH_IMAGE084
需要在第r个快拍(tap)上施加相位因子
Figure 756044DEST_PATH_IMAGE089
和第n个通道上施加相位因子
Figure 697455DEST_PATH_IMAGE155
；
Figure 483009DEST_PATH_IMAGE156
可以被吸收到窗函数的多相分支
Figure 939398DEST_PATH_IMAGE092
中，令
Figure 898127DEST_PATH_IMAGE093
第n个通道上的附加相位因子定义为：
Figure 130263DEST_PATH_IMAGE157
。
S22. 修正窗函数值的计算；
多相分解后，窗函数被修正为：
Figure 211351DEST_PATH_IMAGE095
对于第一级多相分解频谱计算，其窗函数的多相分解形式为
Figure 92719DEST_PATH_IMAGE158
，修正后的窗函数多相分解形式为
Figure 527243DEST_PATH_IMAGE097
，
Figure 505563DEST_PATH_IMAGE159
为频率通道序号。
对于第二级多相分解频谱，其窗函数的多相分解形式为
Figure 695236DEST_PATH_IMAGE160
，修正后的窗函数多相分解式为：
Figure 1584DEST_PATH_IMAGE100
，
Figure 98853DEST_PATH_IMAGE101
为精细频谱通道序号。
采用多相分解计算频谱的方式，能够将实际频谱分辨率尽可能贴近目标频谱分辨率，并且通过并行频谱计算单元降低系统延时。
步骤S3中，频谱计算单元接收来自系统参数计算单元的
Figure 807046DEST_PATH_IMAGE102
参数值和相位权重因子、来自多相寄存器组单元的多相分解采样数据，结合FFT模块共同计算通道频谱和精细频谱，作为功率谱计算单元的输入。频谱计算单元中的FFT模块对相位加权后的多相信号进行频谱计算，频谱计算的结果被输入到功率谱计算单元中，最终得到通道功率和精细频谱功率。
另外根据约束目标的不同，如“最小计算延迟”和“最少计算资源”，频谱计算单元可做灵活伸缩修改，采用全并行频谱计算实现，或者时分复用(TDM: time-domain-multiplexer)计算实现。
系统参数计算单元根据输入观测参数，通过约束最优化算法(最小均方误差)计算得到参数值；根据多相寄存器组中不同抽取索引计算相应的相位因子，然后将相位因子和窗函数对应值结合得到频谱计算单元中的最终相位权重。
本实施例中的观测模式为谱线模式，对时间分辨率的要求为~ 1ms，通道分辨率和精细频谱分辨率决定实际计算单元参数
Figure 964357DEST_PATH_IMAGE161
依然不变。此时可将约束目标设定为“最少计算资源”，以降低观测系统计算成本。当以最少计算资源为目标时，实现方式如图4所示：多相寄存器的输出被相位加权后通过时分复用(TDM: time-domain-multiplexer)模块，被同一FFT模块进行计算，这样在保证通道分辨率和精细频谱分辨率的前提下可节省大量DSP资源。
以FPGA计算平台为例，其通用IP核利用直接FFT实现
Figure 820318DEST_PATH_IMAGE162
点数频谱计算实现采样带宽通道化，在通道频谱分辨率往往与观测目标要求的频谱分辨率有一定偏离。实际中往往选择大点数FFT频谱计算，但这样尽管能满足分辨率要求，但往往会导致不必要的过高频谱分辨实现。由于观测的通道频谱分辨率与灵敏度有关，过高频谱分辨会降低实际灵敏度，从而需要在后续处理中增加积分点数，消耗额外的逻辑资源；过高频谱分辨本身也意味着大点数FFT频谱计算会消耗更多的计算资源。因而通过这种方式得到的通道频谱分辨率不是“最优”的。如果采用多相分解计算频谱的方式，能够将实际频谱分辨率尽可能贴近目标频谱分辨率，并且通过并行频谱计算单元降低系统延时。
参数
Figure 360759DEST_PATH_IMAGE102
的另一种计算方法是计算
Figure 48092DEST_PATH_IMAGE163
，此频谱点数为保证
Figure 579567DEST_PATH_IMAGE164
的最小计算点数，在此初值基础上反复用因子2长除，在所有结果中找到满足
Figure 594928DEST_PATH_IMAGE165
和
Figure 299579DEST_PATH_IMAGE166
约束的
Figure 716785DEST_PATH_IMAGE119
，此算法适合迭代计算实现，并且最终结果按照升序排列。
具体算法伪代码如下所示：
Figure 419161DEST_PATH_IMAGE167
问题II)的计算实现和问题I)类似，只是把目标函数从
Figure 312031DEST_PATH_IMAGE168
改为
Figure 430160DEST_PATH_IMAGE169
在一个具体实施例中，以最小计算延迟为目标，则实现方式如图3所示。
1) 参数
Figure 295348DEST_PATH_IMAGE170
的计算
通过先前的算法可以得到：
Figure 965363DEST_PATH_IMAGE171
2) 多相分解频谱通道相位因子计算
2-1) 计算第一级多相分解频谱通道相位因子
第一级多相分支数为
Figure 459711DEST_PATH_IMAGE172
，每个多相分支上进行
Figure 506164DEST_PATH_IMAGE173
点FFT计算，需要计算的相位因子为：
Figure 225858DEST_PATH_IMAGE174
除去
Figure 207721DEST_PATH_IMAGE175
或
Figure 544024DEST_PATH_IMAGE176
(相位因子为1)，共需计算的相位因子数为:
Figure 331852DEST_PATH_IMAGE177
；
抽取索引
Figure 578156DEST_PATH_IMAGE178
的前10个通道相位因子
Figure 855554DEST_PATH_IMAGE179
如下表所示，其他抽取索引的通道相位因子可按公式计算
Figure 147995DEST_PATH_IMAGE180
2-2) 计算第二级多相分解频谱通道相位因子
第二级多相分支数为
Figure 146038DEST_PATH_IMAGE181
，每个多相分支上进行
Figure 902641DEST_PATH_IMAGE182
点FFT计算，需要计算的相位因子为：
Figure 459262DEST_PATH_IMAGE183
除去
Figure 973420DEST_PATH_IMAGE184
或
Figure 165367DEST_PATH_IMAGE185
（相位因子为1），共需计算的相位因子数为:
Figure 386264DEST_PATH_IMAGE186
；
抽取索引
Figure 739885DEST_PATH_IMAGE187
的前10个通道相位因子(
Figure 741339DEST_PATH_IMAGE188
)如下表所示，其他抽取索引的通道相位因子可按公式计算
Figure 346764DEST_PATH_IMAGE189
3）修正窗函数的多相分解式计算
截断窗函数可根据需要选择，在射电天文中一般选用Hanning窗。
Figure 812380DEST_PATH_IMAGE190
下面给出抽取索引
Figure 540165DEST_PATH_IMAGE178
，频率通道索引
Figure 966598DEST_PATH_IMAGE191
第一级修正窗函数的各多相快拍数(tap:
Figure 500348DEST_PATH_IMAGE192
)计算值，其他值可根据公式算出：
Figure 758154DEST_PATH_IMAGE193
下面给出抽取索引
Figure 827479DEST_PATH_IMAGE194
，频率通道索引
Figure 865842DEST_PATH_IMAGE195
时，第二级修正窗函数的各多相分多相快拍数(tap:
Figure 813069DEST_PATH_IMAGE196
)计算值，其他值可根据公式算出：
Figure 925382DEST_PATH_IMAGE197
系统参数计算单元根据输入观测参数，通过约束最优化算法(最小均方误差)计算得到
Figure 322865DEST_PATH_IMAGE198
，进一步根据多相寄存器组中不同抽取索引计算相应的相位因子，然后将相位因子和窗函数对应值结合得到频谱计算单元中的最终相位权重。对于多相寄存器A和频谱计算单元A，其相位权重为
Figure 192732DEST_PATH_IMAGE199
，对于多相寄存器B和频谱计算单元B，其相位权重为
Figure 68284DEST_PATH_IMAGE200
。
Figure 35103DEST_PATH_IMAGE199
的数量为
Figure 213275DEST_PATH_IMAGE201
，
Figure 960651DEST_PATH_IMAGE200
的数量为
Figure 577577DEST_PATH_IMAGE202
。频谱计算单元中的FFT模块对相位加权后的多相信号进行频谱计算，频谱计算的结果被输入到功率谱计算单元中，最终得到通道功率和精细频谱功率。
本发明的技术优势如下：
1、实现最优通道频率分辨率，其与特殊点数信号的直接DFT频谱分辨率的误差最小；
2.灵活选择所需观测频段，对观测频段覆盖通道进行超精细频谱分辨，并达到最优精细频谱分辨率；
3.可对频谱计算单元进行灵活伸缩修改，在满足最优通道分辨率和精细频谱分辨率的前提下，达到“最小计算延迟”或“最少计算资源”约束目标。
在本说明书的描述中，参考术语“实施例”、“示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不必须针对的是相同的实施例或示例。此外，本领域的技术人员可以在不产生矛盾的情况下，将本说明书中描述的不同实施例或示例以及其中的特征进行结合或组合。
上述内容虽然已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在本发明的范围内可以对上述实施例进行变化、修改、替换和变型等更新操作。

Claims (10)

1.一种非2次幂特殊点数频谱计算方法，其特征在于，所述方法包括以下步骤：
S1.系统参数计算单元接收观测参数，生成计算框架参数值及相位权重因子；
S2.多相寄存器组单元对采样数据进行多相分解和载入缓存；
S3. 频谱计算单元根据计算框架参数值计算通道频谱和精细频谱；
S4. 通道选择单元根据计算框架频段通道参数值选择相应通道；
S5. 功率谱计算单元根据频谱计算单元的频谱输出计算相应的功率谱。
2.根据权利要求1所述的非2次幂特殊点数频谱计算方法，其特征在于，步骤S1中，观测参数包括采样率、选择频段、通道分辨率、精细频谱分辨率、和时间分辨率，系统参数计算单元通过非线性整数规划约束算法计算得到计算框架所需的参数值；并计算多相寄存器组单元所需
Figure 166237DEST_PATH_IMAGE001
的相位权重因子。
3.根据权利要求2所述的非2次幂特殊点数频谱计算方法，其特征在于，步骤S1中包括以下子步骤：
步骤S11.非线性整数规划约束算法计算参数
Figure 810582DEST_PATH_IMAGE002
，其中
Figure 779675DEST_PATH_IMAGE003
为第一级滤波器组多相寄存器组的数量，
Figure 200292DEST_PATH_IMAGE004
为第一级滤波器组频谱计算的点数；类似
Figure 813807DEST_PATH_IMAGE005
为第二级滤波器组多相寄存器组的数量，
Figure 723995DEST_PATH_IMAGE006
为第二级滤波器组频谱计算的点数；具体如下：
计算参数
Figure 55750DEST_PATH_IMAGE007
和采样率
Figure 14479DEST_PATH_IMAGE008
决定通道分辨率和精细频谱分辨率分别为
Figure 138292DEST_PATH_IMAGE009
和
Figure 563589DEST_PATH_IMAGE010
；观测目标给定通道分辨率
Figure 507274DEST_PATH_IMAGE011
，以及精细频谱分辨率
Figure 269694DEST_PATH_IMAGE012
，“最优”通道分辨率和精细频谱分辨率可按最小均方误差定义，其中
Figure 621915DEST_PATH_IMAGE013
分别代表第一级和第二级滤波器组频谱计算的点数；
Figure 608326DEST_PATH_IMAGE014
分别为第一级和第二级滤波器组的多相寄存器组数量，另外在均方误差最小时还需满足频谱分辨率不差于
Figure 976990DEST_PATH_IMAGE011
和
Figure 949625DEST_PATH_IMAGE012
；
因此可以定义非线性整数规划问题如下：
I)
Figure 782452DEST_PATH_IMAGE015
Figure 80710DEST_PATH_IMAGE016
II)
Figure 733408DEST_PATH_IMAGE017
Figure 837630DEST_PATH_IMAGE018
对于I)求解出最优值
Figure 665909DEST_PATH_IMAGE019
，
Figure 994122DEST_PATH_IMAGE020
；
注意到
Figure 337379DEST_PATH_IMAGE021
和
Figure 952949DEST_PATH_IMAGE022
对应的解空间数量最大为
Figure 698051DEST_PATH_IMAGE023
，因此可用穷举法搜索求最优解，以问题I)为例，具体过程包括如下步骤：
S111.构造解空间集合
Figure 72531DEST_PATH_IMAGE024
,
Figure 699822DEST_PATH_IMAGE025
为自然数集；
Figure 349109DEST_PATH_IMAGE026
为整数集;
S112.搜索解空间，找到满足约束的解空间子集
Figure 620821DEST_PATH_IMAGE027
，
Figure 556416DEST_PATH_IMAGE028
为满足I)中
Figure 44904DEST_PATH_IMAGE029
和
Figure 294620DEST_PATH_IMAGE030
的特定解;
S113.对于子集
Figure 545473DEST_PATH_IMAGE031
，计算
Figure 527335DEST_PATH_IMAGE032
，
Figure 129218DEST_PATH_IMAGE033
为取特定值
Figure 589149DEST_PATH_IMAGE034
时,和I)相关的目标函数;
S114.按照升序对目标函数集内的所有元素
Figure 960088DEST_PATH_IMAGE035
进行排序，得到索引向量
Figure 175168DEST_PATH_IMAGE036
Figure 405293DEST_PATH_IMAGE037
即为所求最优解，
Figure 527969DEST_PATH_IMAGE038
表示
Figure 127316DEST_PATH_IMAGE039
在索引向量
Figure 310035DEST_PATH_IMAGE040
中取值时的
Figure 761876DEST_PATH_IMAGE041
和
Figure 891506DEST_PATH_IMAGE042
的特定解;
S115.构造解空间集合
Figure 237037DEST_PATH_IMAGE043
;
S116.搜索解空间，找到满足约束的解空间子集
Figure 731604DEST_PATH_IMAGE044
，
Figure 529795DEST_PATH_IMAGE045
为满足II)中
Figure 869641DEST_PATH_IMAGE005
和
Figure 69678DEST_PATH_IMAGE046
的特定解;
S117.对于子集
Figure 968102DEST_PATH_IMAGE044
，计算目标函数集合
Figure 253589DEST_PATH_IMAGE047
，
Figure 397126DEST_PATH_IMAGE048
为取特定值
Figure 451670DEST_PATH_IMAGE049
时的和II)相关的目标函数;
S118.按照升序对目标函数集内的所有元素
Figure 553618DEST_PATH_IMAGE050
进行排序，得到索引向量
Figure 326402DEST_PATH_IMAGE051
Figure 335946DEST_PATH_IMAGE052
即为
Figure 385942DEST_PATH_IMAGE053
所求的最优解，
Figure 517846DEST_PATH_IMAGE054
表示
Figure 169405DEST_PATH_IMAGE055
在索引向量
Figure 779378DEST_PATH_IMAGE056
中取值时的
Figure 11776DEST_PATH_IMAGE057
和
Figure 455527DEST_PATH_IMAGE046
的特定解；
S12.计算多相寄存器组单元所需的相位权重因子如下：
多相寄存器组单元通道相位权重因子计算
第一级多相分解频谱计算的通道相位因子定义为：
Figure 937324DEST_PATH_IMAGE058
，
类似可定义第二级多相分解频谱计算的通道相位因子：
Figure 554250DEST_PATH_IMAGE059
。
4.根据权利要求1-3任一项所述的非2次幂特殊点数频谱计算方法，其特征在于，步骤S2中，多相寄存器组单元接收系统参数计算单元计算得到的P、Q参数值，形成相应寄存器组，缓存多相分解后的采样数据，作为频谱计算单元的输入；多相寄存器组单元包括第一级多相寄存器组和第二级多相寄存器组，其长度分别为
Figure 313259DEST_PATH_IMAGE060
和
Figure 52545DEST_PATH_IMAGE061
，
Figure 897004DEST_PATH_IMAGE062
满足“最优频谱分辨率”，输入信号采样率为
Figure 848779DEST_PATH_IMAGE063
；将输入信号进行多相分解后载入到第一级多相寄存器组，第一级多相寄存器组后的第一级频谱计算单元将实现频谱通道化计算和通道选择，通道分辨率为
Figure 55770DEST_PATH_IMAGE064
；第二级多相寄存器组及第二级频谱计算单元在通道选择后实现精细频谱计算，精细频谱分辨率为
Figure 339858DEST_PATH_IMAGE065
。
5.根据权利要求4所述的非2次幂特殊点数频谱计算方法，其特征在于，步骤S2中，多相寄存器组单元对采样数据进行多相分解计算频谱的方式包括以下步骤：
S21. 多相分解计算频谱：对长度为M的窗函数截断后的采样序列作频谱计算，将M分解为
Figure 796248DEST_PATH_IMAGE066
，
Figure 754976DEST_PATH_IMAGE067
为每一多相分支上的快拍数tap，
Figure 488577DEST_PATH_IMAGE068
为多相分解后需计算的频谱通道数；对于N点多相分解后频谱计算来讲，输出的N个频率格bin和直接作M点频谱计算后输出的M个bin之间存在以下抽取关系：
对于M个bin，从首个bin开始，每
Figure 38507DEST_PATH_IMAGE069
个bin抽取1个，得到
Figure 123138DEST_PATH_IMAGE068
个抽取后的bin；从bin[1]开始抽取，然后每
Figure 885557DEST_PATH_IMAGE069
个bin抽取1个，得到另一组
Figure 598298DEST_PATH_IMAGE068
个抽取后的bin；依次类推，如果从bin[2],…bin[
Figure 460075DEST_PATH_IMAGE070
]分别开始抽取，最终可得到
Figure 891057DEST_PATH_IMAGE066
个抽取后的bin；
在用多相分解频谱计算时，需要引入额外的相位因子实现抽取初始位置的改变，从而使
Figure 191588DEST_PATH_IMAGE071
个bin和
Figure 398316DEST_PATH_IMAGE072
个bin完全对应；
Figure 555628DEST_PATH_IMAGE073
为采样序列，
Figure 349272DEST_PATH_IMAGE074
为
Figure 453494DEST_PATH_IMAGE075
点截断窗函数，计算
Figure 875248DEST_PATH_IMAGE076
的
Figure 344406DEST_PATH_IMAGE072
点频谱可作多 相分解如下：
Figure 484401DEST_PATH_IMAGE077
Figure 392314DEST_PATH_IMAGE078
令，
Figure 809520DEST_PATH_IMAGE079
,
Figure 574214DEST_PATH_IMAGE080
；仅当
Figure 575405DEST_PATH_IMAGE081
时，
Figure 21430DEST_PATH_IMAGE082
多相分解计算后的
Figure 417777DEST_PATH_IMAGE068
点频谱
Figure 228738DEST_PATH_IMAGE083
是
Figure 343324DEST_PATH_IMAGE084
点原始频谱
Figure 327461DEST_PATH_IMAGE085
的
Figure 984838DEST_PATH_IMAGE067
倍抽取，抽取初始位置为
Figure 825755DEST_PATH_IMAGE086
；
Figure 630900DEST_PATH_IMAGE087
和
Figure 90832DEST_PATH_IMAGE088
对应着
Figure 461770DEST_PATH_IMAGE089
和
Figure 118928DEST_PATH_IMAGE090
；
定义抽取索引数
Figure 942528DEST_PATH_IMAGE091
，
Figure 534046DEST_PATH_IMAGE092
对应着抽取初始位置为
Figure 634857DEST_PATH_IMAGE093
的多相分解计算频谱bin；
Figure 83156DEST_PATH_IMAGE094
相对于
Figure 534997DEST_PATH_IMAGE089
的计算，
Figure 461365DEST_PATH_IMAGE095
需要在第r个快拍tap上施加
Figure 744579DEST_PATH_IMAGE096
和第
Figure 239145DEST_PATH_IMAGE097
个通道上施加
Figure 302916DEST_PATH_IMAGE098
；
Figure 406876DEST_PATH_IMAGE099
可以被吸收到窗函数的多相分支
Figure 544596DEST_PATH_IMAGE100
中，令
Figure 334698DEST_PATH_IMAGE101
第
Figure 761131DEST_PATH_IMAGE097
个通道上的因子定义为：
Figure 29301DEST_PATH_IMAGE102
。
6.根据权利要求5所述的非2次幂特殊点数频谱计算方法，其特征在于，步骤S2中，还包括以下步骤：
S22. 修正窗函数值的计算；
多相分解后，窗函数被修正为：
Figure 287107DEST_PATH_IMAGE103
对于第一级多相分解频谱计算，其窗函数的多相分解形式为
Figure 123476DEST_PATH_IMAGE104
，修正后的窗函数多相分解形式为
Figure 270162DEST_PATH_IMAGE105
，
Figure 217389DEST_PATH_IMAGE106
为频率通道序号；
对于第二级多相分解频谱，其窗函数的多相分解形式为
Figure 392018DEST_PATH_IMAGE107
，修正后的窗函数多相分解式为：
Figure 461606DEST_PATH_IMAGE108
,
Figure 597052DEST_PATH_IMAGE109
为精细频谱通道序号。
7.根据权利要求6所述的非2次幂特殊点数频谱计算方法，其特征在于，步骤S3中，频谱计算单元接收来自系统参数计算单元的
Figure 472604DEST_PATH_IMAGE110
参数值和相位权重因子、来自多相寄存器组单元的多相分解采样数据，频谱计算单元中的FFT模块对相位加权后的多相信号进行频谱计算，频谱计算的结果被输入到功率谱计算单元中，最终得到通道功率和精细频谱功率。
8.根据权利要求7所述的非2次幂特殊点数频谱计算方法，其特征在于，根据约束目标的不同，如“最小计算延迟”和“最少计算资源”，频谱计算单元可做灵活伸缩修改，采用全并行频谱计算实现，或者时分复用计算实现。
9.一种非2次幂特殊点数频谱计算系统，其特征在于，所述系统用于实现根据权利要求1-8任一项所述的非2次幂特殊点数频谱计算方法，系统包括系统参数计算单元、多相寄存器组A、频谱计算单元A、多相寄存器组B、频谱计算单元B、通道选择单元和功率谱计算单元。
10.根据权利要求9所述的非2次幂特殊点数频谱计算系统，其特征在于，采样信号首先被输入到多相寄存器组A中，在信号通路上依次经过频谱计算单元A、通道选择单元、多相寄存器组B、频谱计算单元B，最终由功率谱计算单元输出通道功率和精细频谱功率。
CN202211562557.5A 2022-12-07 2022-12-07 一种非2次幂特殊点数频谱计算方法及系统 Active CN115587288B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211562557.5A CN115587288B (zh) 2022-12-07 2022-12-07 一种非2次幂特殊点数频谱计算方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211562557.5A CN115587288B (zh) 2022-12-07 2022-12-07 一种非2次幂特殊点数频谱计算方法及系统

Publications (2)

Publication Number Publication Date
CN115587288A true CN115587288A (zh) 2023-01-10
CN115587288B CN115587288B (zh) 2023-04-28

Family

ID=84783125

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211562557.5A Active CN115587288B (zh) 2022-12-07 2022-12-07 一种非2次幂特殊点数频谱计算方法及系统

Country Status (1)

Country Link
CN (1) CN115587288B (zh)

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7027942B1 (en) * 2004-10-26 2006-04-11 The Mitre Corporation Multirate spectral analyzer with adjustable time-frequency resolution
CN101707582A (zh) * 2009-11-05 2010-05-12 东南大学 基于多相分解的多天线信道估计方法
CN103020015A (zh) * 2012-11-30 2013-04-03 桂林卡尔曼通信技术有限公司 点数为非2次幂的离散傅里叶变换快速计算的实现方法
CN107976572A (zh) * 2017-12-18 2018-05-01 中国科学院深圳先进技术研究院 一种信道化高精度实时频谱分析方法及其系统
CN111414257A (zh) * 2020-03-30 2020-07-14 中国科学院新疆天文台 一种脉冲星信号多通道滤波方法、装置及存储介质

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7027942B1 (en) * 2004-10-26 2006-04-11 The Mitre Corporation Multirate spectral analyzer with adjustable time-frequency resolution
CN101707582A (zh) * 2009-11-05 2010-05-12 东南大学 基于多相分解的多天线信道估计方法
CN103020015A (zh) * 2012-11-30 2013-04-03 桂林卡尔曼通信技术有限公司 点数为非2次幂的离散傅里叶变换快速计算的实现方法
CN107976572A (zh) * 2017-12-18 2018-05-01 中国科学院深圳先进技术研究院 一种信道化高精度实时频谱分析方法及其系统
CN111414257A (zh) * 2020-03-30 2020-07-14 中国科学院新疆天文台 一种脉冲星信号多通道滤波方法、装置及存储介质

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
HASSAN ABOUSHADY等: "Efficient Polyphase Decomposition of Comb Decimation Filters in ∑△ Analog-to-Digital Converters", 《IEEE TRANSACTIONS ON CIRCUITS AND SYSTEMS—II: ANALOG AND DIGITAL SIGNAL PROCESSING》 *
沈梦萍等: "基于图形处理器的转变边缘传感器读取系统信号处理技术研究", 《北京师范大学学报（自然科学版）》 *

Also Published As

Publication number Publication date
CN115587288B (zh) 2023-04-28

Similar Documents

Publication Publication Date Title
US8107747B2 (en) Codevector trainers and real-time compressor for multi-dimensional data
Potluri et al. Multiplier-free DCT approximations for RF multi-beam digital aperture-array space imaging and directional sensing
CN103973324A (zh) 一种宽带数字接收机及其实时频谱处理方法
CN107015063B (zh) 宽带多信道数字相关接收机及接收方法
Golwala et al. Status of MUSIC, the multiwavelength sub/millimeter inductance camera
Mao et al. R-matrix electron-impact excitation data for the C-like iso-electronic sequence
CN112067129A (zh) 一种高光谱处理方法及波段选择方法
Hamini et al. ORFEES–a radio spectrograph for the study of solar radio bursts and space weather applications
Zhang et al. A new multichannel parallel real-time FFT algorithm for a solar radio observation system based on FPGA
Pereira et al. A low-cost hardware accelerator for ccsds 123 predictor in fpga
CN115587288A (zh) 一种非2次幂特殊点数频谱计算方法及系统
Chen et al. FPGA implementation of neural network accelerator for pulse information extraction in high energy physics
Gunst et al. The LOFAR phased array telescope system
Mao et al. The electron energy loss rate due to radiative recombination
Langevin et al. Image and spectral image compression for four experiments on the ROSETTA and Mars Express missions of ESA
Bruno et al. Rotational and vibrational temperatures of hydrogen nonequilibrium plasmas from Fulcher band emission spectra
Yang et al. FPGA optimization for hyperspectral target detection with collaborative representation
CN102118202A (zh) 一种实现宽带信号分析的软件无线电多信道处理方法
Chikada Correlators for interferometry–Today and Tomorrow
De Vos et al. Cluster computers and grid processing in the first radio-telescope of a new generation
Barry et al. Large-format, transmission-line-coupled kinetic inductance detector arrays for hep at millimeter wavelengths
Bergano et al. Digital complex correlator for a C-band polarimetry survey
Kalinovsky et al. Reversible Digital Filters Total Parametric Sensitivity Optimization using Non-canonical Hypercomplex Number Systems
RU2726281C1 (ru) Активная фазированная антенная решетка
CN219304828U (zh) 一种射电天文相位阵馈源时域波束合成终端

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant