CN111033499B

CN111033499B - 用于在O(n log n)的时间和O(n)的存储中对线性调频Z变换求逆的系统和方法

Info

Publication number: CN111033499B
Application number: CN201880050406.7A
Authority: CN
Inventors: 沃洛迪米尔·苏霍伊; 亚历山大·施托伊切
Original assignee: Iowa State University Research Fund Co ltd
Current assignee: Iowa State University Research Fund Co ltd
Priority date: 2017-07-24
Filing date: 2018-07-24
Publication date: 2022-07-29
Anticipated expiration: 2038-07-24
Also published as: CN111033499A; US20210141856A2; US20210397674A1; KR20200023486A; JP6928208B2; KR102183973B1; KR20200133283A; US20200159808A1; KR102338456B1; EP3659052A4; JP2020526852A; WO2019023220A1; EP3659052A1

Abstract

本公开涉及一种用于在O(n log n)的时间和O(n)的存储中对线性调频Z变换求逆的系统和方法。本公开的实施例描述了一种高效的O(n log n)方法，该方法实现了逆线性调频Z变换(ICZT)。该变换是公知的前向线性调频Z变换(CZT)的逆变换，该变换通过允许采样点落在对数螺旋轮廓上而不是单位圆上来泛化快速傅立叶变换(FFT)。因此，ICZT可以看作是快速傅里叶逆变换(IFFT)的泛化。

Description

用于在O(n log n)的时间和O(n)的存储中对线性调频Z变换求逆的系统和方法

技术领域

本发明总体上涉及逆线性调频Z变换，更具体地涉及用于在O(n log n)的时间和O(n)的存储中计算逆线性调频Z(ICZT)变换的系统和方法。

背景技术

线性调频Z变换(CZT)通过允许样本位于对数螺旋轮廓而不是单位圆上来扩展离散傅里叶变换(DFT)。更具体地说，该变换将样本沿着由公式A W^-k定义的对数螺旋轮廓分布，其中k＝0,1,2,…,M-1。非零复数A和W指定螺旋轮廓的位置和方向，以及轮廓上采样点的间距。更具体地说，给定N-元素的输入向量x，则CZT计算M-元素的输出向量X，其中，X的第k个元素由

给出。

先前已经描述了可以计算正向线性调频Z变换的有效算法。该算法运行时间为O(nlog n)，其中n是变换的大小。因为输入的数量N和输出的数量M可以不同，所以在最一般的情况下，CZT算法的计算复杂度由n＝max(M,N)确定。

可以使用索引替换得出高效的CZT算法，该索引替换最初是由Bluestein在“用于计算离散傅里叶变换的线性滤波方法”(A linear filtering approach to thecomputation of discrete Fourier transform)，IEEE Transactions on Audio andElectroacoustics，18(4)：451-455，1970)中提出的，通过引用将其整体并入本文，以使用快速卷积来表达变换。目前已经为CZT算法提出了各种有用的优化。但是，其计算复杂度保持固定为O(n log n)。

逆线性调频Z变换(ICZT)是线性调频Z变换(CZT)的逆变换。也就是说，ICZT将CZT的输出映射回输入。由于CZT是线性变换，因此可以使用CZT矩阵与输入向量的矩阵乘积表示。可以使用标准算法对该矩阵进行求逆。但是，在算法形式上，此过程可能需要多达O(n³)个操作。

即使有运行比O(n³)快的高效矩阵求逆算法，也需要至少n²个运算才能访问n×n矩阵的每个元素。因此，O(n²)是任何与存储器中以完整的n×n矩阵工作的ICZT算法的计算复杂度的下限。

为提高效率，ICZT算法必须具有与CZT算法相同的计算复杂度，即O(n log n)。该要求排除了需要将转换矩阵存储在存储器中的任何方法。

目前已经进行了得出高效ICZT算法的几种尝试。在一种情况下(描述参见Frickey，“使用逆线性调频Z变换进行仿真雷达信号的时域分析”，技术报告，爱达荷州国家工程实验室，爱达荷州，ID，1995年，通过引用将其全部内容并入本文)，将正向CZT算法的修改版本描述为ICZT算法，其中对数螺旋轮廓沿相反方向移动。但是，这种方法并不能真正对CZT求逆。它仅在某些特殊情况下有效，例如，当A＝1且W＝e^-2πi/n时。也就是说，在CZT简化为DFT的情况下。在一般情况下，即当

时，此方法生成不对CZT求逆的变换。

发明内容

本公开的实施例描述了一种用于计算ICZT的有效的O(n log n)方法。该方法使用生成矢量来在O(n)存储器中更紧凑地表示结构化矩阵(例如对角矩阵、范德蒙(Vandermonde)矩阵、托布里兹(Toeplitz)矩阵、循环矩阵或斜循环矩阵)。ICZT可以看作快速傅里叶逆变换(IFFT)脱离单位圆周的泛化。换句话说，ICZT的复数参数A和W描述了复数平面中的对数螺旋轮廓。与IFFT不同，ICZT可以使用单位圆外的轮廓点，这些轮廓点对应于指数增长或指数衰减的频率分量。本公开的实施例还描述了一种用于在|W|<1的情况下提高CZT或ICZT的数值精度的方法，其在附录A中给出。

据发明人所知，其相信本公开是对有效ICZT算法的首次描述。特别地，工作算法及其推导在此被描述。此外，使用自动测试案例可以验证该算法的准确性。

该算法是通过使用结构化矩阵乘法表达CZT公式然后找到一种将该公式求逆的方法而得出的。ICZT计算的本质简化为对被特殊构建的范德蒙矩阵W进行求逆。此问题又简化为对从W得出的特殊构建的托布里兹矩阵

求逆。

结合附图时，根据以下详细描述，本发明的其他方面、目的和优点将变得更加明显。

附图说明

结合在说明书中并构成说明书一部分的附图示出了本发明的几个方面，并且与说明书一起用于解释本发明的原理。在附图中：

图1A-1B描绘了分别具有32和64个点的对数螺旋轮廓的两个示例，其中每个轮廓的起点用未填充的圆表示；

图2是这样的图表，其提供对用于针对形状类似于图1A-1B中所示但具有不同数量M的点、且针对使用不同浮点精度的计算而执行CZT接着是ICZT的数值精度的说明；

图3是提供对针对不同类型的结构化矩阵的生成矢量和矩阵形状的说明的图表；

图4A-4B给出了一个示例，其分别说明了FFT或IFFT使用的固定幅度复指数与CZT或ICZT可以使用的按指数增长或指数衰减的复指数之间的差异。

图5A-5D示出了对数螺旋轮廓的示例，其点分别位于单位圆外、单位圆内、在单位圆外开始并且在其内结束、以及在单位圆内开始并且在其外结束；

图6A-6C示出了对数螺旋轮廓的示例，其点分别位于跨越90度、180度和360度的单位圆上；

图7A-7B描绘了对数螺旋轮廓的示例，其点分别跨越两个360度旋转和五个360度旋转。

尽管将结合某些优选实施例描述本发明，但是并非旨在于将其限制于那些实施例。相反，本发明的意图是覆盖由所附权利要求书所限定的本发明的精神和范围内所包括的所有替换、修改和等同形式。

具体实施方式

本公开的实施例描述了一种有效的O(n log n)ICZT算法，该算法对于

实现了ICZT。ICZT通过使得可以将采样点分布在复平面中的对数螺旋轮廓上而不是单位圆上，来泛化快速傅里叶逆变换(IFFT)。在图1A-1B中示出了两个示例性的对数螺旋轮廓。只要采样点不在单位圆上，相应的频率分量就会呈指数增长或衰减。图4A-4B示出了由FFT/IFFT和CZT/ICZT使用的增长、衰减和恒定幅度的频率分量的一些示例。

通过将CZT公式表示为结构化矩阵的乘积并找到一种将矩阵方程求逆的方法，得出了ICZT算法。更具体地，计算ICZT简化为对范德蒙矩阵W的特例求逆。这可以通过对从W导出的托布里兹矩阵

的特例求逆来实现。

结构化矩阵是可以用比其元素总数少的参数来描述的矩阵。结构化矩阵的示例包括托布里兹、Hankel、范德蒙和Cauchy矩阵。其他示例包括循环、斜循环和DFT矩阵。对角矩阵也可以视为结构化矩阵。图3示出了这些结构化矩阵中的一些的通用形状。图3中所示的示例为3×3矩阵。此外，图3示出了可用于生成每个矩阵的(多个)生成向量。在大多数情况下，这些用r和c表示，其指定矩阵的第一行和第一列。

Gohberg和Semencul表明，托布里兹矩阵的求逆可以使用托布里兹矩阵的两个乘积之差来表示。这项工作可以在I.Gohberg等人的文章中找到，“关于有限托布里兹矩阵及其连续类似物的求逆”(On the inversion of finite Toeplitz matrices and theircontinuous analogs)，Mat.issled，2(1)：201-233，1972以及William Trench，“有限托布里兹矩阵求逆的算法”(An algorithm for inversion of finite Toeplitz matrices)，工业与应用数学学会学报，12(3)：515-522，1964，这两篇参考文献都通过引用整体并入本文。在文献中，此结果通常称为“Gohberg-Semencul公式”。该公式中的四个矩阵(上三角或下三角托布里兹)是由两个向量u和v生成的。对于ICZT的情况，要求逆的托布里兹矩阵也是对称的。这导致了Gohberg-Semencul公式的简化版本，其仅使用一个生成向量，即仅使用向量u，来表达求逆。

发明人已经确定，在ICZT的情况下，可以导出一个公式，该公式将生成矢量u定义为复数W的函数。该公式导致了高效的ICZT算法。该算法包括将托布里兹矩阵乘以向量的步骤，这可以在O(n log n)中完成，而不将完整的托布里兹矩阵存储在存储器中。附录C、D和E基于可以计算该乘积的这些引用实现了三种不同的算法。这些算法中的每一个都可以用作ICZT算法的构成部分。这三种算法中的每一种都使用附录B中陈述的快速傅立叶变换(FFT)算法。

广泛影响

离散傅里叶变换(DFT)及其使用快速傅里叶变换(FFT)的高效实现被深深地嵌入到各种各样的现代应用中。因为CZT是DFT的泛化，而ICZT是逆DFT的泛化，所以本发明的潜在应用数量非常大。到目前为止，只有CZT算法具有与FFT相同的计算复杂度，即O(n logn)。本公开介绍了一种ICZT算法，其具有与逆FFT(IFFT)相同的复杂度，其也是O(n log n)。

换言之，本发明是变革性的，不仅因为它实现了泛化逆FFT的变换，而且还因为该新算法具有与其泛化的算法相同的运行时间复杂度。

可以与ICZT一同使用的对数螺旋轮廓包括轮廓，这些轮廓的点位于单位圆外部、单位圆内部、在单位圆外部开始并在其内部结束、在单位圆内部开始并在其外部结束。这样的轮廓的示例在图5A-5D中示出。此外，轮廓点可以覆盖单位圆上的部分弧，例如90度弧、180度甚至是完整的360度旋转，在这种情况下该轮廓可以与FFT/IFFT使用的轮廓相同。这样的轮廓的示例在图6A-6C中示出。最后，轮廓点甚至可以覆盖多个360度旋转。图7A-7B示出了具有两个和五个完整旋转的示例。通常，转数甚至可能不是整数。

使用结构化矩阵表达线性调频Z变换

以矩阵形式表示的示例4×4CZT

该示例使用矩阵表示法描述了特殊情况下的CZT，该特殊情况具有四个输入和四个输出。在这种情况下，CZT使用以下公式定义：

令x＝(x₀,x₁,x₂,x₃)^T表示CZT输入向量，X＝(X₀,X₁,X₂,X₃)^T表示输出向量。使用这两个向量，还可以使用以下矩阵公式来表示4×4CZT：

也就是说，示例CZT可以看作是将输入向量x乘以由复数参数A的负幂生成的对角矩阵，并将所得的向量乘以由复数参数W确定的范德蒙矩阵W。

在该示例中，矩阵W具有以下形式：

对矩阵W求逆得到ICZT。

使用范德蒙矩阵W表示CZT

在一般情况下，通常使用以下公式定义CZT：

在该公式中，A和W是用于变换的复数参数，其定义了对数螺旋轮廓和样本在其上的位置。参数N是整数，其指定输入向量x的大小。最后，参数M是指定输出向量X的大小的整数。通常，N可能不等于M。也就是说，输入的维数可能不等于输出的维数。

CZT还可以以矩阵形式表示：

公式(2.5)可以以以下缩短形式来重述：

X＝W diag(A^-0，A^-1，A^-2，...，A^-(N-1))x， (2.6)

其中，W是如下M×N矩阵：

可以看出，W是范德蒙矩阵。即，W的每一行形成几何级数。此外，这些级数中每个级数的共同比等于参数W的相应整数幂。

公式(2.6)中A的负整数幂对范德蒙矩阵W的列进行缩放。因为矩阵W是范德蒙矩阵的特例，所以它可以表示为对角矩阵和托布里兹矩阵

的乘积，如下所述。

从范德蒙矩阵W得出托布里兹矩阵

由于矩阵W是范德蒙矩阵的特例，所以其可以表示为两个对角矩阵和托布里兹矩阵

的乘积。可以使用以下公式表达在矩阵W的每个元素中的参数W的幂：

公式(2.8)表示(2.4)的右侧可以被表示为如下：

可以重新排列最后一个公式中的项，以便可以更轻松地将其映射到矩阵乘积：

项

映射到M×M对角矩阵P。项

映射到两个N×N对角矩阵Q和A的乘积。最后，项

映射到M×N托布里兹矩阵

如下所示：

这意味着可以使用实现该操作的几种O(n log n)算法之一来有效地计算矩阵

与矢量的乘积(请参阅附录C、D和E)。换句话说，公式(2.10)具有以下矩阵形式：

其中所有矩阵-向量乘积都可以用O(n log n)计算，其中n＝max(M,N)。

综上，CZT算法可以看作是以下矩阵方程的有效实现：

其中

A＝diag(A^-0，A^-1，...，A^-(N-1))，W是在(2.7)中定义的范德蒙矩阵的特例，以及

是在(2.11)中定义的托布里兹矩阵的特例。如前所述，x是对于CZT的输入向量，X是来自CZT的输出向量。

表示托布里兹矩阵的求逆

用于4×4情况的Gohberg-Semencul公式

令T为4×4托布里兹矩阵，即

其中a,b,c,d,f,g和h是生成矩阵T的七个复数。

Gohberg-Semencul公式指出，可以使用以下公式对非奇异的4×4托布里兹矩阵T进行求逆：

其中u＝(u₀,u₁,u₂,u₃)是一个四元素矢量，使得u₀≠0，以及v＝(v₀,v₁,v₂,v₃)是另一个四元素矢量，使得v₃≠0。这两个向量由生成矩阵T的数字a,b,c,d,f,g和h确定。但是，很难将它们明确表示为生成T的七个数字的函数。

综上所述，公式(3.2)使用四个结构化矩阵表示4×4托布里兹矩阵T的求逆：1)由向量u生成的下三角托布里兹矩阵

2)由向量v的反转生成的上三角托布里兹矩阵

3)由向量(0，v₀，v₁，v₂)产生的下三角托布里兹矩阵

该向量通过将v向右移动一个元素而获得，4)由向量(0，u₃，u₂，u₁)生成的上三角托布里兹矩阵

该向量通过将u的反转向右移动一个元素而获得。

Gohberg-Semencul公式的一般情况

令T是由向量r＝(r₀,r₁,r₂,...,r_n-1)和向量c＝(c₀,c₁,c₂,...,c_n-1)生成的n×n托布里兹矩阵，向量r＝(r₀,r₁,r₂,...,r_n-1)指定了它的第一行，而向量c＝(c₀,c₁,c₂,...,c_n-1)指定其第一列。假设r₀＝c₀。更规范地，

托布里兹矩阵T的逆矩阵T^-1可以不是托布里兹矩阵。尽管如此，Gohberg-Semencul公式仍可以使用上三角和下三角托布里兹矩阵表达T^-1。

也就是说，Gohberg-Semencul公式指出存在向量u＝(u₀,u₁,u₂,...,u_n-1)和向量v＝(v₀,v₁,v₂,...,v_n-1)满足以下矩阵方程：

换言之，逆矩阵T^-1可以看作是由向量u和向量v生成的结构化矩阵。此类推延伸到将T^-1乘以向量。如果确定了生成向量u和v，则可以通过使用结构化矩阵乘法的实现(3.4)，在O(n log n)中计算矩阵T^-1与向量的乘积。

使用结构化矩阵表达ICZT

公式(2.13)表明，可以将CZT简化为矩阵-向量乘法的序列，其中除

之外的所有矩阵都为对角，

是如下托布里兹矩阵：

这意味着可以将ICZT视为对角矩阵，逆矩阵

和向量的乘积(参见公式(2.13))。

令u和v为使用Gohberg-Semencul公式生成

的两个向量，即，

u＝(u₀,u₁,u₂,…,u_n-1)^T, (4.2)

v＝(v₀,v₁,v₂,…,v_n-1)^T. (4.3)

由于矩阵

是对称的，所以其逆

也是对称的。

的公式是Gohberg-Semencul公式的特例。针对n的前几个值检查逆矩阵

的特例，建议做出三个有用的假设。首先，向量u等于

的第一列。更规范地，

对于每个j∈{1，2，...，n}. (4.4)

其次，向量v等于W^-1的最后一列，即

对于每个j∈{1，2，...，n}. (4.5)

第三，向量v等于向量u的反转，

v_k＝u_n-1-k 对于每个k∈{0，1，2，...，n-1}. (4.6)

换言之，仅需要确定向量u就能使用结构化矩阵乘法将逆矩阵

与向量相乘。

一种推导用于计算向量u的显式公式的方法是从n的前几个值的特例推导出向量u。如果这种归纳成功，那么应该可以对每个有限n使用所得的公式。由于u确定v，所以足以使用Gohberg-Semencul公式和结构化矩阵乘法来公式化ICZT算法。

用于对托布里兹矩阵

求逆的公式的特例

n＝1的特例

假设n＝1。那么，

是一个由等于1的单一数字组成的1×1矩阵，更规范地，

这意味着逆矩阵

也是由单个数字1组成的1×1矩阵，即

的第一行和最后一行是相同的：都由等于1的单个元素组成，更规范地，

这导致Gohberg-Semencul公式的退化特例：

n＝2的特例

假设n＝2。那么，

是由以下公式指定的2×2矩阵：

对于每个

其定义是明确的。

在这种情况下，逆矩阵

也是2×2矩阵，这意味着其第一列u是二元素向量，即u＝(u₀,u₁)^T。通过求解以下线性系统，其元素可以表示为W的函数：

在该公式中，e₀表示一个向量，其中第一个元素设置为1，其余元素设置为0，即

e₀＝(1,0)^T. (4.14)

换言之，

的第一行与向量u的点积必须等于1，而

的第二行与u的点积必须等于0。

线性系统(4.13)可以用

的四个元素和u的两个元素表示。这导致(4.13)的以下扩展形式：

先前的方程式也可以看作是两个标量线性方程式的系统：

只要当W≠0时，这两个方程都是定义明确的。

可以从(4.17)的两侧减去项

这导致了u₁的值作为W和u₀的函数的公式：

先前等式的右侧可以插入到(4.16)中。这样得出的方程式仅取决于u₀和W。更规范地，

通过将前一个方程式的两边乘以

得出u₀值作为W的函数的公式：

将该值插入到(4.18)的右侧，得到将u₁的值表示为W的函数的公式。

更规范地，

总之，在2×2的情况下，指定逆矩阵

的第一列的列向量u可以表示为参数W的以下函数：

向量v可以通过反转u的元素来获得，即

可以通过组合u和v的表达式来得出逆矩阵

的封闭型的表达式。

在该特例下，可以使用以下矩阵变换序列来验证Gohberg-Semencul公式。它从(3.4)中所述公式的左侧和右侧之间的差开始，并表明该差等于零：

n＝3的特例

在此特例下，

是3×3托布里兹矩阵。

的每个元素都是变换参数W的函数。更规范地，

逆矩阵

也是3×3矩阵，但是与2×2情况相比，它不再是托布里兹。设u为

的第一列，即

与2×2情况类似，向量u是以下线性系统的解：

然而，与先前的情况相反，在这种情况下，e₀表示由三个元素而不是两个元素组成的列向量。e₀的第一个元素等于1。e₀的其余两个元素等于0。

换言之，

e₀＝(1,0,0)^T. (4.29)

线性系统(4.28)可以用构成

的九个元素和构成u的三个元素来表示。展开的结果如下所示：

该线性系统还可以使用具有三个未知数的三个方程的系统以标量形式表示。这三个方程中每个方程的左侧等于

的行与向量u之间的点积。右侧等于e₀的相应元素。更规范地，

求解该线性方程系统得到一个公式，该公式将向量u表示为参数W的函数，如下所示：

逆矩阵

具有以下形式

注意，向量v指定

的最后一列，并且其等于向量u的反转，向量u指定矩阵的第一列。换言之，先前的公式确认了三个假设(4.4)、(4.5)和(4.6)。

对n的其他值重复该过程允许我们导出确定u的通用公式。下一部分将说明此公式。

用于对托布里兹矩阵

求逆的通用公式

在一般情况下，可以使用以下公式来计算指定

的第一列的向量u的元素：

对于每个k，可以使用出现在(4.36)右边分母中的多项式乘积的预先计算值，在O(1)中计算u_k的值。这些值在一个单独的循环中计算，该循环以O(n)运行。

可以将u的结果值及其等于v的反转插入Gohberg-Semencul公式(3.4)中，以计算逆矩阵

其也可以用作计算托布里兹矩阵和向量的乘积的有效算法的生成向量。可以按实现将矩阵

与向量相乘而不实际将

存储在存储器中的顺序来调用这些算法。这导致ICZT算法以O(n log n)运行。

矩阵形式的ICZT

以下等式以矩阵形式显示了ICZT：

其中A^-1＝diag(A⁰，A¹，...，A^N-1),

以及

即，该公式示出了如何从其输出矢量X计算CZT输入矢量x。该方程式是通过将前向变换公式(2.13)中的对角矩阵求逆并将矩阵

替换为其逆

得出的。

公式(4.37)中的每个矩阵都是对角矩阵，或者可以表示为托布里兹矩阵的乘积之差。实施该公式将得出用于计算ICZT的O(n log n)算法，该算法在下面的算法5.2中给出。该算法通过使用生成向量以计算所有矩阵-向量运算的结果而不在存储器中存储任何中间矩阵，从而高效地执行所有矩阵计算。其需要O(n)存储器来运行。

示例

在3×3的情况下，CZT可以用以下矩阵方程表示：

在那种情况下，用于ICZT的扩展矩阵方程看来为：

在(4.35)中给出逆矩阵

但是，该算法未明确构造此矩阵。如上所述，

可以表示如下：

其中向量u＝(u₀,u₁,u₂)由(4.34)给出，向量v＝(v₀,v₁,v₂)＝(u₂,u₁,u₀)是u的反转。对于u和v的这些值，方程式(4.40)具有以下形式：

另外，矩阵

等于矩阵

的转置，并且矩阵

等于矩阵

的转置。因此，逆矩阵

可以用以下简化形式表示：

此外，

和

都可以从向量u的元素生成。

将(4.42)代入(4.37)得到ICZT的以下矩阵方程式：

相反，通过利用这些矩阵的结构可以有效地计算方程式(4.43)。如上所述，这四个托布里兹矩阵并未明确构建。

CZT算法和ICZT算法

算法5.1给出了用于CZT算法的伪代码。该算法使用附录E中所述的基于卷积的算法TOEPLITZMULTIPLYC，以将托布里兹矩阵与矢量相乘。可以替代地使用其他算法来计算托布里兹矩阵和向量的乘积。例如，附录C中描述的TOEPLITZMULTIPLYE和附录D中描述的TOEPLITZMULTIPLYP都可以替换算法5.1中的TOEPLITZMULTIPLYC，而不更改其计算复杂性。

算法5.2给出了ICZT算法的伪代码。其实现了逆公式(4.42)，而没有实际在存储器中存储任何矩阵。该算法以O(n log n)运行，其中n＝max(M,N)。该实现仅支持n＝M＝N的情况。如上所述，在第23-26行使用的函数TOEPLITZMULTIPLYC，可以用TOEPLITZMULTIPLYE或TOEPLITZMULTIPLYP代替，而不更改算法的计算复杂性。

图2给出了ICZT的数值精度，其被作为变换大小M和在计算过程中浮点数使用的位数的函数而给出。图1A-1B中的两个对数螺旋轮廓对应于图2中的M＝32和M＝64。对于图2的所有行，将A的值设置为1.1。使用公式

确定W的值。通过执行CZT然后执行ICZT，并计算CZT的输入与ICZT的输出之间的欧几里德距离，计算出数值精度。最终精度对100个随机输入向量求平均。图2中的第二列示出了变换矩阵的条件数。尽管条件数很高，但是当使用高精度浮点数时，即使对于较大的M值，也可以解决此问题。

算法5.3给出了ICZT-RECT算法的伪代码，该算法为N≤M的情况实现了逆算法。该算法调用算法5.2，其中，将第二参数设置为M然后仅返回其输出向量的前N个元素。

已经描述了ICZT算法，现在提供该算法的应用。在这方面，以下只是ICZT算法的一些有用应用的部分列表。在实施例中，ICZT算法在用于生成或分析从基于雷达的传感器接收或发送到基于雷达的传感器的数据的应用中使用。在其他实施例中，ICZT算法在用于生成或分析从声纳传感器接收或发送到声纳传感器的数据的应用中使用。在另外其他实施例中，ICZT算法在用于生成或分析从其他基于范围的传感器接收或发送到其他基于范围的传感器的数据的应用中使用。

在其他实施例中，ICZT算法用于ICZT算法的硬件实现。在实施例中，硬件实现是各种系统，包括存储设备和处理器，例如个人计算机、芯片(例如，片上系统)、包括专用集成电路(ASIC)的系统、包括图形处理单元(GPU)或现场可编程门阵列(FPGA)的系统。

在进一步的实施例中，ICZT算法用于通过替换较慢的方法来提高计算ICZT需要的应用程序的速度，这些较慢的方法对矩阵显式求逆或求解显式公式化的线性系统。另外，ICZT算法可以全部或部分地替代常规使用的FFT和IFFT。这种替换可以扩展计算系统的范围，使其覆盖可能位于单位圆内或外的对数螺旋轮廓上的点。更进一步，在计算ICZT之前或期间，可以修改频域向量，以使得由ICZT返回的时域向量不同于用于生成所述频域向量的时域向量。例如，可以对频域矢量执行操作，诸如从矢量中滤除某些元素或频率。

此外，在实施例中，ICZT算法用于医学成像应用中，例如CT、PET或MRI。此外，ICZT算法还用于信号处理、信号分析、信号合成、语音和音频处理、以及图像处理的应用中。在信号处理应用中，CZT或FFT的输入向量可以是从音频或图像信号(例如通过对信号进行采样)得出的时域向量。然后，CZT或FFT的输出向量可以是那些信号的频域向量表示。

下列附录A-E提供有关ICZT算法的实现的其他信息。

附录A：反转对数螺旋轮廓的方向以提高CZT和ICZT的数值精度

本附录描述了CZT和ICZT算法的替代版本，其提高了计算出的变换的数值精度。这是通过在|W|<1时反转线性调频轮廓的方向而在|W|≥1时保持原始方向来实现的。反转轮廓的参数为W′＝W^-1和A′＝AW^-(M-1)。取决于变换参数的具体值，数值精度的提高可能超过几个数量级。

算法A.1给出了CZT算法的替代版本的伪代码，该伪代码在|W|<1时执行轮廓反转。算法A.2给出了ICZT算法的替代版本的伪代码，在这种情况下，它也执行轮廓反转。

使用这些算法，对于作为增长着的对数螺旋的线性调频轮廓，即，当|W|<1时，CZT和ICZT两者的数值精度都可以得到改善。在那些情况下，算法A.1和算法A.2反转线性调频轮廓的方向。此外，这种反转还伴随着向量X中元素顺序的反转，向量X是算法A.1中的CZT输出向量和算法A.2中的ICZT输入向量。

附录B：FFT和逆FFT算法

该附录陈述了FFT算法和逆FFT算法。算法B.1给出了用于快速傅立叶变换(FFT)的伪代码。算法B.2给出了快速傅里叶逆变换(IFFT)的伪代码。FFT和IFFT都以O(n log n)运行。

附录C：将托布里兹矩阵乘以向量

该附录描述了一种将托布里兹矩阵乘以矢量的流行方法。令A为n×n托布里兹矩阵，令x为向量。矩阵A由其第一行r和其第一列c指定，其中假定r₀等于c₀。本附录中的示例假定n是2的幂。但是，可以修改该方法以适用于n的其他值。更规范地，

为了使用基于FFT的算法高效地计算乘积y＝Ax，可以将A嵌入2n×2n循环矩阵

如下所示：

向量x的末尾填充有n个零，这导致长度为2n的向量

和

的乘积可以使用以O(n log n)运行的基于FFT的算法来计算。也就是说，目标是计算长度为2n的矢量

其等于

更规范地，

的前n个元素等于y的对应元素，即，对于每个k∈{0,1,2,…,n-1}，都有

剩下的只是示出如何使用FFT和IFFT将(2n)×(2n)循环矩阵B与长度为2n的向量v相乘。矩阵B由其第一列向量b生成。更规范地，令

可以使用DFT和逆DFT来计算该乘积，这可以根据循环卷积定理得出。换言之，乘积的每个元素(B v)_k等于b和v的循环卷积中的第k个元素。

即

DFT(B v)＝DFT(b)×DFT(v), (C.6)

其中×表示元素对元素乘法。因此，

B v＝DFT^-1(DFT(b)×DFT(v)). (C.7)

可以使用FFT和逆FFT来计算DFT和逆DFT，即

B v＝FFT^-1(FFT(b)×FFT(v)). (C.8)

总之，可以使用以下六个步骤来计算n×n托布里兹矩阵与向量的乘积，其中n被假定为2的幂：1)通过连接第一列、单个零、以及第一行中最后n-1个元素的反转，从托布里兹变为循环矩阵；2)在向量的末尾填充n个零；3)计算循环矩阵第一列的DFT；4)计算填充向量的DFT；5)将两个DFT以元素对元素的方式相乘；6)计算元素对元素乘积的逆DFT；7)返回由逆DFT计算的结果向量的前n个元素。使用FFT算法计算DFT，使用IFFT算法计算逆DFT。请注意，在步骤1)中，既不明确计算托布里兹矩阵也不计算循环矩阵。相反，仅使用它们的生成向量来执行必要的计算。

算法C.1给出了上述算法的伪代码。该算法的计算复杂度为O(n log n)。该函数的名称为TOEPLITZMULTIPLYE，其中的“E”是“嵌入(embeddding)”的缩写。

代替在该算法中包括FFT计算，可以使用FCIRCULANTMULTIPLY(请参阅下面的算法D.2)。但是，FCIRCULANTMULTIPLY比FFT计算序列执行多出O(n)的工作，因为其需要对于每个k∈{0,1,…,n-1}计算

和

附录D：使用Pustylnikov分解的托布里兹向量积

可以使用Pustylnikov分解来制定将托布里兹矩阵乘以向量的算法。此算法没有按附录C中所述将托布里兹矩阵嵌入循环矩阵中。相反，本附录描述了TOEPLITZMULTIPLY的不同版本。

下面示出了说明该算法的示例。令A为下面的4×4托布里兹矩阵，该矩阵由其第一行r和第一列c生成。

然后，可以如下将A分解为两个矩阵的和：

A＝A'+A", (D.2)

其中A'是循环矩阵，而A"是f＝-1的f-循环矩阵(即，斜循环矩阵)。对于此示例，这两个矩阵具有以下形式：

下面提供了用于计算A'和A"的值的公式。令a'＝(a'₀,a'₁,a'₂,…,a'_n-1)为A'的第一列，即

类似地，令a"＝(a"₀,a"₁,a"₂,…,a"_n-1)是A"的第一列。即

然后，可以使用以下公式来计算a'和a"的元素：

其中r是A的第一行，c是其第一列。

算法D.1给出了上述算法的伪代码。该函数的名称为TOEPLITZMULTIPLYP，其中“P”是“Pustylnikov”的缩写。第27和28行调用算法D.2，如下所述。该算法以O(n log n)的时间来运行。

算法D.2给出了函数FCIRCULANTMULTIPLY(用于将f-循环矩阵与向量相乘)的伪代码，该伪代码在算法D.1中使用。其以O(n log n)运行。

附录E：利用基于FFT的卷积的托布里兹向量积

算法E.1示出了用于将托布里兹矩阵乘以向量的又一种算法的伪代码。在这种情况下的函数名称为TOEPLITZMULTIPLYC，其中“C”是卷积的缩写。该算法的计算复杂度为O(nlog n)。与附录C中描述的方法类似，该算法可以解释为循环嵌入的一种形式。

算法E.2给出了基于FFT的卷积算法的伪代码，其是用于计算离散卷积的O(n logn)方法。其在算法E.1中使用。

本文所引用的所有参考文献，包括出版物、专利申请和专利，均通过引用并入本文，其程度如同每个参考文献被单独地且具体地指示为通过引用并入本文并在此完整阐述。

在描述本发明的上下文中(特别是在所附权利要求的上下文中)术语“个”(不定冠词a和an)和“该”(定冠词the)和类似提及应解释为涵盖单数形式和复数形式两者，除非本文另有说明或与上下文明显矛盾。除非另有说明，否则术语“包含”、“具有”、“包括”和“包含”应解释为开放式术语(即，意思是“包括但不限于”)。除非在此另外指出，否则本文中值范围的陈述仅旨在用作独立提及落入该范围内的每个单独值的简写方法，并且每个单独值都被并入说明书中，如同其在本文中被单独叙述一样。除非本文另外指出或与上下文明显矛盾，否则本文描述的所有方法可以以任何合适的顺序执行。除非另外要求，否则本文提供的任何和所有示例、或示例性语言(例如“诸如”)的使用仅旨在更好地阐明本发明，并且不对本发明的范围构成限制。说明书中的任何语言都不应解释为指示任何未要求保护的要素对于实施本发明必不可少。

本文描述了本发明的优选实施方式，包括发明人已知的用于实施本发明的最佳方式。通过阅读前述说明，那些优选实施例的变型对于本领域普通技术人员而言将变得显而易见。发明人预期本领域技术人员适当地采用这样的变型，并且发明人意图让本发明以不同于本文具体描述的方式来实践。因此，本发明包括适用法律所允许的所附权利要求书中所述主题的所有修改和等同物。而且，除非本文另外指出或与上下文明显矛盾，否则本发明涵盖上述元件在其所有可能的变化中的任何组合。

Claims

1.一种对信号进行逆线性调频Z变换(CZT)的方法，其中，所述CZT具有参数A和W通过公式A W^-k来定义复平面中的对数螺旋轮廓，k＝0,1,2,…,M-1，其中所述CZT能够被表示为X＝W A x，其中W是由参数W及其幂定义且具有维度M×N的范德蒙矩阵，A为由参数A及其幂定义的第一对角矩阵，x为第一向量，并且X是通过处理器使用所述CZT从所述第一向量x计算出的第二向量，使得X的第k个元素由

给出，其中A和W是非零复数，所述方法包括以下步骤：

将矩阵W表示为第二对角矩阵P、托布里兹矩阵

和第三对角矩阵Q的乘积，从而所述CZT被表示为

将逆线性调频Z变换(ICZT)表示为

其中A^-1、Q^-1、

和P^-1分别是A、Q、

和P的逆矩阵，其中，所述第二向量X存储在存储器中，而矩阵A^-1、Q^-1、

和P^-1中没有一个完全存储在所述存储器中；

通过使用所述处理器从右向左在所述乘积

中执行乘法以计算第三向量

来计算所述ICZT。

2.根据权利要求1所述的方法，其中所述第一向量x的每个元素基本上等于所述第三向量

的对应元素。

3.根据权利要求1所述的方法，其中所述第一向量x的至少一个元素不等于所述第三向量

的对应元素。

4.根据权利要求3所述的方法，其中在计算步骤之前，所述方法还包括对所述第二向量X执行修改其元素中的至少一个的操作的步骤。

5.根据权利要求4所述的方法，其中所述操作是滤波操作。

6.根据权利要求1所述的方法，其中所述方法具有小于或等于O(n²)的计算复杂度，其中n＝max(M,N)。

7.根据权利要求6所述的方法，其中所述计算复杂度为O(n log n)。

8.根据权利要求1所述的方法，其中所述矩阵

能够被表示为

其中

是第一下三角托布里兹矩阵，

是等于所述

的转置的第一上三角托布里兹矩阵，

是第二上三角托布里兹矩阵，

是等于所述

的转置的第二下三角托布里兹矩阵，并且u₀等于

其中n基于M或N中的至少一个。

9.根据权利要求8所述的方法，其中用于

和

四个矩阵中的每个矩阵的生成向量能够从以下给出的向量u＝(u₀,u₁,...,u_n-1)得出：

其中k＝0,1,...,n-1.

10.根据权利要求8所述的方法，其中在计算

与向量的所述乘积的所述步骤中，从右到左执行所述乘法，使得没有矩阵与另一个矩阵相乘。

11.根据权利要求10所述的方法，其中托布里兹矩阵

或

中的至少一个与向量的乘积还包括以下步骤：

将所述托布里兹矩阵嵌入到循环矩阵中，其中所述循环矩阵由其生成向量表示，并且其中所述循环矩阵的维度基于M或N中的至少一个；以及

在相乘之前用零填充所述向量以产生填充向量，其具有与所述循环矩阵的列的数目相同的长度；以及

通过将所述循环矩阵与所述填充向量相乘来计算乘积向量；以及

提取结果向量，其包括所述乘积向量中对应于嵌入了所述托布里兹矩阵的所述循环矩阵的行的元素。

12.根据权利要求11所述的方法，其中所述循环矩阵是方矩阵，并且其中所述循环矩阵的所述行的数目是2的幂。

13.根据权利要求10所述的方法，其中，托布里兹矩阵

或

中的至少一个与向量的乘积包括在与所述向量相乘之前对所述矩阵执行Pustylnikov分解。

14.根据权利要求10所述的方法，其中托布里兹矩阵

或

中的至少一个与向量的所述乘积还包括以下步骤：

通过用零填充其生成向量，将所述托布里兹矩阵嵌入具有基于M或N中的至少一个的维度的扩展托布里兹矩阵；

用零填充将所述托布里兹矩阵与之相乘的所述向量以产生填充向量，所述填充向量具有与所述扩展托布里兹矩阵中的列的数目相同的长度；

使用Pustylnikov分解将所述扩展托布里兹矩阵表示为循环矩阵和斜循环矩阵的和；

将所述循环矩阵与所述填充向量相乘，将所述斜循环矩阵与所述填充向量相乘，并将两个结果向量相加以产生填充结果向量；

通过保留所述填充结果向量中与嵌入所述托布里兹矩阵的所述扩展托布里兹矩阵的行相对应的元素，从所述填充结果向量中提取结果向量。

15.根据权利要求14所述的方法，其中，所述扩展托布里兹矩阵中的行数和列数是2的幂。

16.根据权利要求8所述的方法，其中在执行所述方法期间，矩阵

和

中没有一个完全存储在所述存储器中。

17.根据权利要求1所述的方法，其中执行所述方法所需的存储器不超过O(n)，并且其中n＝max(M,N)。

18.根据权利要求1所述的方法，其中M等于N。

19.根据权利要求1所述的方法，其中M不等于N。

20.根据权利要求1所述的方法，其中所述第一向量x从音频信号中得出。

21.根据权利要求20所述的方法，其中所述音频信号包括语音信号。

22.根据权利要求20所述的方法，其中所述音频信号包括声纳信号或超声信号中的至少一个。

23.根据权利要求1所述的方法，其中所述第一向量x从图像信号中得出。

24.根据权利要求23所述的方法，其中所述图像信号包括计算机断层摄影(CT)信号、正电子发射断层摄影(PET)信号或磁共振成像(MRI)信号中的至少一种。

25.根据权利要求1所述的方法，其中所述第一向量x从由基于雷达的传感器接收的信号中得出。

26.根据权利要求1所述的方法，其中所述第一向量x用于生成被发送到基于雷达的传感器的信号。

27.根据权利要求1所述的方法，其中所述第二向量X不是使用具有相同的参数A和W的所述CZT从特定的第一向量x中计算的。

28.一种扩展计算系统的范围以覆盖公式AW^-k的对数螺旋轮廓的方法，该计算系统对信号进行快速傅里叶变换(FFT)或逆FFT(IFFT)中的至少一个，其中A和W是非零复数，k＝0,1,…,M–1，其中所述计算系统包括处理器和存储器，所述处理器使用所述FFT计算频域向量或使用所述IFFT计算时域向量，所述存储器存储所述频域向量、所述时域向量、或者所述频域向量和所述时域向量二者，而不存储任何中间矩阵，所述方法包括使用所述处理器通过利用逆线性调频Z变换(ICZT)代替IFFT来计算所述时域向量的步骤，

其中，所述ICZT是线性调频Z变换(CZT)的逆变换，所述CZT能够被表示为X＝W A x，其中X是长度为M的频域向量，x是长度为N的时域向量，A是大小为N×N的对角矩阵，并且W是具有M×N维度的范德蒙矩阵，其中所述CZT能够由两个非零复数A和W参数化，所述A和W通过公式A W^-k定义对数螺旋轮廓，其中k＝0,1,…,M–1，并且其中在所述计算系统上利用所述ICZT的计算复杂度小于或等于O(n²)，其中n＝max(M,N)。

29.根据权利要求28所述的方法，其中所述计算复杂度为O(n log n)。

30.根据权利要求29所述的方法，其中所述计算系统在所述ICZT的执行期间使用不超过O(n)的存储器。

31.一种提高计算系统的效率的方法，所述方法针对信号对能够被表示为X＝W A x的线性调频Z变换(CZT)进行求逆，其中X是长度为M的频域向量，x是长度为N的时域向量，A是大小为N×N的对角矩阵，并且W是具有M×N维度的范德蒙矩阵，其中所述计算系统包括处理器和存储器，所述方法包括利用一种对所述CZT求逆的方法以生成时域向量的步骤，其中，所述时域向量存储在所述存储器中，而不在所述存储器中存储对所述CZT求逆的所述方法的任何中间矩阵，当在所述处理器上执行时，所述方法具有小于或等于O(n²)的计算复杂度，其中n＝max(M,N)，并且其中所述CZT能够由两个非零复数A和W参数化，所述A和W通过公式A W^-k定义对数螺旋轮廓，其中k＝0,1,…,M–1。

32.根据权利要求31所述的方法，其中所述计算系统在计算逆CZT时需要不超过O(n)的存储器。

33.根据权利要求31所述的方法，其中所述用于对所述CZT求逆的所述方法具有O(nlog n)的计算复杂度。

34.一种提高计算系统的数值精度的方法，该计算系统针对信号使用所述计算系统的处理器计算线性调频Z变换(CZT)或逆线性调频Z逆变换(ICZT)中的至少一种，其中，所述线性调频Z变换能够被表示为X＝W A x，其中X是长度为M的频域向量，x是长度为N的时域向量，A是大小为N×N的对角矩阵，并且W是具有M×N维度的范德蒙矩阵，其中，所述时域向量和所述频域向量中的一个或两者存储在所述计算系统的存储器中，而不在所述存储器中存储所述方法中涉及的任何中间矩阵，所述CZT、ICZT、或者所述CZT和ICZT两者能够由两个非零复数A和W参数化，所述A和W通过公式A W^-k定义对数螺旋轮廓，其中k＝0,1,…,M–1，其中当|W|<1时，所述方法包括对所述对数螺旋轮廓的方向反转使得该轮廓的起点和终点互换的步骤，并且其中，所述CZT的输出向量、所述ICZT的输入向量、或所述CZT的所述输出向量和所述ICZT的所述输入向量两者中元素的次序也被反转。

35.一种针对信号将长度为M的频域向量X转换为长度为N的时域向量x的方法，其中所述频域向量X和所述时域向量x存储在计算系统的存储器中，包括以下步骤：使用所述计算系统的处理器通过从右向左执行乘法来计算能够被表示为