一种基于振动信号的高铁钢轨伤损探测方法及装置的制作方法

文档序号:6005113阅读:281来源:国知局
专利名称:一种基于振动信号的高铁钢轨伤损探测方法及装置的制作方法
技术领域
本发明涉及一种基于振动信号的高铁钢轨伤损探测方法及装置,属于铁路安全监测与防护技术领域。
背景技术
铁路是交通运输的大动脉,对国民经济的发展起着十分重要的作用。随着列车运行速度的不断提高及密度和载重的加大,特别是高铁的飞速发展,对现有的铁路钢轨探伤方式提出了新的挑战。现有的铁路钢轨探伤是以大型钢轨探伤车为主、小型钢轨探伤仪为辅的探伤模式。而目前国内钢轨探伤车的速度很低(一般为50-80公里/小时),但高铁的时速却在 250-350公里/小时之间,当这些设备对钢轨进行检测时将占用高铁线路,严重地影响了列车运行的效率。同时,小型钢轨探伤仪多为人推探伤,只能进行局部范围内有限的低速检测。显然,现有的探伤方法、探伤设备、探伤速度远远不能满足高铁发展的需要。

发明内容
本发明目的是为了解决现有的测量设备检测高铁钢轨伤损情况时会严重影响列车运行效率的问题,提供了一种基于振动信号的高铁钢轨伤损探测方法及装置。本发明所述一种基于振动信号的高铁钢轨伤损探测方法包括以下步骤步骤一、通过振动加速度传感器采集高铁钢轨的振动信号,并对所述振动信号进行经验模态分解,获取η个IMF分量和一个残差;IMF为内固模态函数;步骤二、计算步骤一获取的η个IMF分量与分解之前采集到的振动信号的相关性系数Ui(i = 1,2, ...,n),选取相关性系数1^大于阈值λ的IMF分量作为伤损特征信息 IMF ;步骤三、分别计算步骤二获取的伤损特征信息IMF的功率谱密度,步骤四、根据步骤三获取的伤损特征信息IMF的功率谱密度确定高铁钢轨伤损部位,完成高铁钢轨伤损的探测。基于上述一种基于振动信号的高铁钢轨伤损探测方法的装置它包括振动加速度传感器、信号转换电路、经验模态分解模块、伤损特征信息提取模块、伤损特征信息的功率密度获取模块和高铁钢轨的伤损部位判定模块,设置在高铁钢轨上的振动加速度传感器检测高铁钢轨的伤损信息,并输出给信号转换电路,信号转换电路将接收的模拟信号转换为数字形式的输入伤损振动信号为x(t)并输出给经验模态分解模块,经验模态分解模块将输入伤损振动信号为X (t)进行处理,将获取的η个IMF分量输出给伤损特征信息提取模块,伤损特征信息提取模块将满足伤损条件的IMF分量输出给伤损特征信息的功率密度获取模块,伤损特征信息的功率密度获取模块计算出每个包含伤损特征的IMF分量的功率密度,并输出给高铁钢轨的伤损部位判定模块,高铁钢轨的伤损部位判定模块根据各IMF分量的功率密度判断被测高铁钢轨的伤损部位。
本发明的优点1)针对现有钢轨探伤方式效率较低的情况,在基于铁路沿线振动传感器的基础上,有效的利用了实时振动信号并挖掘出特征信息,提出了一种详细的伤损判断方法及装置,可实现对整个路段的钢轨进行实时探测的目的。2)利用EMD对信号时间尺度的自适应能力,较好地处理了非平稳、非线性的钢轨振动信号;并且通过相应的PSD计算,从中振动信号中准确地获得钢轨伤损特征频率。其中,利用相关性分析提取出包含主要伤损特征信息的IMF的方法,防止了 EMD分解中低频范围内引入无关IMF的干扰,以进一步精确分析伤损特征。


图1为本发明的流程图;图2为钢轨振动信号检测装置结构示意图;图3为经验模态分解(EMD)流程图;图4为钢轨振动传感器测点布置图;图5为三种不同的伤损在钢轨上的位置;图6为无损钢轨振动信号的各阶IMF和残差;图7为无损钢轨振动信号的第一阶IMF及其相应PSD ;图8为轨头伤损时钢轨振动信号的第一阶IMF及其相应PSD ;图9为轨底伤损时钢轨振动信号的第一阶IMF及其相应PSD ;图10为轨腰伤损时钢轨振动信号的第一阶IMF及其相应PSD ;图11为无损以及三种不同位置损伤振动信号的PSD图。
具体实施例方式具体实施方式
一下面结合图1、图2和图5说明本实施方式,本发明针对现有的探伤技术和探伤机制无法满足高速铁路的探伤需求,采用一种新的基于振动信号特征提取的高速列车伤损探测方法。通过在铁路沿线建立相应的轨道伤损探测传感器网络,来测量钢轨振动信号,并对信号进行预处理、发送给探伤车或特定的信息中心;探伤车或信息中心接收到信号后,利用经验模态分解(EMD)对接收到的振动信号进行分解,从而得到相应的内固模态函数(IMF),所得到的各阶IMF分量包含了钢轨伤损的特征信息(该特征信息包含无损和有损,而有损又包含伤损位置、程度等信息)。在EMD分解过程中,往往部分IMF分量明显反映伤损特征,而其它IMF分量基本不反映伤损特征或反映不明显,为此需要选取与分解之前信号相关性较大的IMF对其特征进行分析并提取相应的特征。通过计算包含伤损特征信息IMF的功率谱密度(PSD),研究PSD和钢轨损伤的关系,从中得出伤损的频率范围和特征。最后,建立钢轨损伤的存在及其位置的判断准则。本发明通过以上方法提出一种针对高铁伤损钢轨的实时探测的具体方法和装置, 也就是利用铁路沿线的振动传感器获得其振动信号,通过对该信号进行EMD分解和相关 IMF的PSD计算,得到相应的钢轨伤损特征,然后对伤损特征进行分析并建立钢轨伤损判断准则。本发明改变了以往传统的铁路钢轨探伤模式,不仅提高了探伤的效率,而且能做到对整个路段钢轨的实时探测。
本实施方式所述一种基于振动信号的高铁钢轨伤损探测方法包括以下步骤步骤一、通过振动加速度传感器采集高铁钢轨的振动信号,并对所述振动信号进行经验模态分解,获取η个IMF分量和一个残差;步骤二、计算步骤一获取的η个IMF分量与分解之前采集到的振动信号的相关性系数Ui(i = 1,2,. . .,η),选取相关性系数Ui大于阈值λ的IMF分量作为包含伤损特征信息的IMF ;步骤三、分别计算步骤二获取的伤损特征信息IMF的功率谱密度,功率谱密度简写为PSD ;步骤四、根据步骤三获取的伤损特征信息IMF的功率谱密度确定高铁钢轨伤损部位,完成高铁钢轨伤损的探测。高铁钢轨的伤损可以发生在不同的位置,比如轨头、轨腰和轨底,具体参见图5所示。步骤一中通过振动加速度传感器采集它的振动信号,分析这些振动信号中的异常情况。步骤二中大于阈值λ的相关性系数Ui可以是一个或多个,如为一个,则将该IMF 分量作为包含伤损特征信息的IMF分量,如为多个,则将这些IMF分量都作为包含伤损特征信息的IMF分量。步骤三中选择其中包含伤损特征信息的IMF并计算其PSD,然后,研究其PSD和钢轨损伤的关系,从中得出不同位置伤损的频率范围和特征,建立损伤的存在及其位置的判断准则。参见图2所示,与被测试钢轨相连的振动传感器装置,与该传感器相连的信号转换电路,以及连接于该信号转换电路的微型计算机,然后通过计算机进行伤损的匹配,确立伤损准则。其中,传感器测量的振动信号经过转换后,可以通过有线直接连接到附近检测中心的计算机,也可以通过无线网络发送到检测中心。
具体实施方式
二 下面结合图4说明本实施方式,本实施方式对实施方式一作进一步说明步骤一中通过振动加速度传感器采集高铁钢轨的振动信号,将振动加速度传感器设置在高铁钢轨的轨底上表面,且距离轨底外边缘d距离处,所述d的取值范围是 18mm 20mmo具体实施方式
三下面结合图3说明本实施方式,本实施方式对实施方式一作进一步说明步骤一中获取η个IMF分量和一个残差的过程为设定输入伤损振动信号为χ⑴,t = 1,2,· · ·,N,步骤11、IMF分解过程初始化n = 1,且满足关系式IV1U) = x(t)成立,其中 ^(t)为第(n-1)次分解后趋势函数;步骤12、筛选过程初始化,k = 1,且满足关系式Iin0ri) (t) = (t)成立,其中 K(^1) (t)为第η次经验模态分解中经过第(k-Ι)次筛选后的剩余函数;步骤13、根据筛选程序获取经过第k次筛选后的剩余函数hnk(t);步骤14、采用标准偏差准则判断步骤13获得的剩余函数hnk(t)是否满足本
征模态函数IMF的条件,即(/^“_υ0;)- "-(Ο)2 /是否小于阈值&D,
0. 2 彡 Hsd 彡 0. 3 ;判断结果为是,执行步骤15,判断结果为否,则k = k+Ι,然后执行步骤13,步骤15、提取一个本征模态函数分量Cn (t)=、⑴;
步骤16、获取输入伤损振动信号x(t)经过第η次经验模态分解的趋势函数rn(t) =rn_! (t) -cn (t);步骤17、判断趋势函数rn(t)是否为单调函数;判断结果为否,则n = η+l,然后执行步骤12,判断结果为是,完成提取过程,获得 η 个 IMF 分量{Cl(t), c2(t)-cn(t)};和 1 个残差 RES :rn(t)。步骤13中获取剩余函数hnk (t)的过程为步骤131、利用三次样条函数获取输入伤损振动信号x(t)经过第η次经验模态分解中经过第k-Ι次筛选后的剩余函数Kari) (t)的上、下包络,步骤132、计算所述剩余函数Ilna^a)上、下包络曲线在各个t的均值,Ο),步骤133、获取输入伤损振动信号χ (t)经过第η次经验模态分解中经过第k次筛选后的剩余函数/^⑴=Hnik^t)-mn{k-i){t)。
具体实施方式
四本实施方式与实施方式三的不同之处在于,步骤14的中Hsd = 0. 25,其它与实施方式三相同。
具体实施方式
五本实施方式对实施方式一作进一步说明步骤二中的相关性系数Ui按如下公式获取
IMFnUi =-
,其中,IMFn为振动信号进行经验模态分解获取的η个IMF分量,x(t)为分解出该 η个IMF的分解之前的振动信号。步骤二中的阈值λ = 0. 4 0. 6。
具体实施方式
六本实施方式与实施方式一的不同之处在于,步骤二中阈值λ = 0.5,其它与实施方式一相同。
具体实施方式
七本实施方式对实施方式一作进一步说明步骤三中的包含伤损特征信息IMF的功率谱密度按如下公式获取S(f) =「R{z)e—2mfT τ,
‘-OO其中,Rh)为包含伤损特征信息IMF的自相关函数。步骤四中确定高铁钢轨伤损部位的判断准则为若包含伤损特征信息IMF的功率谱密度曲线中存在频率为3880Hz 士 IOOHz的频率点,且该频率点幅值处于无损钢轨特征频率3360Hz的幅值上下浮动5%范围内,则判定为高铁钢轨伤损部位为轨头;若包含伤损特征信息IMF的功率谱密度曲线中存在频率为3880Hz 士 IOOHz的频率点,且该频率点幅值处于无损钢轨特征频率3360Hz的幅值的300% 350%范围内,则判定为高铁钢轨伤损部位为轨底;若包含伤损特征信息IMF的功率谱密度曲线中存在频率为3600Hz 士 IOOHz的频率点且该频率点幅值处于无损钢轨特征频率3360Hz的幅值的1200% 1300%范围内,则判定为高铁钢轨伤损部位为轨腰。
具体实施方式
八下面结合图2说明本实施方式,实现实施方式一所述的一种基于振动信号的高铁钢轨伤损探测方法的装置,它包括振动加速度传感器1、信号转换电路 2、经验模态分解模块3、伤损特征信息提取模块4、伤损特征信息的功率密度获取模块5和高铁钢轨的伤损部位判定模块6,设置在高铁钢轨上的振动加速度传感器1检测高铁钢轨的伤损信息,并输出给信号转换电路2,信号转换电路2将接收的模拟信号转换为数字形式的输入伤损振动信号为x(t)并输出给经验模态分解模块3,经验模态分解模块3将输入伤损振动信号为χ (t)进行处理,将获取的η个IMF分量输出给伤损特征信息提取模块4,伤损特征信息提取模块4将满足伤损条件的IMF分量输出给伤损特征信息的功率密度获取模块 5,伤损特征信息的功率密度获取模块5计算出每个包含伤损特征的IMF分量的功率密度, 并输出给高铁钢轨的伤损部位判定模块6,高铁钢轨的伤损部位判定模块6根据各IMF分量的功率密度判断被测高铁钢轨的伤损部位。
具体实施方式
九下面结合图2、图6至图11说明本实施方式,本实施方式给出一个具体的实施例本发明提出的一种基于振动信号的高铁钢轨伤损探测装置的实施例,如图2所示,包括以下几个部分振动加速度传感器,信号转换电路,以及微型计算机及其相应的接口电路。各部分具体实现结构及工作原理说明如下首先,钢轨结构振动是机车车辆对钢轨动态作用的结果,采用振动加速度传感器对其进行测量。传感器的型号为ULT20系列的ULT2003/V压电加速度传感器;钢轨的型号为60Kg/m钢轨。本实施例中传感器的布置如图4所示,传感器位于轨底上表面距离轨底外边缘20mm处。经过多次的实验,该点能很好的测得不同钢轨伤损部位的振动信号(伤损分别位于轨头、轨腰和轨底如图5)并区分出它们之间的伤损特征信息。针对无损伤以及伤损发生在轨头、轨腰和轨底三种不同位置情况下钢轨的振动信号进行了测量,得到了 4组测量数据用于分析。其次,测量后的振动信号,通过信号转换电路进行采样,转换芯片的型号为 AD9071,是一个10位,100MSPS模数转换芯片。采用100M高速采样,可获得更为详细的信
肩、ο最后,将转换后的信号通过有线或者无线方式发送到微型计算机进行伤损信号特征提取和伤损匹配,建立伤损判断准则。微型计算机型号为ThinkpadR400笔记本。下面结合实施例和

本发明的检测方法
具体实施例方式执行步骤一将测得的钢轨不同伤损部位(伤损位于轨头、轨腰和轨底三种位置) 的振动信号以及无损伤振动信号χ⑴分别进行EMD分解,得到各个信号的各阶IMF和残差,记为Cl,c2,...CnW&r。以无损伤钢轨的振动信号为例,分解后的各阶IMF和残差如图6所示。同理,可以得到其他三种信号分解后的各阶IMF和残差。执行步骤二 分别计算三种伤损情况以及无损情况振动信号分解后的η个IMF同其各自的未分解之前信号的相关性系数^。同样,以无损伤钢轨的振动信号为例,其η个 IMF同其未分解之前信号的相关性系数见表1。同理,可以得到其他三种信号各自的η个 IMF同其各自未分解之前信号的相关性系数。表1无损钢轨的振动信号各阶IMF同其未分解之前信号的相关性系数IMFlIMF2IMF3IMF4IMF5IMF6IMF7相关性系数0. 87980. 38660. 05250. 02020. 01260. 01040. 0075选取λ = 0. 5,这样就可以选择相关性较强包含伤损特征信息的IMF进行分析。 从表1中可以看出满足Ui > λ = 0. 5的为0. 8798,其他IMF的相关性较弱,因此选择IMFl 进行分析。同理,对其他三种伤损信号进行IMF相关性分析且选择相关性较大的IMF。本实例中通过选择后都为各自信号分解后的IMF1。执行步骤三利用功率谱密度对步骤二中所得的四种不同伤损情况的IMFl进行计算,得到不同损伤情况下(轨头、轨腰和轨底的伤损以及无损伤)振动信号的PSD,从不同伤损情况下PSD图的对比中得出钢轨伤损的频率范围和特征,最后得出判断准则。以无伤损的钢轨振动信号为例,其分解后的IMFl及其PSD如图7。从图7中可以看出,无伤损钢轨振动信号的频率范围是0-4000HZ。振幅较高的频率是1840Hz和3360Hz, 这是无损钢轨的固有振动频率。同理,轨头、轨底和轨腰三种不同伤损位置钢轨振动信号的 IMFl的PSD分别如图8、图9和图10。最后,将四种情况的PSD合并在一幅图中,如图11。在图11中,从上到下分别是无损,伤在轨头、轨底和轨腰的PSD图,从图中的PSD 对比可以看出,这些信号的PSD中都包含了原无损钢轨的固有振动频率。同时,从伤在轨头、轨底和轨腰的PSD图中,看出伤损信号的频率范围是3580Hz-3920Hz。3880Hz是轨头的有伤的特征频率,其幅值接近于无损钢轨特征频率3360Hz的幅值。当伤损发生在轨底时,其特征频率同样也是接近于3880Hz,但其幅值是伤损发生在轨头的特征频率的3倍。 3600Hz是伤损发生在轨腰时的特征频率,其幅值是伤发生在轨底的频率的4倍。从以上规律中可以看出,无伤以及伤损发生在三种不同位置时的信号的特征频率及其特征规律。因此,通过以上规律我们能建立相应的伤损的存在及其位置的判断准则。即通过对振动信号的EMD分解,提取包含伤损特征信息的IMF并对其进行PSD计算,然后从PSD图观察其特征频率并确定是否发生了伤损及其发生的位置。从该实例的PSD图可以看出,很好的识别出了不同位置伤损的特征频率,确定了伤损发生的位置。综上,对于高铁钢轨的伤损探伤,本发明给出了利用振动信号建立伤损判断准则的详细方法及其装置,能较好的实现对整个路段的实时探测,确定伤损的存在及其位置,提高了探伤效率,以保证高速列车的正常运行。
权利要求
1.一种基于振动信号的高铁钢轨伤损探测方法,其特征在于,该方法包括以下步骤 步骤一、通过振动加速度传感器采集高铁钢轨的振动信号,并对所述振动信号进行经验模态分解,获取η个IMF分量和一个残差;IMF为内固模态函数;步骤二、计算步骤一获取的η个IMF分量与分解之前采集到的振动信号的相关性系数 Ui (i = 1,2,...,11),选取相关性系数1^大于阈值λ的IMF分量作为伤损特征信息IMF ; 步骤三、分别计算步骤二获取的伤损特征信息IMF的功率谱密度, 步骤四、根据步骤三获取的伤损特征信息IMF的功率谱密度确定高铁钢轨伤损部位, 完成高铁钢轨伤损的探测。
2.根据权利要求1所述的一种基于振动信号的高铁钢轨伤损探测方法,其特征在于, 步骤一中通过振动加速度传感器采集高铁钢轨的振动信号,将振动加速度传感器设置在高铁钢轨的轨底上表面,且距离轨底外边缘d距离处,所述d的取值范围是18mm 20mm。
3.根据权利要求1所述的一种基于振动信号的高铁钢轨伤损探测方法,其特征在于, 步骤一中获取η个IMF分量和一个残差的过程为设定输入伤损振动信号为x(t),t = 1,2,· · ·,N,步骤11、IMF分解过程初始化η= 1,且满足关系式IV1U) =x(t)成立,其中ivJt) 为第(n-1)次分解后趋势函数;步骤12、筛选过程初始化,k= 1,且满足关系式hn(k_D (t) =Tn^1 (t)成立,其中hndt) 为第η次经验模态分解中经过第(k-Ι)次筛选后的剩余函数;步骤13、根据筛选程序获取经过第k次筛选后的剩余函数、⑴;步骤14、采用标准偏差准则判断步骤13获得的剩余函数hnk(t)是否满足本征模态函数IMF的条件,即(hn(k_O (O - hnk (O)2 / H2nik^ 0)是否小于阈值Hsd, 0. 2 彡 Hsd 彡 0. 3 ;判断结果为是,执行步骤15,判断结果为否,则k = k+Ι,然后执行步骤13, 步骤15、提取一个本征模态函数分量Cn (t) = hnk(t);步骤16、获取输入伤损振动信号x(t)经过第η次经验模态分解的趋势函数rn(t)= Iv1 (t)-cn (t);步骤17、判断趋势函数rn(t)是否为单调函数;判断结果为否,则η = η+1,然后执行步骤12,判断结果为是,完成提取过程,获得η个 IMF 分量{Cl(t), c2(t)-cn(t)};和 1 个残差 RES:rn(t)。
4.根据权利要求3所述的一种基于振动信号的高铁钢轨伤损探测方法,其特征在于, 步骤13中获取剩余函数hnk(t)的过程为步骤131、利用三次样条函数获取输入伤损振动信号x(t)经过第η次经验模态分解中经过第k-Ι次筛选后的剩余函数Kari) (t)的上、下包络,步骤132、计算所述剩余函数hn(k_D (t)上、下包络曲线在各个t的均值‘(^0), 步骤133、获取输入伤损振动信号x(t)经过第η次经验模态分解中经过第k次筛选后的剩余函数Κ ) = (0-^-1)(0。
5.根据权利要求3所述的一种基于振动信号的高铁钢轨伤损探测方法,其特征在于, 步骤14的中Hsd = 0. 25。
6.根据权利要求1所述的一种基于振动信号的高铁钢轨伤损探测方法,其特征在于, 步骤二中的相关性系数Ui按如下公式获取
7.根据权利要求1所述的一种基于振动信号的高铁钢轨伤损探测方法,其特征在于, 步骤二中的阈值λ = 0. 4 0. 6。
8.根据权利要求1所述的一种基于振动信号的高铁钢轨伤损探测方法,其特征在于, 步骤三中的伤损特征信息IMF的功率谱密度按如下公式获取
9.根据权利要求1所述的一种基于振动信号的高铁钢轨伤损探测方法,其特征在于, 步骤四中确定高铁钢轨伤损部位的判断准则为若包含伤损特征信息IMF的功率谱密度曲线中存在频率为3880Hz 士 IOOHz的频率点, 且该频率点幅值处于无损钢轨特征频率3360Hz的幅值上下浮动5%范围内,则判定为高铁钢轨伤损部位为轨头;若包含伤损特征信息IMF的功率谱密度曲线中存在频率为3880Hz士 IOOHz的频率点, 且该频率点幅值处于无损钢轨特征频率3360Hz的幅值的300% 350%范围内,则判定为高铁钢轨伤损部位为轨底;若包含伤损特征信息IMF的功率谱密度曲线中存在频率为3600Hz士 IOOHz的频率点且该频率点幅值处于无损钢轨特征频率3360Hz的幅值的1200% 1300%范围内,则判定为高铁钢轨伤损部位为轨腰。
10.实现权利要求1所述的一种基于振动信号的高铁钢轨伤损探测方法的装置,其特征在于,它包括振动加速度传感器(1)、信号转换电路O)、经验模态分解模块(3)、伤损特征信息提取模块G)、伤损特征信息的功率密度获取模块(5)和高铁钢轨的伤损部位判定模块(6),设置在高铁钢轨上的振动加速度传感器⑴检测高铁钢轨的伤损信息,并输出给信号转换电路O),信号转换电路( 将接收的模拟信号转换为数字形式的输入伤损振动信号为x(t)并输出给经验模态分解模块(3),经验模态分解模块C3)将输入伤损振动信号为χ (t)进行处理,将获取的η个IMF分量输出给伤损特征信息提取模块(4),伤损特征信息提取模块(4)将满足伤损条件的IMF分量输出给伤损特征信息的功率密度获取模块(5), 伤损特征信息的功率密度获取模块( 计算出每个包含伤损特征的IMF分量的功率密度, 并输出给高铁钢轨的伤损部位判定模块(6),高铁钢轨的伤损部位判定模块(6)根据各IMF 分量的功率密度判断被测高铁钢轨的伤损部位。
全文摘要
一种基于振动信号的高铁钢轨伤损探测方法及装置,属于铁路安全监测与防护技术领域,本发明为解决现有的测量设备检测高铁钢轨伤损情况时会严重影响列车运行效率的问题。本发明方法包括步骤一、通过振动加速度传感器采集高铁钢轨的振动信号,并对所述振动信号进行经验模态分解,获取n个IMF分量和一个残差;步骤二、计算步骤一获取的n个IMF分量与分解之前采集到的振动信号的相关性系数ui(i=1,2,...,n),选取相关性系数ui大于阈值λ的IMF分量作为伤损特征信息IMF;步骤三、分别计算步骤二获取的伤损特征信息IMF的功率谱密度,步骤四、根据步骤三获取的伤损特征信息IMF的功率谱密度确定高铁钢轨伤损部位,完成高铁钢轨伤损的探测。
文档编号G01N29/04GK102175768SQ201110042898
公开日2011年9月7日 申请日期2011年2月22日 优先权日2011年2月22日
发明者沈毅, 王艳, 章欣 申请人:哈尔滨工业大学
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1