一种基于Sobel梯度的直方图均衡算法及其在红外图像上的应用

万昕, 刘坤, 崔昌浩

万昕, 刘坤, 崔昌浩. 一种基于Sobel梯度的直方图均衡算法及其在红外图像上的应用[J]. 红外技术, 2024, 46(4): 452-459.
引用本文: 万昕, 刘坤, 崔昌浩. 一种基于Sobel梯度的直方图均衡算法及其在红外图像上的应用[J]. 红外技术, 2024, 46(4): 452-459.
WAN Xin, LIU Kun, CUI Changhao. Histogram Equalization Algorithm Based on Sobel Gradient and Its Application on Infrared Images[J]. Infrared Technology , 2024, 46(4): 452-459.
Citation: WAN Xin, LIU Kun, CUI Changhao. Histogram Equalization Algorithm Based on Sobel Gradient and Its Application on Infrared Images[J]. Infrared Technology , 2024, 46(4): 452-459.

一种基于Sobel梯度的直方图均衡算法及其在红外图像上的应用

详细信息
    作者简介:

    万昕(1993-),男,湖北广水人,硕士,主管技师,主要从事图像算法的研究

    通讯作者:

    刘坤(1992-),男,湖北武汉人,硕士,算法工程师,主要从事图像算法的研究。E-mail: 271581673@qq.com

  • 中图分类号: TN911.73

Histogram Equalization Algorithm Based on Sobel Gradient and Its Application on Infrared Images

  • 摘要: 为了能在动态范围压缩的同时增强红外图像的对比度,提出了一种基于Sobel梯度直方图均衡算法(gradient histogram equalization,GHE)。与以往的直方图均衡化(histogram equalization,HE)方法不同,该方法自适应地为图像强梯度的灰阶分配高对比度,保留并增强16 bit图像中更多的细节。随后使用双Gamma映射对映射曲线进行调整,有效地抑制图像亮部的过曝现象,同时提高暗部的细节。该方法相比于传统的直方图均衡化算法在暗区细节处理、过曝抑制、对比度增强等方面都有较好的效果。
    Abstract: A Sobel gradient histogram equalization(GHE) algorithm is proposed to enhance the contrast of infrared images during dynamic range compression. In contrast to previous histogram equalization(HE) methods, this method adaptively assigns a high contrast to the strongly graded parts of the image, preserving and enhancing more details in the 16-bit image. Dual Gamma mapping is then used to adjust the mapping curve to effectively suppress overexposure in the bright parts of the image while improving the detail in the shadows. Compared with the traditional histogram equalization algorithm, this method has better effects on dark area detail processing, overexposure suppression, and contrast enhancement.
  • 区别于传统的波纹管自调形式,记忆合金型自调式制冷器采用了另外一条自调技术路线,即利用记忆材料随温度变化产生的形状记忆效应来实现制冷器流量的调节,具有结构简单、装调方便、寿命长等优点[1]。记忆合金制冷器自调后,流量会趋于稳定,但在一些因素如振动、温度变化、制冷器记忆合金的不稳定等的影响下,稳定流量可能会发生较大变化,造成制冷器不到温或工作时间缩短。关于记忆合金自调式制冷器流量稳定性的研究,学者多从提升记忆合金弹簧的疲劳稳定性和设计预紧力来维持制冷器流量的稳定[2-5],很少有学者研究节流孔大小对制冷器流量稳定性的影响。

    由于焦耳汤姆逊效应存在,高压气体流经节流孔,压力明显降低,相应地温度大幅下降直至相变为液体,液态制冷工质蒸发后吸热来实现被冷却对象的制冷。

    节流孔直径的设计对制冷器性能有着重要的影响:直径偏小,制冷流量也偏小,制冷量不足;直径偏大,制冷工质液化率低,制冷效率不高,且探测器工作时长会受到影响。为满足制冷量和探测器工作时长的要求,节流孔直径将被限制在一个范围内。对于对流量稳定性要求较高的自调式制冷器,考虑节流孔直径对流量稳定性的影响,可以进一步寻找到一个更佳的节流孔直径。本文对一款典型的记忆合金自调式制冷器,自调机构见图 1,对其节流孔孔径在0.10~0.25 mm之间变化时的流量稳定性进行了理论和实验研究,探讨了节流孔孔径对制冷器流量稳定性的影响[6]

    图  1  记忆合金自调式制冷器结构图
    Figure  1.  The sketch of a SMA self-regulated cryocooler

    图 1是一款典型的记忆合金自调式制冷器结构,主要包括主动弹簧、形状记忆合金调节器、补偿块、平衡弹簧和阀针。制冷器通气后,高压气体先大流量流经节流孔,经节流后温度骤降,形状记忆合金弹簧被冷却收缩,主动弹簧和平衡弹簧相应地伸长,带动阀针运动关小节流孔实现流量的自动调节。

    现从理论层面分析,改变节流孔孔径大小是否有助于提高自调制冷器流量稳定性。图 2为节流孔和阀针结构示意图,气流由1断面流向2断面。

    图  2  节流孔和阀针结构示意图
    Figure  2.  The sketch of the orifice and needle

    为简化计算,现作两点假设:

    1)气流由1断面流到2断面的沿程阻力损失忽略不计;

    2)气流相变发生在2断面之后,1-2断面制冷工质保持气体状态。

    现通过能量守恒来计算制冷工质流量。由于阀针阻碍导致的局部阻力损失系数为[7]

    $$\zeta {\rm{ = }}0.5\left( {1 - \frac{{{A_2}}}{{{A_1}}}} \right)$$ (1)

    式中:A1A2分别为1、2断面的过流面积。

    由能量守恒知:

    $${P_1} = 0.5(1 - \frac{{{A_2}}}{{{A_1}}})\frac{{\rho v_2^2}}{2} + {P_{{\rm{s2}}}} + {P_{{\rm{d2}}}}$$ (2)

    式中:P1为断面1的全压;等式右侧第一项为断面1~2之间的局部压力损失,ρ为气流的密度;v2为断面2处的流速;Ps2为断面2的静压,与外界大气相通,接近于大气压,取值为0;Pd2为断面2处的动压,其表达式为:

    $${P_{{\rm{d}}2}} = \frac{{\rho v_2^2}}{2}$$ (3)

    流过断面2的流量为:

    $$ Q = {A_2}{v_2} $$ (4)

    式(2)、(3)、(4)联立:

    $${P_1} = 0.5(1 - \frac{{{A_2}}}{{{A_1}}})\frac{{\rho {Q^2}}}{{2A_2^2}} + \frac{{\rho {Q^2}}}{{2A_2^2}}$$ (5)
    $$Q = \sqrt {\frac{{A_2^2{P_1}}}{{\frac{\rho }{2} + \frac{\rho }{4}(1 - \frac{{{A_2}}}{{{A_1}}})}}} $$ (6)

    断面2面积A2可用节流孔孔径R1、阀针进入节流孔的深度l、阀针的角度α表示:

    $$Q = \sqrt {\frac{{{{\left[ {{\rm{ \mathsf{ π} }}R_1^2 - {\rm{ \mathsf{ π} }}{{(l*\tan \frac{\alpha }{2})}^2}} \right]}^2}{P_1}}}{{\frac{\rho }{2} + \frac{{\rho {{(l*\tan \frac{\alpha }{2})}^2}}}{{4R_1^2}})}}} $$ (7)

    由流量公式可知,影响流量的主要因素是节流孔孔径、阀针进入节流孔的距离、阀针的角度、进气压力、气流密度。

    本研究通过计算制冷器流量变化量来评估不同节流孔孔径制冷器的流量稳定性。对于一种特定的自调式制冷器,其调试流量一般都会设定在某一区间,计算中控制不同直径的制冷器调试流量均相同,为14.61 g/min;设扰动因素会导致阀针进入节流孔中的距离减少0.01 mm;制冷工质为氮气,节流前密度为506.25 kg/m3;1断面处的压力为27 MPa,阀针角度为30°,计算结果如表 1图 3所示。

    表  1  节流孔孔径对流量稳定性影响算例
    Table  1.  Example of the influence of orifice diameter on flow stability
    Orifice diameter/mm Setup flow rate/(g/min) Distance of needle into orifice at setup flow rate/mm Distance of needle into orifice after disturbance/mm Flow rate after disturbance/(g/min) Flow rate variation/(g/min)
    0.10 14.61 0.1634 0.1534 21.07 6.46
    0.11 14.61 0.1840 0.1740 21.80 7.19
    0.12 14.61 0.2042 0.1942 22.52 7.91
    0.13 14.61 0.2242 0.2142 23.23 8.62
    0.14 14.61 0.2439 0.2339 23.94 9.33
    0.15 14.61 0.2635 0.2535 24.65 10.04
    0.16 14.61 0.2829 0.2729 25.35 10.74
    0.17 14.61 0.3022 0.2922 26.05 11.44
    0.18 14.61 0.3214 0.3114 26.74 12.13
    0.19 14.61 0.3406 0.3306 27.44 12.83
    0.20 14.61 0.3597 0.3497 28.14 13.53
    0.21 14.61 0.3787 0.3687 28.83 14.22
    0.22 14.61 0.3977 0.3877 29.52 14.91
    0.23 14.61 0.4166 0.4066 30.22 15.61
    0.24 14.61 0.4356 0.4256 30.91 16.30
    0.25 14.61 0.4544 0.4444 31.60 16.99
    下载: 导出CSV 
    | 显示表格
    图  3  流量变化量和节流孔孔径的关系
    Figure  3.  The correlation between flow rate variation and orifice diameter

    表 1数据显示,在14.61 g/min的调试流量下,扰动因素导致阀针进入节流孔的距离减少0.01 mm时,节流孔孔径为0.25 mm的制冷器流量增加了16.99 g/min,而节流孔孔径为0.10 mm的制冷器流量增加了6.46 g/min,小于0.25 mm孔径制冷器10.53 g/min;图 3显示制冷器流量变化量随节流孔孔径增大呈线性增大趋势;说明节流孔孔径越小,制冷器流量越不容易发生变化,流量稳定性越好。

    制冷器调试流量既要满足制冷量要求,又不能超差,要处于一个合理的范围内。现分析调试流量大小对流量稳定性的影响。节流孔孔径设定为0.15 mm;调试流量在14.61~23.65 g/min之间变化;干扰因素相同,均使得阀针在节流孔中的距离减少0.01 mm;计算结果见表 2图 4

    表  2  调试流量对流量稳定性影响算例
    Table  2.  Example of the influence of adjusting flow rate on flow stability
    Orifice diameter/mm Setup flow rate/(g/min) Distance of needle into orifice at setup flow rate/mm Distance of needle into orifice after disturbance/mm Flow rate after disturbance/(g/min) Flow rate variation/(g/min)
    0.15 14.61 0.2635 0.2535 24.65 10.04
    0.15 15.62 0.2625 0.2525 25.64 10.02
    0.15 16.63 0.2615 0.2515 26.63 10.00
    0.15 17.64 0.2605 0.2505 27.63 9.99
    0.15 18.64 0.2595 0.2495 28.62 9.97
    0.15 19.65 0.2585 0.2485 29.60 9.96
    0.15 20.65 0.2575 0.2475 30.59 9.94
    0.15 21.65 0.2565 0.2465 31.58 9.93
    0.15 22.65 0.2555 0.2455 32.56 9.91
    0.15 23.65 0.2545 0.2445 33.54 9.89
    下载: 导出CSV 
    | 显示表格
    图  4  流量变化量和调试流量的关系
    Figure  4.  The correlation between flow rate variation and adjusting flow rate

    表 2数据显示,在14.61 g/min的调试流量下,受扰动因素影响,制冷器流量变化量为10.04 g/min,当调试流量增大至23.65 g/min时,制冷器流量变化量为9.89 g/min,流量变化量减少了0.15 g/min;图 4中显示流量变化量随调试流量的增大呈线性减少趋势。整体来看,调试流量增加,制冷器的流量变化量会减小,但减小幅度不大,而调试流量调的过大,很容易造成制冷器流量超差,增大调试流量对提高制冷器的流量稳定性的作用较为有限。

    理论分析表明:减小节流孔直径有助于提高制冷器的流量稳定性,改变调试流量对制冷器流量稳定性的影响较小。

    制冷器在受到自身或外界因素变化的影响下,流量会发生变化。为引入扰动因素,实验中对制冷器进行了疲劳测试和振动测试。记忆合金弹簧由于其自身材料的特点存在疲劳稳定性的问题,其低温下收缩量不稳定会影响阀针进入节流孔的距离,导致制冷器流量的变化,疲劳测试可以反映记忆合金弹簧不稳定对制冷器流量稳定性的影响;振动测试是为了模拟制冷器机动过程中受到的加速度冲击,制冷器受加速度冲击后,自调机构之间的相对位置会发生一定的变化,从而导致制冷器流量发生变化。

    疲劳测试台如图 5所示,疲劳测试设备一端与气源连接,另一端连有10个接口,可同时供10只制冷器测试。用户在控制台的可视化界面中输入工作时间、停机时间及运行次数,程序根据输入参数控制气源的输送和切断。测试在恒温、恒湿的洁净间中进行,温度为22℃、湿度为46%、净化等级为10万级。测试中,通过调节减压阀,将供气压力调节至29 MPa,工作时间设为5 min,停机时间设为15 min,运行次数设为100次。在工作时间内,程序打开气源开关,向制冷器输送高压气体,气体节流制冷,记忆合金弹簧被冷却后收缩;在停机时间内,程序关闭气源开关,停止向制冷器输送高压气体,制冷器无冷量输出,记忆合金弹簧逐渐恢复至原长,在100次的运行次数下,记忆合金弹簧经历100次的疲劳变形。

    图  5  疲劳测试实验台
    Figure  5.  Experimental platform for fatigue test

    振动测试台如图 6所示,功能振动功率谱密度如图 7所示。试验中,制冷器固定于制冷器装卡夹具中,振动过程中制冷器全程通气,气源压力29 MPa。振动频率为20~2000 Hz,最大功率谱密度为0.04 g2/Hz,总加速度均方根值为7.68 g,振动方向为振动台轴向,振动时间为10 min。

    图  6  振动测试实验台
    Figure  6.  Experimental platform for vibration test
    图  7  功能振动功率谱密度
    Figure  7.  Functional vibration power spectral density

    记忆合金自调式制冷器节流孔直径通常在0.10~0.25 mm之间。为使对比明显,应选择直径跨度较大的节流孔。考虑到0.10 mm附近的节流孔加工难度大,精度较难保证,因此选择直径为0.15 mm和0.25 mm两种节流孔用实验研究。实验中,0.15 mm和0.25 mm两种节流孔规格的制冷器各制作9只,由于将制冷器的流量调至完全相同是非常困难的,实验中将制冷器在29 MPa进行调试,调试流量保持在15~18 g/min的小区间变化。对制冷器调试后、疲劳测试后、振动后的流量进行测试,测试压力为29 MPa和22 MPa,共有6种不同的工况。所用流量计为质量流量计,如图 8所示。测得的流量数据如表 3表 4所示。

    图  8  质量流量计
    Figure  8.  Mass flow meter
    表  3  0.15 mm节流孔制冷器流量数据
    Table  3.  Flow rate data of cryocoolers with 0.15 mm orifice
    Cryocooler number Setup flow rate/(g/min) Flow rate after fatigue test/(g/min) Flow rate after vibration test /(g/min)
    29 MPa 22 MPa 29 MPa 22 MPa 29 MPa 22 MPa
    A1 17.35 14.16 15.61 11.4 16.75 15.07
    A2 16.47 12.32 16.9 14.3 15.73 14.57
    A3 16.5 13.25 16.18 14.6 15.89 13.87
    A4 16.15 12.38 14.32 11.83 11.92 8.13
    A5 16.45 13.66 17.17 15.7 18.04 15.1
    A6 15.16 11.85 16.45 13.95 16.19 13.56
    A7 15.84 13.35 15.04 13.11 16.07 12.9
    A8 16.67 12.36 18.07 14.75 17.99 14.05
    A9 15.44 13.68 15.19 13.01 16.06 15.34
    下载: 导出CSV 
    | 显示表格
    表  4  0.25 mm节流孔制冷器流量数据
    Table  4.  Flow rate data of cryocoolers with 0.25 mm orifice
    Cryocooler number Setup flow rate/(g/min) Flow rate after fatigue test/(g/min) Flow rate after vibration test /(g/min)
    29 MPa 22 MPa 29 MPa 22 MPa 29 MPa 22 MPa
    B1 16.63 12.33 15.1 12.1 17.13 13.56
    B2 15.28 12.97 15.9 10.54 15.84 9.89
    B3 17.8 14.43 17.91 11.07 18.33 12.32
    B4 15.35 11.67 16.4 13.34 16.17 13.57
    B5 16.92 15.3 14.31 11.61 21.42 18.87
    B6 16.32 14.6 17.2 13.41 21.13 17.23
    B7 17.87 15.33 16.66 14.12 16.42 13.89
    B8 15.19 13.3 19.4 17.4 12.41 5.77
    B9 16.77 14.62 16.77 14.62 18.04 15.46
    下载: 导出CSV 
    | 显示表格

    为比较两种不同规格制冷器的流量稳定性,现将调试后流量作为基准流量,记为Q0,疲劳测试后流量和振动测试后流量记为Q1Q2,计算制冷器的流量方差σ和方差均值$\bar \sigma $,如式(8)和式(9)所示。计算制冷器在29 MPa和22 MPa下的流量方差,结果如表 5所示。

    $$\sigma = \frac{{{{({Q_1} - {Q_0})}^2} + {{({Q_2} - {Q_0})}^2}}}{2}$$ (8)
    $$\bar \sigma = \frac{{\sum\limits_{i = 1}^9 {{\sigma _i}} }}{9}$$ (9)
    表  5  0.15 mm/0.25 mm孔径制冷器流量方差
    Table  5.  Flow rate variance of cryocoolers with 0.15 mm and 0.25 mm orifice
    Cryocooler number Flow rate variance under 29 MPa Flow rate variance under 22 MPa Cryocooler number Flow rate variance under 29 MPa Flow rate variance under 22 MPa
    A1 1.6938 4.2229 B1 1.29545 0.7829
    A2 0.3662 4.4915 B2 0.349 7.69565
    A3 0.2373 1.1035 B3 0.1465 7.87085
    A4 10.6209 9.1825 B4 0.88745 3.19945
    A5 1.5233 3.1176 B5 13.53105 13.1805
    A6 1.3625 3.6671 B6 11.95525 4.1665
    A7 0.3465 0.1301 B7 1.7833 1.76885
    A8 1.8512 4.2841 B8 12.72625 36.75545
    A9 0.2235 1.6023 B9 0.80645 0.3528
    Avarage 2.6932 3.5335 Avarage 4.8311 8.4192
    下载: 导出CSV 
    | 显示表格

    从方差均值来看:0.15 mm节流孔的制冷器在29 MPa和22 MPa下的流量方差分别为2.6932和3.5335,而0.25 mm节流孔的制冷器在两种压力下的方差均值分别为4.8311和8.4192,明显大于节流孔直径为0.15 mm的情况。实验结果表明:0.15 mm节流孔孔径的制冷器在29 MPa和22 MPa下的流量方差波动均较小,不容易发生较大的流量变化。

    流量稳定性是评价自调型制冷器性能的一个重要指标。制冷器流量增大会导致工作时长缩短,流量减小会导致制冷量不足。本文从理论分析和实验研究两种方法出发,探究了节流孔孔径对记忆合金自调式制冷器的流量稳定性的影响。研究表明:在满足其他使用要求的情况下,节流孔孔径设计地越小,制冷器的流量稳定性越好。本文研究内容有助于记忆合金型自调式制冷器设计优化。

  • 图  1   直方图与梯度直方图对比

    Figure  1.   Comparison of histogram with gradient histogram

    图  2   室内梯度直方图效果

    Figure  2.   Result of GHE

    图  3   不同gamma值下的映射曲线

    Figure  3.   Mapping curves with different gamma values

    图  4   直方图均衡化和梯度直方图均衡化效果对比

    Figure  4.   Comparison of HE result with GHE result

    图  5   噪声抑制效果对比

    Figure  5.   Comparison of noise suppression result

    图  6   室外双Gamma映射效果对比

    Figure  6.   Comparison of dual Gamma mapping

    图  7   室内双Gamma映射效果对比

    Figure  7.   Comparison of dual Gamma mapping

    表  1   两种算法对比度对比

    Table  1   Comparison of HE with GHE

    Algorithms Sky Architecture Picture
    HE 3.74 54.33 74.14
    GHE 0.19 123.96 98.19
    下载: 导出CSV
  • [1]

    LIU C, SUI X, KUANG X, et al. Optimized contrast enhancement for infrared images based on global and local histogram specification[J]. Remote Sensing, 2019, 11(7): 849. DOI: 10.3390/rs11070849

    [2] 李凌杰, 陈菲菲. 基于改进直方图的红外图像增强方法[J]. 航空兵器, 2022, 29(2): 101-105. https://www.cnki.com.cn/Article/CJFDTOTAL-HKBQ202202015.htm

    LI Lingjie, CHEN Feifei. Infrared image enhancement method based on improved histogram[J]. Aero Weaponary, 2022, 29(2): 101-105. https://www.cnki.com.cn/Article/CJFDTOTAL-HKBQ202202015.htm

    [3]

    CHEN X, LV L. A compositive contrast enhancement algorithm of IR image[C]//International Conference on Information Technology & Applications of IEEE, 2013, DOI: 10.1109/ITA.2013.20

    [4]

    LIU N, CHEN X. Infrared image detail enhancement approach based on improved joint bilateral filter[J]. Infrared Physics & Technology, 2016, 77: 405-413.

    [5]

    HUANG Z, ZHANG T, LI Q, et al. Adaptive gamma correction based on cumulative histogram for enhancing near-infrared images[J]. Infrared Physics & Technology, 2016, 79: 205-215.

    [6]

    JUNG C, LU W. Tone-preserving contrast enhancement in images using rational tone mapping and constrained optimization[C/OL]. 2016 IEEE Visual Communications and Image Processing, https://ieeexplore.ieee.org/document/7805525.

    [7]

    LI S, JIN W, LI L, et al. An improved contrast enhancement algorithm for infrared images based on adaptive double plateaus histogram equalization[J]. Infrared Physics & Technology, 2018, 90: 164-174.

    [8]

    Hum Y C, Tee Y K, W. -S. Y, et al. A contrast enhancement framework under uncontrolled environments based on just noticeable difference[J]. Signal Processing, 2022, 103: 116657

    [9]

    Bae J, Kim K, Yun Y J, et al. Adaptive tone-mapping operator for HDR images based on image statistics[C]//Proc. of IEEE, 2011, DOI: 10.1109/TENCON.2011.6129047

    [10]

    DENG Weitao, LIU Lei, CHEN Huateng, et al. Infrared image contrast enhancement using adaptive histogram correction framework[J]. Optik, 2022, 271: 170114. DOI: 10.1016/j.ijleo.2022.170114

    [11]

    Arici T, Dikbas S, Altunbasak Y. A histogram modification framework and its application for image contrast enhancement[J]. IEEE Trans Image Process, 2009, 18(9): 1921-1935. DOI: 10.1109/TIP.2009.2021548

    [12]

    Ok J, Lee C. HDR tone mapping algorithm based on difference compression with adaptive reference values[J]. Journal of Visual Communication and Image Representation, 2016, 43: 61-76.

    [13]

    Arici T, Dikbas S, Altunbasak Y. A histogram modification framework and its application for image contrast enhancement[J]. IEEE Trans Image Process, 2009, 18(9): 1921-1935. DOI: 10.1109/TIP.2009.2021548

    [14]

    Lee S. An efficient content-based image enhancement in the compressed domain using Retinex theory[J]. IEEE Transactions on Circuits & Systems for Video Technology, 2007, 17(2): 199-213.

  • 期刊类型引用(1)

    1. 张巍,苏鹏翼,高涛,贾宝军. 弹簧参数对记忆合金式J-T制冷器流量的影响. 低温与超导. 2022(08): 93-99 . 百度学术

    其他类型引用(0)

图(7)  /  表(1)
计量
  • 文章访问数:  162
  • HTML全文浏览量:  39
  • PDF下载量:  65
  • 被引次数: 1
出版历程
  • 收稿日期:  2023-04-08
  • 修回日期:  2024-04-07
  • 刊出日期:  2024-04-19

目录

/

返回文章
返回
x 关闭 永久关闭

尊敬的专家、作者、读者:

端午节期间因系统维护,《红外技术》网站(hwjs.nvir.cn)将于2024年6月7日20:00-6月10日关闭。关闭期间,您将暂时无法访问《红外技术》网站和登录投审稿系统,给您带来不便敬请谅解!

预计6月11日正常恢复《红外技术》网站及投审稿系统的服务。您如有任何问题,可发送邮件至编辑部邮箱(irtek@china.com)与我们联系。

感谢您对本刊的支持!

《红外技术》编辑部

2024年6月6日