红外偏振成像系统性能评估模型

王霞, 赵家碧, 孙琪扬, 金伟其

王霞, 赵家碧, 孙琪扬, 金伟其. 红外偏振成像系统性能评估模型[J]. 红外技术, 2023, 45(5): 437-445.
引用本文: 王霞, 赵家碧, 孙琪扬, 金伟其. 红外偏振成像系统性能评估模型[J]. 红外技术, 2023, 45(5): 437-445.
WANG Xia, ZHAO Jiabi, SUN Qiyang, JIN Weiqi. Performance Evaluation Model for Infrared Polarization Imaging System[J]. Infrared Technology , 2023, 45(5): 437-445.
Citation: WANG Xia, ZHAO Jiabi, SUN Qiyang, JIN Weiqi. Performance Evaluation Model for Infrared Polarization Imaging System[J]. Infrared Technology , 2023, 45(5): 437-445.

红外偏振成像系统性能评估模型

基金项目: 

国家自然科学基金资助项目 62171024

详细信息
    作者简介:

    王霞(1972-),女,副教授,博士,主要从事图像处理、红外偏振成像、光电探测等方向的教学和研究工作。E-mail: angelniuniu@bit.edu.cn

  • 中图分类号: TN219

Performance Evaluation Model for Infrared Polarization Imaging System

  • 摘要: 红外偏振成像系统快速发展且应用广泛,但评估其性能的成像系统性能模型发展不足。迫切需要能够与先进的偏振成像系统相匹配的性能模型。利用深度学习网络的训练过程与人脑提取认知信息过程的相似性,本文首次将深度学习方法引入系统性能模型领域,提出了一种基于二维图像的可自动评估系统性能的红外偏振成像系统性能模型。该模型主要包含两个主要模块:退化模块、性能感知模块。在评估一个新的系统时,需要输入高质量的原始图像,并根据系统的硬件参数量身定制成像系统退化模块,退化完成后输入性能感知模块,从而得到最终的目标获取性能。为验证模型有效性,本文基于红外辐射理论自建了面向海面场景的红外偏振数据集,训练网络并进行测试。应用该模型对红外偏振成像系统的性能进行评估,评估结果与主观感知具有较好的一致性。
    Abstract: Although infrared polarization imaging systems have been developed rapidly and widely, a model for evaluating their performance has not been sufficiently developed. Performance models that can match advanced polarization imaging systems are urgently required. Regarding the similarity between the training process of a deep learning network and the process of extracting cognitive information from the human brain, this paper introduces a deep learning method in the field of system performance modeling for the first time and proposes a performance model for infrared polarization imaging systems that can automatically evaluate system performance based on two-dimensional images. The model includes two main modules: a degradation module and a performance awareness module. When evaluating a new system, high-quality original images are input and sequentially passed through an imaging system degradation module, customized according to the hardware parameters of the system, and input into a performance awareness module to obtain the final target acquisition performance. Moreover, to verify the effectiveness of the model, we realized a self-built infrared polarization dataset for sea surface scenes based on infrared radiation theory, and trained and tested the networks. The results obtained when the model was applied to evaluate the performance of infrared polarization imaging systems showed good agreement with subjective perception.
  • 自适应光学的基本思想是在光学系统中引入一个表面形状可变的光学元件(波前校正器)和一个波前误差传感器,用波前传感器测量出不断变化的波前误差,根据测量结果由控制系统控制波前校正器对波前误差进行校正。如果这一过程足够快,就可以用不断变化的波前校正量来校正不断变化的动态波前误差。

    夏克-哈特曼波前传感器是目前自适应光学系统中应用最为广泛的实时波前探测器。自适应光学系统将夏克-哈特曼波前传感器所测到的波前畸变信息转化成校正器的控制信号,以实现自适应光学系统的闭环控制。用于补偿大气湍流引入像差的天文自适应光学系统的空间分辨力和时间分辨力应分别与大气相干长度和时间常数相匹配。优秀观测站址的大气相干长度r0和时间常数分别约十几厘米和几毫秒,随着望远镜口径的增大,要求自适应光学系统有更多的子孔径数。以直径8m的望远镜为例,若r0取15cm,则要求夏克-哈特曼传感器在直径方向上至少有53个子孔径。由此给波前处理系统带来的问题是计算量大、实时性要求高。为提高自适应光学系统的控制带宽,必须提高自适应光学系统的波前采样频率和波前校正速度。在系统采样频率一定的情况下,波前处理系统的计算延时会直接影响系统的控制带宽。在当前图像帧结束的情况下,越早给出波前校正量越能更好地补偿波前像差,实现准实时的波前校正。

    2002年中科院光电所研制的2900Hz、61单元的波前处理机[1]共用17片DSP,计算延时为340μs;2011年中科院电子技术研究所研制的基于FPGA的22×22子孔径波前斜率处理器[2]完成一帧所有子孔径斜率计算的延迟时间为0.33μs;2015年中科院电子技术研究所研制的基于PC机的949个子孔径的系统[3]进行2000 Hz的处理,处理延迟低于240μs;2016年中科院长春光机所研制的基于GPU的液晶自适应光学系统的波前处理[4]斜率延时为18μs,其相机像素为240×240,帧率2000帧/s,40×40个子孔径。以上夏克-哈特曼传感器所产生的图像数据大致在每秒100M像素数左右,现代的FPGA都可以在单周期内完成波前斜率计算,因此只要FPGA运行在100MHz的时钟频率下都可以进行实时处理。面对自适应光学对传感器分辨率、帧率要求的不断提高,夏克-哈特曼传感器所产生的图像数据大致在500M/s以上的像素数,传感器为达到传输带宽的要求,都采用多通道(8~10通道)同时传输的方式,FPGA难以运行在如此高的时钟频率。

    本文针对具有高分辨率(1020×1020),高帧速(560帧/s),大规模子孔径数(56×56单元)的夏克-哈特曼传感器,提出一种基于现场可编程门阵列(field-programmable gate array,FPGA)的实时波前处理机,在有效降低硬件资源的前提下,可在当前子孔径数据进入FPGA后延迟10行完成当前子孔径波前处理运算,提高了系统的波前处理速度。

    夏克-哈特曼波前传感器由微透镜阵列和探测器组成,每个微透镜对应探测器上的一个子探测区。当以参考波前(近似平面波)入射到微透镜阵列时,在探测区域形成近似等间距排列的光斑阵列,此时记下光斑质心位置作为参考。当以畸变波前入射时,每个子透镜对应的光斑质心位置相对参考光斑可能产生XY两个方向的移动,移动量即对应于畸变波前在子孔径内的平均斜率。获得畸变波前在子孔径内的平均斜率后,通过波前重构算法(模式法或区域法)即可获得波前相位。

    图 1为夏克-哈特曼波前传感器的工作原理。光波场被阵列透镜分解成许多小孔径光斑,这些小孔径光斑汇聚到夏克-哈特曼探测器焦平面上。假定不同的子孔径的光成像在探测器靶面的不同位置上,任意两个子孔径所产生的邻近像点没有重叠,这样每个子孔径在像平面上都对应的一块专用的成像面积。进入阵列透镜的光束在探测器上形成一阵列衍射光斑。

    图  1  夏克-哈特曼波前传感器结构
    Figure  1.  Structure of Shaker-Hartman sensor

    根据光斑质心的定义可写出离散采样情况下光斑质心的计算公式为:

    $$\begin{array}{*{20}{c}} {{x_c} = \sum\limits_{i,j}^{L,M} {{x_i}{I_{i,j}}} /\sum\limits_{i,j}^{L,M} {{I_{i,j}}} }\\ {{y_c} = \sum\limits_{i,j}^{L,M} {{y_i}{I_{i,j}}} /\sum\limits_{i,j}^{L,M} {{I_{i,j}}} } \end{array}$$ (1)

    式中:xiyi分别为探测器各子孔径位置坐标;Ii, j为第(i, j)个探测器像元灰度值。

    在每一个子孔径内进行光斑质心计算,可以得到参考波前质心(Xc0Yc0)和畸变波前质心(Xc, Yc)。波前斜率gxgy定义为信标光波在xy方向的偏导数在该子孔径上的平均,波前斜率和光斑质心变化量满足如下关系:

    $$ {g_x} = \frac{{({x_c} - {x_0}) \times p}}{f}, {g_y} = \frac{{({y_c} - {y_0}) \times p}}{f} $$ (2)

    式中:p为探测器像元尺寸;f为子透镜焦距。

    综上所述,夏克-哈特曼波前传感器波前处理的核心为子孔径内光斑质心的计算。

    基于Kintex7-XC7K325T的FPGA实时波前处理机组成如图 2所示,夏克-哈特曼传感器将实时采集的图像数据发送给以FPGA为核心的自适应光学波前处理机,FPGA处理完后将数据发送给上位机,同时上位机可以发送控制指令对处理机的参数进行控制。

    图  2  处理机系统组成
    Figure  2.  Composition of processor system

    夏克-哈特曼传感器的FPGA处理硬件结构如图 3所示。采用具有高帧速、大面阵CMOS探测器的EoSens CL型相机,在Camera Link Extended-Full传输模式下以1020×1020分辨率可实现560帧/s的帧率。为了传输高达4.4 Gbit/s的数据,需要3对Camera Link数据链,因此处理机数据接收端采用3片DS90CR288MTD芯片将Camera Link数据转成LVCMOS信号传输给FPGA。

    图  3  处理机硬件结构图
    Figure  3.  Structure of hardware

    本夏克-哈特曼传感器微透镜整列数为56×56,每个子孔径在xy方向上用8bit表示整数坐标位置,16bit表示小数坐标位置,则每秒至少需要传输560×56×56×2×24=84295680bit数据,大大超过了通信常用的异步串行通信(RS232/422)的传输极限。因此在传输波前质心数据时,本文采用Cypress公司CYUSB2014-BZXI与上位机间进行USB3.0数据传输(USB3.0的最大速率5Gbit),为了降低USB数据传输开发复杂性,本文将USB3.0传输设置为单向传输,上位机与波前处理质心计算系统的控制由RS422完成。

    夏克-哈特曼传感器采用的EoSens CL相机工作在1020×1020分辨率,帧频为560帧/s,采用Camera Link Extended-Full模式传输,像素时钟为75MHz。每时钟同时输出10个像素,每个像素为8bit,共80bit,行有效周期为102个像素时钟。其行时序如图 4所示(图中每一个方块代表一个像素,方块中的数字表示像素的行方向坐标,Tx代表第x个时钟周期)。

    图  4  EoSens CMOS探测器一行像素输出时序
    Figure  4.  Pixel output timing of EoSens CMOS detector

    如果实时将并行输入的10个像素转成串行的单个像素进行处理FPGA需要750MHz的处理时钟,大大超过了FPGA能够达到的最高时钟,因此本文设计了一个行缓冲将高速多像素并行输入的图像数据变成低速像素为单位的串行的数据(如图 5所示),FPGA程序设计中将每一行数据按10个像素(80bit)75MHz时钟频率写入两个FIFO(64bit位宽和16bit位宽两个FIFO组成一个80bit位宽的FIFO)中,再以1个像素(8bit)频率FFPGA(75MHz~200MHz,FPGA一般时钟频率)从FIFO中读出,FIFO深度为128。每行数据的处理时间为1020个FFPGA时钟。

    图  5  行缓冲将10像素80bit输入转成串行1像素8bit
    Figure  5.  Line buffer converts 10-pixel 80-bit parallel input into serial 1-pixel 8-bit

    每一个质心计算核心模块如图 6所示,由一个如图 6的行缓冲FIFO和两个乘法累加器分别进行子孔径内xy方向的乘法运算及累加运算组成,质心计算核心模块将FIFO读出点在对应子孔径参数进行运算。该运算为两周期流水运算。质心计算核心模块将计算后该行对应在子孔径的3个乘积累计结果进行输出。图中YPos为当前质心计算核心模块处理行在整个图像中的行位置,XPos为当前处理像素在整个图像中的列位置,Sum输出为质心计算核心模块处理完当前行后输出的行像素值总和,IntgX为每一个子块在行方向的乘法累加值,IntgY为每一个子块在列方向的乘法累加值。

    图  6  质心计算核心模块
    Figure  6.  Module of centroid computing

    质心计算核心模块工作在相机的像素时钟频率下(75MHz),相机行有效时间为Tlv时,质心计算核心模块需要10Tlv的时间来处理一行的数据,因此一个质心计算核心模块是显然不够的,但是如果每一行都有一个质心计算核心模块会浪费大量的FPGA资源,甚至导致FPGA由于资源不够而无法完成布局布线。为了满足功能的需要又不造成FPGA资源的浪费,如图 7所示,本文采用核心模块复用的质心计算方法,设计了共10个质心计算核心模块,通过核心模块选择器分时复用使用这10个质心计算核心模块,实现质心计算核心模块的复用,减少FPGA资源的使用。

    图  7  质心计算核心模块复用
    Figure  7.  Centroid computing core module reuse

    最终FPGA资源使用情况如图 8所示,FPGA资源使用中最重要的几个指标LUT,FF,BRAM的使用量都在30%以下。

    图  8  FPGA使用资源情况
    Figure  8.  Usage of FPGA resources

    使用QuestaSim对图 9实际采集的夏克-哈特曼传感器成像图进行了仿真,仿真时序如图 10所示,每一行子孔径质心位置结果可以在当前子孔径所有行均输入后延时10行时间输出,采用的EoSens CL型相机工作在1020×1020分辨率560帧时,行周期小于1/(560×1020)=1.75μs,因此本系统的计算输出延迟小于17.5μs。同时将QuestaSim仿真结果与Matlab计算的结果进行了对比,结果一致。

    图  9  夏克-哈特曼传感器成像图
    Figure  9.  Shack-Hartman sensor imaging
    图  10  时序仿真结果
    Figure  10.  Timing simulation results

    表 1将本文所采用的系统和方法与所查到文献中其他处理系统与方法进行了比较,比较主要包括子孔径的数量、传感器的数据量和处理延时,本文所采用的系统和方法在子孔径数量最多,数据量最大的情况下,也能够将延时控制得最小。

    表  1  处理能力比较
    Table  1.  Comparison of processing power
    Processor Number of sub-apertures The data volume of the sensor Processing delay
    2002 Institute of Optics and Electronics, Chinese Academy of Sciences[1] 17DSP 61 7.8MB/s 340μs
    2011 Institute of Optics and Electronics, Chinese Academy of Sciences[5] FPGA 177 No data 338μs
    2014 Technische Universität Ilmenau[6] FPGA 120 43.3MB/s 1050μs
    2015 Institute of Optics and Electronics, Chinese Academy of Sciences[2] Computer 949 88.2MB/s 163μs
    2016 Changchun Institute of Optics and Mechanics[4] GK104 GPU+CPU 1600 115.2MB/s 18μs
    2018 Changchun Institute of Optics and Mechanics[7] FPGA 349 113.5MB/s 198μs
    2018 Indian Institute of Astrophysics[8] FPGA 2500 No data 24-26μs
    System and method of this paper FPGA 3136 580MB/s μs
    下载: 导出CSV 
    | 显示表格

    实验系统光路如图 11所示。光纤耦合输出激光光源(635nm)经准直后通过分束镜1垂直入射到变形镜上。变形镜用于产生畸变波前。畸变波前在分束镜2处分为两路,一路用于夏克-哈特曼传感器进行波前探测,一路用于检测远场光斑的变化。由于变形镜的通光口径一般大于哈特曼传感器的探测区域尺寸,实验系统采用透镜组L3和L4将变形镜出射的光进行缩束,焦距之比等于变形镜口径与哈特曼探测区域直径之比。夏克-哈特曼传感器获得的光斑阵列图像经Cameralink接口发送给FPGA波前处理机进行实时质心计算,质心计算结果由上位机进行波前重构后得到畸变波前分布。

    图  11  实验光路图
    Figure  11.  Experimental light path

    目前夏克-哈特曼传感器所用的高分辨率、高帧率探测器大多采用多通道并行传输的方式才能满足自适应光学所要求的高分辨率、高帧频的需求,本文提出的一种基于FPGA的夏克-哈特曼波前处理机,采用核心模块复用的质心计算方法,兼顾了处理效率和FPGA资源的平衡,可以实时处理1020×1020分辨率、560帧/s图像的成像器件和56×56子孔径数微透镜整列组成的夏克-哈特曼传感器,数据输出延迟小于17.5μs(10个行周期),而且FPGA资源占用率不到30%,由此推算用该FPGA至少可以处理图像分辨率或图像帧率增加1倍的哈特曼传感器。

  • 图  1   基本系统性能模型结构(绿色)。灰色:两个主要模块; 蓝色:性能感知模块

    Figure  1.   Basic system performance model structure(green). Gray: two main components; Blue: performance awareness module

    图  2   PRI-YOLOv5的结构及输出定义示意图

    Figure  2.   Illustration of the PRI-YOLOv5 structure and output definition

    图  3   预测网络示意图

    Figure  3.   Illustration of our prediction network

    图  4   场景仿真示意图

    Figure  4.   Illustration of scene simulation

    图  5   仿真结果示例

    Figure  5.   An example of the simulation results

    图  6   PRI-YOLOv5训练结果

    Figure  6.   Training results of PRI-YOLOv5

    图  7   预测网络测试结果

    Figure  7.   Test results of the prediction network

    图  8   两组数据预测差值与其相似性的关系

    Figure  8.   The difference of the predicted results v.s. NMI

    图  9   几款待评估偏振系统的退化效果示意图

    Figure  9.   Schematic diagram of degradation effects of several polarization systems

    图  10   不同距离处的目标获取概率(三款系统)

    Figure  10.   Target acquisition probability as a function of range(three systems)

    表  1   红外偏振成像系统性能模型研究现状

    Table  1   Research status of performance models for infrared polarization imaging systems

    Name Year Principle Illustration
    Edson Guimaraes[9] 1999 MRTD, Johnson Criterion Based on the MRTD, this model further consider the transfer function and transmittance of the polarizer.
    Mehmet Yildirim[10] 2000 Designed for the second generation forward looking infrared sensor.
    Zhou Chenghao[11] 2013 Emphatically analyzed the models of different targets such as point source targets and extended source targets.
    Xia Runqiu[12] 2016 Ignore registration errors and the impact of the polarizer on MTF, the MRTD model of the polarization system is calculated.
    Liang Jianan[13] 2019 Modify the MRTD model and Johnson criterion based on interference factors such as background clutter.
    下载: 导出CSV

    表  2   船模型及类别

    Table  2   Ship models and classes

    Target class
    Frigate Destroyer Patrol
    Target type ship1 Size
    X: 0.280m
    Y: 1.948m
    ship3 Size
    X: 0.436m
    Y: 3.919m
    ship5 Size
    X: 0.367m
    Y: 2.001m
    ship2 Size
    X: 0.344m
    Y: 2.830m
    ship4 Size
    X: 0.395m
    Y: 3.108m
    ship6 Size
    X: 0.482m
    Y: 2.000m
    Note: All 3D models are downloaded from https://www.3d66.com/
    下载: 导出CSV

    表  3   仿真中变量和常量参数设置

    Table  3   Variables and constant settings during simulation

    Variables
    Group 1: Model ship Group 2: Real ship
    Name Number of values Values(unit) Number of values Values(unit)
    Focal length(f) 2 60, 85(mm) 1 105(mm)
    Pixel size(p) 3 14, 17, 20(μm) 3 14, 17, 20(μm)
    View radius(r) 2 100, 150(m) 4 2, 3, 4, 5(km)
    View zenith angle 4 30, 45, 60, 75(°) 4 30, 45, 60, 75(°)
    Wind speed 1 2(m/s) 1 13(m/s)
    Sea surface size 1 8 m×8 m 1 338 m×338 m
    Target class 6/3 see Table 1 6/3 see Table 1
    Constants
    Name Values Name Values
    Solar zenith angle/° 45 Solar azimuth angle/° 60
    Max reflected number 3 Detected wavelength/μm 3.8
    Image size 300×300 Sea temperature/ K 300
    Deck temperature/K 318.15 Hull temperature/K 303.15
    Hot part temperature/K 328.15 Others temperature/K 313.15
    下载: 导出CSV

    表  4   待评估红外偏振成像系统主要参数

    Table  4   The main parameters of infrared polarization imaging systems to be evaluated

    System Focal length/mm Wave-length/μm F# Resolution Pixel size/μm
    System A 105 3~5.2 3.5 160×120 20
    System B 135 3~5.2 1.5 640×480 17
    System C 105 3~5.2 2 320×240 20
    下载: 导出CSV
  • [1]

    Scott L B, Condiff L R. C2NVEO advanced FLIR systems performance model[C]//Proc. SPIE, 1990, 1309: 168-180.

    [2]

    Scott L B, D'Agostino J A. NVEOD FLIR92 thermal imaging systems performance model[C]//Proc. SPIE, 1992, 1689: 194-203.

    [3]

    Maurer T, Driggers R G, Vollmerhausen R H. 2002 NVTherm improve-ments[C]//Proc. SPIE, 2002, 4719: 15-23.

    [4]

    Edwards T C, Vollmerhausen R H, Driggers R G, et al. NVESD time-limited search model[C]//Proc. SPIE, 2003, 5076: 53-59.

    [5]

    Bijl P, Valeton J M. Triangle orientation discrimination: the alternative to minimum resolvable temperature difference and minimum resolvable contrast[J]. Optical Engineering, 1998, 37: 1976-1983. DOI: 10.1117/1.601904

    [6]

    Vollmerhausen R H, Jacobs E and Driggers R G. New Metric for Predicting target acquisition Performance [J]. Optical Engineering, 2004, 43(11): 2806-2818. DOI: 10.1117/1.1799111

    [7]

    Greif H J, Weiss A R, Wittenstein W. pcSitoS: A new tool for image-based IR system simulation [C]//Proc. SPIE, 2009, 7481: 748107.

    [8]

    Hogervorst M A, Bijl P, Valeton J M. Capturing the sampling effects: a TOD sensor performance model [C]//Proc. SPIE, 2001, 4372: 62-73.

    [9]

    Guimaraes E F. Investigation of minimum resolvable temperature difference formulation for polarized thermal imaging range prediction[D]. America: Naval Postgraduate School, 1999.

    [10]

    Yildirim M. Modeling second generation FLIR sensor detection recognition and identification range with polarization filtering[D]. America: Naval Postgraduate School, 2000.

    [11] 周程灏. 红外偏振成像系统作用距离建模与分析[D]. 哈尔滨: 哈尔滨工业大学, 2013.

    ZHOU Chenghao. Modeling and Analysis of the Operating Range of Infrared Polarization Imaging System[D]. Harbin: Harbin Institute of Technology, 2013.

    [12] 夏润秋, 王霞, 金伟其, 等. 海面环境中红外偏振成像系统作用距离模型[J]. 红外与激光工程, 2016, 45(3): 74-78. https://www.cnki.com.cn/Article/CJFDTOTAL-HWYJ201603010.htm

    XIA Runqiu, WANG Xia, JIN Weiqi, et al. Distance model of infrared polarization imaging system used in sea-surface environment[J]. Infrared and Laser Engineering, 2016, 45(3): 74-78. https://www.cnki.com.cn/Article/CJFDTOTAL-HWYJ201603010.htm

    [13] 梁建安. 基于红外偏振成像技术的水面杂波抑制方法研究[D]. 北京: 北京理工大学, 2019.

    LIANG Jianan. Research on Water Surface Clutter Suppression Method Based on Infrared Polarization Imaging Technology[D]. Beijing: Beijing Institute of Technology, 2019.

    [14]

    ZHANG J Q, WANG X R. Modeling and Performance Evaluation Theory of Photoelectric Imaging System[M]. Xi'an: Xidian University Publishing House, 2010.

    [15]

    Redmon J, Divvala S, Girshick R. You only look once: Unified, real-time object detection[C]//IEEE Conference on Computer Vision and Pattern Recognition, 2016: 779-788.

    [16]

    Jocher G. YOLOv5 release v4.0[EB/OL]. [2023-01-20]. https://github.com/ultralytics/yoloV5.

    [17]

    SU Shaolin, YAN Qingsen, ZHU Yu, et al. Blindly assess image quality in the wild guided by a self-adaptive hyper network[C]//IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2020: 3664-3673.

    [18]

    HE Si, WANG Xia, XIA Runqiu, et al. Polarimetric infrared imaging simulation of a synthetic sea surface with Mie scattering[J]. Applied Optics, 2018, 57(7): B150-B159. DOI: 10.1364/AO.57.00B150

  • 期刊类型引用(1)

    1. 朱强,周维虎,陈晓梅,石俊凯,李冠楠. 高速实时近红外弱信号检测系统. 光学精密工程. 2022(24): 3116-3127 . 百度学术

    其他类型引用(2)

图(10)  /  表(4)
计量
  • 文章访问数:  244
  • HTML全文浏览量:  98
  • PDF下载量:  100
  • 被引次数: 3
出版历程
  • 收稿日期:  2023-04-13
  • 修回日期:  2023-04-27
  • 刊出日期:  2023-05-19

目录

/

返回文章
返回