CN102165172A - 频率成分分析装置 - Google Patents

Abstract

本发明提供一种与内燃机的旋转同步来进行针对内燃机的运转参数的检测值的频率成分分析的频率成分分析装置。以预定时间间隔对运转参数进行采样，将采样值转换为数字值。针对预定数的采样值计算与检测值内包含的多个频率成分对应的第1要素的强度和相对于第1要素相位偏移了90度的第2要素的强度。使用第1要素强度和第2要素强度，与内燃机的旋转同步地计算多个频率的频率成分强度。当内燃机的旋转速度在设定阈值以上时，将第1要素强度的累计值的一部分和第2要素强度的累计值的一部分分别置换为上次计算值来进行频率成分强度的计算。

Description

频率成分分析装置

技术领域

本发明涉及将针对内燃机的运转参数的检测值的频率成分分析与内燃机的旋转同步来进行的频率成分分析装置。

背景技术

在专利文献1中公开了这样的信号处理装置：使用离散傅立叶变换(Discrete Fourier Transform)算法(以下称为“DFT算法”)和高速傅立叶变换(Fast Fourier Transform)算法(以下称为“FFT算法”)，进行安装在内燃机上的爆震传感器的输出信号的频率成分分析。FFT算法为了提高运算速度而采用了特殊的算法，因而存在不能获得需要的频率成分(具体地说是爆震传感器输出信号的中心频率成分)的强度的情况。因此，在上述信号处理装置中，针对爆震传感器输出信号的中心频率成分(强度最大的成分)，使用DFT算法来计算强度，针对中心频率成分以外的频率成分，使用FFT算法来计算强度。

专利文献1：日本特开平11-303673号公报

爆震传感器输出信号的中心频率成分分析有必要基本上与内燃机的旋转同步进行，因而在内燃机旋转速度高的状态下，需要进一步缩短频率成分分析的运算时间。

然而，由于在专利文献1公开的装置中没有考虑到这一点，因而可确保用于使用通过频率成分分析获得的频率成分强度的爆震判别处理的时间缩短，准确的判定很有可能变得困难。

发明内容

本发明是着眼于上述这一点而作成的，本发明的目的是提供这样的频率成分分析装置：在与内燃机的旋转同步进行检测参数的频率成分分析时，可更适当地执行频率成分强度运算，特别是可提高在内燃机旋转速度高的状态下的运算速度。

为了达到上述目的，本发明提供一种与所述内燃机的旋转同步来进行针对内燃机的运转参数的检测值(VKNK)的频率成分分析的频率成分分析装置。该频率成分分析装置具有：采样单元，要素强度计算单元，以及频率成分强度计算单元。所述采样单元以预定时间(TSMP)间隔对所述运转参数进行采样，将采样值转换为数字值。所述要素强度计算单元针对预定数(ND)的采样值计算与所述检测值(VKNK)内包含的多个频率成分对应的第1要素的强度(DMFTS)和相对于所述第1要素相位偏移了90度的第2要素的强度(DMFTC)。频率成分强度计算单元使用所述第1要素强度(DMFTS)和所述第2要素强度(DMFTC)，与所述内燃机的旋转同步地计算所述多个频率的频率成分强度。而且，当所述内燃机的旋转速度(NE)在设定阈值(NETH)以上时，所述频率成分强度计算单元将所述第1要素强度的累计值的一部分和所述第2要素强度的累计值的一部分分别置换为上次计算值(TRMS1，TRMC1)来进行所述频率成分强度(STFT)的计算。

根据该结构，以预定时间间隔对内燃机的运转参数进行采样，将采样值转换为数字值。根据获得的数字值，针对预定数的采样值计算与多个频率成分对应的第1要素的强度和相对于第1要素相位偏移了90度的第2要素的强度，使用第1要素强度和第2要素强度，与内燃机的旋转同步地计算多个频率的频率成分强度。当内燃机的旋转速度在设定阈值以上时，将第1要素强度的累计值的一部分和第2要素强度的累计值的一部分分别置换为上次计算值来进行频率成分强度的计算。因此，在计算频率成分强度时，无需再次进行第1和第2要素强度的累计运算的一部分，可提高运算速度。其结果，能充分确保使用计算出的频率成分强度的后处理的时间。

期望的是，所述设定阈值(NETH)是根据获得所述预定数(ND)的采样值所需要的时间(TS)来设定的。

根据该结构，内燃机旋转速度的设定阈值是根据获得预定数的采样值所需要的时间来设定的。具体地说，当内燃机旋转速度增高且频率成分强度的运算周期比获得预定数的采样值所需要的时间短时，在频率成分强度的计算中应用的第1和第2要素的强度的上次计算值和本次计算值的一部分重复，因而能针对该重复的强度使用上次计算值。

期望的是，所述设定阈值(NE)是根据所述频率成分强度(STFT)的计算所需要的每单位时间的处理步骤数(NS)来设定的。

根据该结构，内燃机旋转速度的设定阈值是根据频率成分强度的计算所需要的每单位时间的处理步骤数来设定的。如上所述，当频率成分强度的运算周期比获得预定数的采样值所需要的时间短时，在频率成分强度的计算中应用的第1和第2要素的强度的上次计算值和本次计算值的一部分重复，因而能针对该重复的强度使用上次计算值。然而，在实际的运算装置中，存在这样的情况：通过进行将要素强度的上次计算值应用于频率成分强度的运算的处理，处理步骤数增加，反而招致运算速度的降低。因此，通过将即使进行应用上次计算值的处理、每单位时间的步骤数也不增加的内燃机旋转速度作为设定阈值，可进行与运算装置的性能(特性)对应的适当的强度运算。

期望的是，所述频率成分强度计算单元使用在以与所述内燃机的旋转同步的触发信号(CRK中断)的产生定时为中心的预定期间(TS)内获得的采样值来进行所述频率成分强度(STFT)的计算。

根据该结构，由于使用在以与内燃机的旋转同步的触发信号的产生定时为中心的预定期间内所采样的数字值来进行频率成分强度的计算，因而可进行以成为对象的内燃机的旋转角度为中心的频率成分分析。其结果，与采用使用在以产生触发信号的时刻为起点的预定期间内或者在以产生触发信号的时刻为终点的预定期间内所采样的数字值的方法相比，可获得难以受到内燃机旋转变动的影响的效果。

附图说明

图1是示出本发明的一实施方式涉及的内燃机及其控制装置的结构的图。

图2是用于说明爆震传感器输出的采样和频率成分分析的图。

图3是用于说明计算频率成分强度(STFT)的顺序和定时关系的时序图。

图4是用于说明频率成分强度(STFT)的计算期间(TS)和CRK中断的产生周期(CRME)的关系的时序图。

图5是示出发动机转速(NE)和每单位时间的运算步骤数的关系的图。

图6是用于说明计算频率成分强度(STFT)的顺序的流程图。

图7是将通过频率成分分析获得的强度数据表示为频谱时序映射图和二值化频谱时序映射图的图。

图8是示出进气门的落位噪声的二值化频谱时序映射图的图。

图9是为了说明噪声去除处理而示出二值化频谱时序映射图的图。

图10是为了说明使用主模式映射图的爆震判定而示出二值化频谱时序映射图的图。

图11是示出加权映射图的一例的图。

图12是爆震判定处理的流程图。

图13是示出在图12的处理中参照的映射图的图。

图14是在图12所示的处理中执行的二值化数据映射图计算处理的流程图。

图15是在图14所示的处理中执行的二值化处理的流程图。

图16是示出在图15的处理中参照的映射图的图。

图17是在图12所示的处理中执行的噪声去除处理的流程图。

图18是在图12所示的处理中执行的适合率计算处理的流程图。

图19是用于说明在图18的处理中参照的映射图的图。

图20是在图12所示的处理中执行的噪声学习处理的流程图。

图21是在图20所示的处理中执行的噪声映射图更新处理的流程图。

图22是示出由爆震传感器检测出的振动的频率成分强度的图。

具体实施方式

以下，参照附图说明本发明的实施方式。

图1是本发明的一实施方式涉及的内燃机(以下称为“发动机”)及其控制装置的整体结构图，例如在4气缸的发动机1的进气管2的中途配设有节气门3。节气门3与检测节气门开度TH的节气门开度传感器4连接，传感器4的检测信号被提供给电子控制单元(以下称为“ECU”)5。

针对每个气缸，在位于发动机1与节气门3之间、且处于进气管2中的未图示的进气门的稍稍上游侧位置设置有燃料喷射阀6，各喷射阀与未图示的燃料泵连接并与ECU5电连接，根据来自该ECU5的信号来控制燃料喷射阀6的开启时间。在发动机1的各气缸内设有火花塞7，火花塞7与ECU5连接。ECU5将点火信号提供给火花塞7。

在节气门3的下游侧设有检测进气压PBA的进气压传感器8以及检测进气温度TA的进气温度传感器9。在发动机1的主体上安装有检测发动机冷却水温度TW的冷却水温度传感器10和非共振型的爆震传感器11。传感器8～11的检测信号被提供给ECU5。作为爆震传感器11，例如使用可检测从5kHz到25kHz的频带的振动的传感器。

在进气管2的节气门3的上游侧设有检测吸入空气流量GA的吸入空气流量传感器13，其检测信号被提供给ECU5。

ECU5与检测发动机1的曲轴(未图示)的旋转角度的曲轴角度位置传感器12连接，与曲轴的旋转角度对应的信号被提供给ECU5。曲轴角度位置传感器12由以下传感器构成：在发动机1的特定气缸的预定曲轴角度位置输出脉冲(以下称为“CYL脉冲”)的气缸判别传感器；针对各气缸的吸入行程开始时的上止点(TDC)，在预定曲轴角度之前的曲轴角度位置(在4气缸发动机中是每180度曲轴角度)输出TDC脉冲的TDC传感器；以及以比TDC脉冲短的一定曲轴角度周期(例如6度的周期)产生1个脉冲(以下称为“CRK脉冲”)的CRK传感器，CYL脉冲、TDC脉冲以及CRK脉冲被提供给ECU5。这些脉冲用于燃料喷射正时、点火正时等的各种正时控制、发动机转速(发动机旋转速度)NE的检测。

发动机1具有气门工作特性改变装置20，气门工作特性改变装置20具有第1气门工作特性改变机构和第2气门工作特性改变机构，第1气门工作特性改变机构连续变更进气门(未图示)的气门升程量和开启角(气门开启期间)，第2气门工作特性改变机构连续变更驱动进气门的凸轮的以曲轴旋转角度为基准的工作相位。ECU5将升程量控制信号和工作相位控制信号提供给气门工作特性改变装置20，进行进气门的工作控制。第1和第2气门工作特性改变机构的结构分别在例如日本特开2008-25418号公报和日本特开2000-227013号公报中示出。

ECU5由以下电路等构成：具有对来自各种传感器的输入信号波形进行整形，将电压电平修正为预定电平，将模拟信号值转换为数字信号值等的功能的输入电路；中央运算处理单元(以下称为“CPU”)；存储由该CPU执行的各种运算程序和运算结果等的存储电路(存储器)；以及将驱动信号提供给燃料喷射阀6和火花塞7的输出电路。

在本实施方式中，进行爆震传感器11的输出信号的频率成分分析，根据该分析结果获得的频率成分强度进行爆震判定。因此首先，说明频率成分分析的概要。图2(a)示出爆震传感器11的输出信号波形，图2(b)是将该图2(a)的期间TS的波形进行放大示出的图。在本实施方式中，将采样周期TSMP设定为20微秒，以连续检测出的50个数据作为对象，进行基于离散傅立叶变换(DFT)的频率成分分析。该频率成分分析的结果如图2(c)所示。图2(c)的纵轴是频率成分强度STFT，在本实施方式中，与从5kHz到25kHz的频带中的每1kHz的频率(5，6，7，…，24，25kHz)对应的频率成分强度STFT是针对每6度曲轴角度计算出的(发动机1的曲轴每旋转6度)。

图3是用于说明上述频率成分分析的时序图，图3(a)示出与依次存储有通过每20微秒对爆震传感器11的输出信号进行采样获得的50个检测数据VKNK(数字值)的存储器的地址对应的地址编号，图3(b)示出使用检测数据VKNK所计算的正弦波成分强度DMFTS和余弦波成分强度DMFTC对应于检测数据VKNK被存储在存储器内的状态。地址编号与图3(a)相同，对应于1个检测数据VKNK的采样值，每1毫秒计算与从频率5kHz到25kHz的21个频率对应的正弦波成分强度DMFTS和余弦波成分强度DMFTC，将它们存储在存储器内。

图3(c)示出从曲轴角度位置传感器12输出的CRK脉冲，在本实施方式中，CRK脉冲的下降定时被用作运算执行定时的基准(以下称为“CRK中断”)。在本实施方式中，存储有发生CRK中断的定时，使用与在以该定时为中心的采样期间TS中获得的50个检测数据对应的正弦波成分强度DMFTS和余弦波成分强度DMFTC来计算各频率(5～25kHz)的频率成分强度STFT(j，i)(图3(d))。这里，“j”是表示频率的指数参数(以下称为“频率指数”)，j＝0、1、2、…、20对应于频率5、6、7、…、25kHz。“i”是表示发生CRK中断的曲轴角度CA(将活塞位于上止点时的角度位置设定为0度)的指数参数(以下称为“曲轴角度指数”)，i＝0、1、2、…、14对应于曲轴角度6、12、18、…、90度。

在图3所示的例子中，指数参数m(＝1～50)表示在频率成分强度STFT(j，i)的计算中应用的数据。在最开始的CRK中断定时，使用以地址编号“5”的检测数据(对应于m＝25的检测数据)为中心的50个数据来计算频率成分强度STFT(j，i)，在下一个CRK中断定时，使用以地址编号“45”的检测数据(对应于m＝25的检测数据)为中心的50个数据来计算频率成分强度STFT(j，i+1)。另外，频率成分强度STFT(j，i)是表示相对强度的无量纲量。

通过将在频率成分强度STFT(j，i)的计算中应用的数据作为以与m＝25对应的检测数据为中心的50个数据，换句话说，在以CRK中断发生定时为中心的采样期间TS中获得的50个检测数据，可进行以作为对象的曲轴角度为中心的频率成分分析。其结果，与采用使用在以发生CRK中断的时刻为起点的采样期间内或者在以发生CRK中断的时刻为终点的采样期间内所采样的检测数据的方法相比，可获得难以受到内燃机旋转变动的影响的效果。

图4是示出CRK中断发生周期CRME和在频率成分强度STFT的计算中应用的检测数据的采样期间(表示为“STFT计算”)TS的关系的图。图4(a)示出发动机转速NE是1000rpm的例子，采样期间TS和中断产生周期CRME一致。因此，当发动机转速NE超过1000rpm时，采样期间TS一部分重复。图4(b)示出发动机转速NE稍高于2000rpm的例子，采样期间TS的重复部分变长。

在本实施方式中，在这样采样期间TS重复的状态下，通过使用上次计算值作为与重复部分相关的正弦波成分强度DMFTS的累计值和余弦波成分强度DMFTC的累计值，缩短频率成分强度STFT的计算所需要的运算时间。

不过，在使用上次计算值的运算中，当可使用的重复数据数的比率小时，依赖于运算装置的特性，每单位时间的处理步骤数增加。图5所示的实线L11表示在本实施方式中使用的CPU中的在使用了上次计算值的情况下的每单位时间的处理步骤数NS和发动机转速NE的关系，在达到转速NESAT(以下称为“步骤数饱和转速NESAT”)之前，处理步骤数NS增加，在NE≥NESAT的范围内为一定值NSSAT。因此在本实施方式中，将阈值NETH设定为步骤数饱和转速NESAT，当发动机转速NE在阈值NETH以上时，进行使用上次计算值的运算。

图5所示的虚线L12表示在使用FFT(高速傅立叶变换)算法的情况下的发动机转速NE和处理步骤数NS的关系。从图5可以看出，在使用FFT算法的情况下，相对于发动机转速NE的增加，处理步骤数NS线性增加，与此相对，通过进行使用DFT算法、且在发动机转速NE高的区域(NE＞NESAT)使用上次计算值的运算，能抑制处理步骤数NS的增加。

图6是示出频率成分强度STFT计算顺序的流程图。在步骤S1中，取得检测数据VKNK，在步骤S2中，根据下述式(1)和(2)计算正弦波成分强度DMFTS(j，k)和余弦波成分强度DMFTC(j，k)。这里，指数参数k是图3(b)所示的地址编号，取从“1”到ND(在本实施方式中是“50”)的值。并且Δt相当于(采样周期×ND)，在本实施方式中是1毫秒。

DMFTS(j，k)＝VKNK(k)×sin{2π×(j+5)×1000×Δt×k/ND}  (1)

DMFTC(j，k)＝VKNK(k)×cos{2π×(j+5)×1000×Δt×k/ND}  (2)

每个采样周期(20微秒)执行步骤S1和S2的处理。

在步骤S3中，判定是否发生了CRK中断，当未发生时回到步骤S1，在发生了CRK中断的定时，判别发动机转速NE是否在阈值NETH以上。当NE＜NETH时，根据不使用上次使用数据的下述第1式(A)计算频率成分强度STFT，另一方面，当NE≥NETH时，根据使用上次使用数据的下述第2式(B)计算频率成分强度STFT(j，i)。指数参数m是被设定成使发生了CRK中断时的地址编号k为“25”的修正地址编号(参照图3(d))。指数参数mp是在上次计算时应用的修正地址编号，mx是上次边界指数，my是本次边界指数，分别由下述式(3)和(4)给出。

[算式1]

$\mathrm{STFT}(j,i)=\sqrt{{\left(\underset{m=1}{\overset{\mathrm{ND}}{\mathrm{\Σ}}}\mathrm{DMFTS}(j,m)\right)}^2+{\left(\underset{m=1}{\overset{\mathrm{ND}}{\mathrm{\Σ}}}\mathrm{DMFTC}(j,m)\right)}^2}---\left(A\right)$

mx＝ND-NV+1    (3)

my＝NV+1       (4)

这里，NV是在上次的STFT计算中使用的数据中的可在本次的STFT计算中使用的数据的数量，发动机转速NE越高，数据数NV就越增加。

第2式(B)的第1项TRMS1和第3项TRMC1相当于累计值的上次计算值，第2项TRMS2和第3项TRMC2相当于新计算的累计值。

这样通过将正弦波成分强度的累计值的一部分和余弦波成分强度的累计值的一部分分别置换为上次计算值(TRMS1，TRMC1)来进行频率成分强度STFT(j，i)的计算，从而在计算频率成分强度STFT时，无需再次进行正弦波成分强度DMFTS和余弦波成分强度DMFTC的累计运算的一部分，可提高运算速度。其结果，能充分确保使用计算出的频率成分强度的后处理的时间。

下面参照图7～图11，说明使用频率成分强度STFT的爆震判定方法的概要。

如上所述，频率成分强度STFT作为时序数据如图7(a)所示被计算为二维矩阵(以下称为“频谱时序映射图”)。在频谱时序映射图中，纵方向是频率f[kHz]，横方向是曲轴角度(燃烧行程开始的上止点后的曲轴角度)CA[度]。这里，作为纵方向和横方向的指标分别使用频率指标“j”(j＝0～20)和曲轴角度指标“i”(i＝0～14)，将频谱时序映射图的要素表示为强度参数KMAP(j，i)(强度参数KMAP在本实施方式中是与频率成分强度STFT(j，i)相同的参数，由于是频谱时序映射图上的值，因而使用KMAP这样的标签)。例如强度参数KMAP(0，0)相当于左下端的“34”，强度参数KMAP(0，14)相当于右下端的“31”，强度参数KMAP(20，0)相当于左上端的“56”，强度参数KMAP(20，14)相当于右上端的“30”。

对图7(a)的频谱时序映射图进行二值化获得图7(b)的二值化频谱时序映射图。当强度参数KMAP(j，i)在“50”以上时，作为图7(b)的映射图的要素的二值化强度参数NKMAP(j，i)取“1”，当强度参数KMAP(j，i)小于“50”时，二值化强度参数NKMAP(j，i)取“0”。

图7(b)所示的二值化频谱时序映射图是与当发生图22的虚线表示的爆震时获得的频率成分强度对应的映射图。另一方面，与图22的实线所示的进气门的落位噪声(由当发动机1的进气门到达全闭位置时产生的振动引起的噪声)对应的二值化频谱时序映射图如图8所示，由于与图7(b)所示的爆震发生时的映射图明确不同，因而不会将落位噪声误判定为发生爆震，可进行准确的爆震判定。

在本实施方式中还构成为，根据当判定为未发生爆震时的二值化频谱时序映射图，通过学习运算生成与噪声对应的二值化频谱时序映射图(以下简称为“噪声映射图”)，从根据检测数据计算出的二值化频谱时序映射图上的对应的值中减去噪声映射图上的噪声学习值，从而去除噪声的影响。

图9(a)是示出通过学习运算获得的噪声映射图的一例的图，示出在曲轴角度6度的定时、在6～25kHz的频带内的噪声成分。图9(b)示出通过从图7(b)所示的二值化频谱时序映射图中减去图9(a)所示的噪声成分来进行了噪声去除处理后的二值化频谱时序映射图。这样通过进行噪声去除处理，能消除噪声的影响，能进行更准确的判定。

并且在本实施方式中，通过将根据检测数据获得的二值化频谱时序映射图与图10(a)所示的主模式(pattern)映射图进行比较，判定是否发生爆震。主模式映射图相当于当发生爆震时获得的典型的二值化频谱时序映射图。

上述比较具体按以下方式进行。通过将二值化频谱时序映射图上的二值化强度参数NKMAP(j，i)与主模式映射图上的主参数MMAP(j，i)的积进行累计，计算强度累计值SUMK，计算主参数MMAP(j，i)自身的累计值即基准累计值SUMM，计算强度累计值SUMK对基准累计值SUMM的比率即适合率PFIT(＝SUMK/SUMM)。然后，当适合率PFIT超过判定阈值SLVL时，判定为发生了爆震。

图10(b)示出图9(b)所示的二值化强度参数NKMAP(j，i)与图10(a)所示的主模式映射图上的主参数MMAP(j，i)的积(NKMAP(j，i)×MMAP(j，i))的映射图。在该例子中，图10(b)的积映射图与图9(b)所示的映射图完全相同，适合率PFIT为“0.863”。

另外，在后述的实际的爆震判定处理中，在计算上述强度累计值SUMK和基准累计值SUMM时进行与发动机运转状态对应的加权。图11示出设定了用于进行该加权的加权参数WMAP(j，i)的加权映射图的一例。图11所示的加权映射图被设定为：对于频率6kHz的成分针对从曲轴角度12度到72度的参数值进行加权，对于频率10kHz的成分针对从曲轴角度36度到66度的参数值进行加权，以及对于频率11kHz、13kHz、15kHz和20kHz的成分针对从曲轴角度12度到30度的参数值进行加权。

基于图11所示的加权映射图的加权是为了补偿二值化频谱时序映射图的针对频率的特性依赖于发动机运转状态而变化这种情况而进行的。通过根据发动机运转状态进行加权，不管发动机运转状态是否变化，都能进行准确判定。

图12是使用上述方法进行爆震判定的处理的流程图，由ECU5的CPU与TDC脉冲的发生同步地执行该处理。

在步骤S11中，执行图14所示的二值化数据映射图计算处理，进行上述的二值化频谱时序映射图的计算。在步骤S12中，执行图17所示的噪声去除处理，使用图9(a)例示的噪声映射图来进行去除噪声成分的处理。

在步骤S13中，执行图18所示的适合率计算处理，使用去除了噪声成分后的二值化频谱时序映射图和主模式映射图来计算适合率PFIT。

在步骤S14中，根据发动机转速NE和进气压力PBA检索例如按图13所示设定的SLVL映射图，计算判定阈值SLVL。针对SLVL映射图中设定的格子点以外的发动机转速NE和进气压力PBA，通过插值运算来计算判定阈值SLVL。

在步骤S15中，判别在步骤S13计算出的适合率PFIT是否大于判定阈值SLVL，当其答案是肯定(YES)时，判定为发生爆震，将爆震标志FKNOCK设定为“1”(步骤S16)。

当在步骤S15中PFIT≤SLVL时，判定为未发生爆震，将爆震标志FKNOCK设定为“0”(步骤S17)。然后执行图20所示的噪声学习处理，进行噪声映射图(参照图9(a))的更新。

图14是在图12的步骤S11中执行的二值化数据映射图计算处理的流程图。

在步骤S21中，将曲轴角度指标i和频率指标j均初始化为“0”。在步骤S22中，判别频率指标j是否大于从频率数据数JN(在本实施方式中是21)中减去“1”后的值。由于最开始其答案是否定(NO)，因而进到步骤S23，判别曲轴角度指标i是否大于从曲轴角度数据数IN(在本实施方式中是15)中减去“1”后的值。

由于最开始步骤S23的答案也是否定(NO)，因而将强度参数KMAP(j，i)设定为存储在存储器内的频率成分强度STFT(j，i)(步骤S25)，然后执行图15所示的二值化处理(步骤S26)。在步骤S27中，使曲轴角度指标i递增“1”。

在图15的步骤S31中，根据发动机转速NE和进气压力PBA检索图16所示的BLVL映射图，计算二值化阈值BLVL。BLVL映射图被设定成：发动机转速NE越增大，则二值化阈值BLVL就越增大，而且进气压力PBA越增大，则二值化阈值BLVL就越增大。由线L1表示的设定值适用于从第1预定进气压力PBA1(例如53kPa(400mmHg))至第2预定进气压力PBA2(例如80kPa(600mmHg))的范围。线L2和L3分别对应于第3预定进气压力PBA3(例如93kPa(700mmHg))和第4预定进气压力PBA4(例如107kPa(800mmHg))。

在步骤S32中，判别强度参数KMAP(j，i)是否大于二值化阈值BLVL，当其答案是肯定(YES)时，将二值化强度参数NKMAP(j，i)设定为“1”(步骤S33)。另一方面，当在步骤S32中KMAP(j，i)≤BLVL时，将二值化强度参数NKMAP(j，i)设定为“0”(步骤S34)。

回到图14，在步骤S23的答案是否定(NO)的期间，重复执行步骤S24～S27，当曲轴角度指标i超过(IN-1)时，进到步骤S28，使曲轴角度指标i恢复为“0”，并使频率指标j递增“1”，回到步骤S22。

在步骤S22的答案是否定(NO)的期间，重复执行步骤S23～S28，当频率指标j超过(JN-1)时，结束本处理。

通过图14的处理，将表示图7(a)所示的频谱时序映射图上的值的强度参数KMAP(j，i)设定为存储在存储器内的作为频率成分分析结果获得的频率成分强度STFT(j，i)，并进行强度参数KMAP(j，i)的二值化，计算二值化强度参数NKMAP(j，i)。即，生成图7(b)所示的二值化频谱时序映射图。

图17是在图12的步骤S12中执行的噪声去除处理的流程图。

在步骤S41中，将曲轴角度指标i和频率指标j都初始化为“0”。在步骤S42中，判别频率指标j是否大于从频率数据数JN中减去“1”后的值。由于最开始其答案是否定(NO)，因而进到步骤S43，判别曲轴角度指标i是否大于从曲轴角度数据数IN中减去“1”后的值。

由于最开始步骤S43的答案也是否定(NO)，因而根据下述式(6)校正二值化强度参数NKMAP(j，i)，计算校正二值化强度参数JKMAP(j，i)。式(6)的NNMAP(j，i)是通过学习处理而更新的噪声映射图上的二值化噪声参数。

JKMAP(j，i)＝NKMAP(j，i)-NNMAP(j，i)    (6)

在步骤S45中，判别校正二值化强度参数JKMAP(j，i)是否是负值，当其答案是否定(NO)时，立即进到步骤S47。当JKMAP(j，i)＜0时，将JKMAP(j，i)设定为“0”(步骤S46)，进到步骤S47。

在步骤S47中，使曲轴角度指标i递增“1”，回到步骤S43。在步骤S43的答案是否定(NO)的期间，重复执行步骤S44～S47，当曲轴角度指标i超过(IN-1)时，进到步骤S48，使曲轴角度指标i恢复为“0”，并使频率指标j递增“1”，回到步骤S42。在步骤S42的答案是否定(NO)的期间，重复执行步骤S43～S48，当频率指标j超过(JN-1)时，结束本处理。

通过图17的处理，获得去除了噪声后的校正二值化强度参数JKMAP(j，i)(参照图9(b))。

图18是在图12的步骤S13中执行的适合率计算处理的流程图。

在步骤S51中，根据发动机转速NE和进气压力PBA选择主模式映射图，在步骤S52中，根据发动机转速NE和进气压力PBA选择加权映射图。加权映射图是为了补偿二值化频谱时序映射图的针对频率的特性依赖于发动机运转状态而变化这种情况而设置的。当发动机转速NE或进气压力PBA(发动机负荷)变化时，燃烧室内的温度变化，二值化频谱时序映射图变化。因此，通过根据发动机转速NE和进气压力PBA选择主模式映射图和加权映射图，不管发动机运转状态是否变化，都能进行准确的判定。

在本实施方式中，如图19所示对应于根据发动机转速NE和进气压力PBA而定义的9个发动机运转区域，预先设定了9个主模式映射图和9个加权映射图，在步骤S51中选择9个主模式映射图中的1个，在步骤S52中，选择9个加权映射图中的1个。在图19中，低转速区域被设定为例如发动机转速NE在2000rpm以下的区域，中转速区域被设定为从2000rpm到4000rpm的区域，高转速区域被设定为超过4000rpm的区域。并且，低负荷区域被设定为例如进气压力PBA在67kPa(500mmHg)以下的区域，中负荷区域被设定为从67kPa到93kPa(700mmHg)的区域，高负荷区域被设定为超过93kPa的区域。

在步骤S53中，将曲轴角度指标i和频率指标j都初始化为“0”，并将强度累计值SUMK和基准累计值SUMM初始化为“0”。强度累计值SUMK和基准累计值SUMM在后述的步骤S57中被更新，在步骤S60中被应用于适合率PFIT的计算。

在步骤S54中，判别频率指标j是否大于从频率数据数JN中减去“1”后的值。由于最开始其答案是否定(NO)，因而进到步骤S55，判别曲轴角度指标i是否大于从曲轴角度数据数IN中减去“1”后的值。

由于最开始步骤S55的答案也是否定(NO)，因而进到步骤S56，根据下述式(7)和(8)计算加权主参数MMW和加权积参数KMW。下述式的WMAP(j，i)是在加权映射图中设定的加权参数。加权积参数KMW是使用加权参数WMAP(j，i)将主参数MMAP(j，i)和校正二值化强度参数JKMAP(j，i)的积进行了加权后的参数。

MMW＝MMAP(j，i)×WMAP(j，i)                 (7)

KMW＝MMAP(j，i)×JKMAP(j，i)×WMAP(j，i)    (8)

在步骤S57中，根据下述式(9)和(10)，将加权主参数MMW和加权积参数KMW进行累计，计算基准累计值SUMM和强度累计值SUMK。

SUMM＝SUMM+MMW    (9)

SUMK＝SUMK+KMW    (10)

在步骤S58中，使曲轴角度指标i递增“1”，回到步骤S55。在步骤S55的答案是否定(NO)的期间，重复执行步骤S56～S58，当曲轴角度指标i超过(IN-1)时，进到步骤S59，使曲轴角度指标i恢复为“0”，并使频率指标j递增“1”，回到步骤S54。在步骤S54的答案是否定(NO)的期间，重复执行步骤S55～S59，当频率指标j超过(JN-1)时，进到步骤S60，根据下述式(11)计算适合率PFIT。

PFIT＝SUMK/SUMM    (11)

图20是在图12的步骤S18中执行的噪声学习处理的流程图。

在步骤S71中，将曲轴角度指标i、频率指标j、加法学习参数LK以及减法学习参数LM均初始化为“0”。在步骤S72中，判别频率指标j是否大于从频率数据数JN中减去“1”后的值。由于最开始其答案是否定(NO)，因而进到步骤S73，判别曲轴角度指标i是否大于从曲轴角度数据数IN中减去“1”后的值。

由于最开始步骤S73的答案也是否定(NO)，因而进到步骤S74，判别二值化强度参数NKMAP(j，i)是否等于二值化噪声参数NNMAP(j，i)。当其答案是肯定(YES)时，立即进到步骤S80。

当步骤S74的答案是否定(NO)、二值化强度参数NKMAP(j，i)与二值化噪声参数NNMAP(j，i)不同时，判别二值化强度参数NKMAP(j，i)是否大于二值化噪声参数NNMAP(j，i)(步骤S75)。当其答案是肯定(YES)时，将加法学习参数LK设定为“1”，将减法学习参数LM设定为“0”(步骤S76)。另一方面，当NKMAP(j，i)＜NNMAP(j，i)时，将加法学习参数LK设定为“0”，将减法学习参数LM设定为“1”(步骤S77)。

在步骤S78中，根据下述式(12)和(13)，修正加法学习参数LK和减法学习参数LM。式(12)、(13)的DSNOISE是被设定为例如0.1的噪声学习系数。

LK＝DSNOISE×LK    (12)

LM＝DSNOISE×LM    (13)

在步骤S79中，将加法学习参数LK和减法学习参数LM应用于下述式(14)，来更新噪声参数NMAP(j，i)。在下面说明的图21的处理中，通过将噪声参数NMAP(j，i)进行二值化，计算出二值化噪声参数NNMAP(j，i)。

NMAP(j，i)＝NMAP(j，i)+LK-LM    (14)

在步骤S80中，使曲轴角度指标i递增“1”，回到步骤S73。在步骤S73的答案是否定(NO)的期间，重复执行步骤S74～S80，当曲轴角度指标i超过(IN-1)时，进到步骤S81，使曲轴角度指标i恢复为“0”，并使频率指标j递增“1”，回到步骤S72。在步骤S72的答案是否定(NO)的期间，重复执行步骤S73～S81，当频率指标j超过(JN-1)时，进到步骤S82，执行图21所示的噪声映射图更新处理。

在图21的步骤S91中，将曲轴角度指标i和频率指标j都初始化为“0”。在步骤S92中，判别频率指标j是否大于从频率数据数JN中减去“1”后的值。由于最开始其答案是否定(NO)，因而进到步骤S93，判别曲轴角度指标i是否大于从曲轴角度数据数IN中减去“1”后的值。

由于最开始步骤S93的答案也是否定(NO)，因而进到步骤S94，判别噪声参数NMAP(j，i)是否大于噪声二值化阈值NLVL(例如0.8)。当其答案是肯定(YES)时，将二值化噪声参数NNMAP(j，i)设定为“1”(步骤S95)。另一方面，当NMAP(j，i)≤NLVL时，将二值化噪声参数NNMAP(j，i)设定为“0”(步骤S96)。

在步骤S97中，使曲轴角度指标i递增“1”，回到步骤S93。在步骤S93的答案是否定(NO)的期间，重复执行步骤S94～S97，当曲轴角度指标i超过(IN-1)时，进到步骤S98，使曲轴角度指标i恢复为“0”，并使频率指标j递增“1”，回到步骤S92。在步骤S92的答案是否定(NO)的期间，重复执行步骤S93～S98，当频率指标j超过(JN-1)时，结束本处理。

当在图20的步骤S75中NKMAP(j，i)＞NNMAP(j，i)时，在步骤S78中为LK＝0.1，LM＝0，根据式(14)，使噪声参数NMAP(j，i)增加“0.1”。另一方面，当在步骤S75中NKMAP(j，i)＜NNMAP(j，i)时，在步骤S78中LK＝0，LM＝0.1，根据式(14)，使噪声参数NMAP(j，i)减小“0.1”。然后，当在图21的处理中噪声参数NMAP(j，i)大于噪声二值化阈值NLVL时，将二值化噪声参数NNMAP(j，i)设定为“1”，另一方面，当噪声参数NMAP(j，i)在噪声二值化阈值NLVL以下时，将二值化噪声参数NNMAP(j，i)设定为“0”。

通过图20和图21的处理，根据当判定为未发生爆震时的二值化强度参数NKMAP(j，i)来更新噪声映射图，例如把如进气门的落位噪声那样稳定发生的噪声反映到噪声映射图上。其结果，能去除噪声的影响来进行高精度的判定。

另外，也可以不进行噪声去除处理(图12，步骤S12)，直接使用二值化强度参数NKMAP来进行适合率PFIT的计算。这是因为，在不进行噪声去除处理的情况下，受到噪声影响的可能性高，然而可以如参照19所说明那样，区别于噪声来进行爆震的判定。

如以上详述，在本实施方式中，以6度曲轴角度的间隔进行爆震传感器11的输出信号的频率成分分析，生成作为其结果而获得的5kHz～25kHz的频率成分强度的时序数据即频谱时序映射图。即，频谱时序映射图的要素作为二维排列数据即强度参数KMAP(j，i)被存储在存储器内。然后，通过对强度参数KMAP(j，i)进行二值化，计算出二值化强度参数NKMAP(j，i)，根据该二值化强度参数NKMAP(j，i)判定是否发生了爆震。由于在二值化强度参数NKMAP(j，i)上反映了伴随发动机旋转的频率成分分布的变化，因而通过将二值化强度参数NKMAP(j，i)与和爆震发生时特有的变化模式对应的主模式映射图上的主参数MMAP(j，i)进行比较，可准确判定爆震的发生。并且通过对作为频率成分分析的结果而获得的强度参数KMAP(j，i)进行二值化，数据量减少，并且时序数据的变化模式简化，因而可减少存储器容量，并可提高运算速度。

并且当二值化强度参数NKMAP(j，i)(JKMAP(j，i))表现出接近主参数MMAP(j，i)的变化模式时，发生了爆震的可能性高，因而通过计算表示二值化强度参数NKMAP(j，i)(JKMAP(j，i))和主参数MMAP(j，i)的类似性(相关性)的参数，可进行准确判定。

在本实施方式中，作为表示该类似性(相关性)的参数使用适合率PFIT，当适合率PFIT超过判定阈值SLVL时，判定为发生了爆震。通过使用适合率PFIT，可利用比较简单的运算可靠地评价二值化强度参数NKMAP(j，i)(JKMAP(j，i))和主参数MMAP(j，i)的类似性(相关性)，可进行准确判定。

并且，根据二值化强度参数NKMAP(j，i)计算噪声成分的时序数据NNMAP(j，i)，利用噪声成分的时序数据即噪声参数NNMAP(j，i)校正二值化强度参数NKMAP(j，i)，根据校正二值化强度参数JKMAP(j，i)进行爆震判定。因此，能去除如上述的落位噪声那样表现稳定的噪声成分来进行更准确的判定。

并且，对二值化强度参数NKMAP(j，i)(JNKMAP(j，i))和主参数MMAP(j，i)乘以根据频率而设定的加权参数WMAP(j，i)，计算适合率PFIT。由于预先判明了在发生爆震时增大的频率成分，因而通过对与该频率附近的频率对应的参数进行大的加权，可提高判定精度。

在本实施方式中，ECU5构成采样单元、要素强度计算单元以及频率成分强度计算单元。

另外，本发明不限于上述的实施方式，能进行各种变形。例如，爆震传感器输出的采样周期和进行频率成分分析的曲轴角度间隔不限于上述的采样周期和曲轴角度间隔(20微秒，6度)，能在可达到本发明目的的范围内进行变更。并且，二值化频谱时序映射图(在上述实施方式中由21行×15列的矩阵构成)也同样能变更。

并且在上述的实施方式中，示出了进行爆震传感器的检测值的频率分析的装置，然而本发明不限于此，能应用于与发动机旋转(CRK中断)同步进行发动机运转参数检测值的频率成分分析的情况。例如本发明也能应用于从CRK传感器输出的CRK脉冲(中断)的产生周期CRME或者作为其倒数获得的发动机1的瞬时旋转速度的频率成分分析等。

并且在上述的实施方式中，将发动机转速NE的阈值NETH设定为每单位时间的运算步骤数饱和的步骤数饱和转速NESAT，然而也可以将在频率成分强度的计算中使用的数据的采样期间TS比CRK中断发生周期CRME短的转速NETS、或者比转速NETS稍高(例如高100rpm左右)的转速设定为阈值NETH。在上述的实施方式中，转速NETS是1000rpm。

并且在上述的实施方式中，对二值化强度参数NKMAP(j，i)和主参数MMAP(j，i)乘以加权参数WMAP(j，i)，计算出适合率PFIT，然而也可以不乘以加权参数WMAP(j，i)，即不进行加权来计算。

并且在上述的实施方式中，根据发动机转速NE和进气压力PBA计算或选择判定阈值SLVL、二值化阈值BLVL、主模式映射图以及加权映射图，然而也可以固定为预先设定的值或1个映射图。

标号说明

1：内燃机；5：电子控制单元(采样单元、要素强度计算单元、频率成分强度计算单元)；11：爆震传感器。

Claims (8)

1.一种频率成分分析装置，该频率成分分析装置与内燃机的旋转同步地进行针对所述内燃机的运转参数的检测值的频率成分分析，所述频率成分分析装置具有：

采样单元，其以预定时间间隔对所述运转参数进行采样，将采样值转换为数字值；

要素强度计算单元，其针对预定数目的采样值计算与所述检测值中包含的多个频率成分对应的第1要素的强度以及相对于所述第1要素相位偏移了90度的第2要素的强度；以及

频率成分强度计算单元，其使用所述第1要素强度和所述第2要素强度，与所述内燃机的旋转同步地计算所述多个频率的频率成分强度，

当所述内燃机的旋转速度大于等于设定阈值时，所述频率成分强度计算单元将所述第1要素强度的累计值的一部分和所述第2要素强度的累计值的一部分分别置换为上次的计算值来进行所述频率成分强度的计算。

2.根据权利要求1所述的频率成分分析装置，其中，所述设定阈值是根据获得所述预定数目的采样值所需要的时间而设定的。

3.根据权利要求1所述的频率成分分析装置，其中，所述设定阈值是根据所述频率成分强度的计算所需要的每单位时间的处理步骤数而设定的。

4.根据权利要求1至3中任一项所述的频率成分分析装置，其中，所述频率成分强度计算单元使用在以与所述内燃机的旋转同步的触发信号的产生定时为中心的预定期间内获得的采样值来进行所述频率成分强度的计算。

5.一种频率成分分析方法，该频率成分分析方法与内燃机的旋转同步地进行针对所述内燃机的运转参数的检测值的频率成分分析，所述频率成分分析方法具有以下步骤：

a)以预定时间间隔对所述运转参数进行采样，将采样值转换为数字值；

b)针对预定数目的采样值计算与所述检测值中包含的多个频率成分对应的第1要素的强度以及相对于所述第1要素相位偏移了90度的第2要素的强度；以及

c)使用所述第1要素强度和所述第2要素强度，与所述内燃机的旋转同步地计算所述多个频率的频率成分强度，

当所述内燃机的旋转速度大于等于设定阈值时，将所述第1要素强度的累计值的一部分和所述第2要素强度的累计值的一部分分别置换为上次的计算值来进行所述频率成分强度的计算。

6.根据权利要求5所述的频率成分分析方法，其中，所述设定阈值是根据获得所述预定数目的采样值所需要的时间而设定的。

7.根据权利要求5所述的频率成分分析方法，其中，所述设定阈值是根据所述频率成分强度的计算所需要的每单位时间的处理步骤数而设定的。

8.根据权利要求5至7中任一项所述的频率成分分析方法，其中，所述频率成分强度计算单元使用在以与所述内燃机的旋转同步的触发信号的产生定时为中心的预定期间内获得的采样值来进行所述频率成分强度的计算。

Images