CN112380723B - 一种菲涅尔衍射分段传播的快速仿真方法 - Google Patents

一种菲涅尔衍射分段传播的快速仿真方法 Download PDF

Info

Publication number
CN112380723B
CN112380723B CN202011341503.7A CN202011341503A CN112380723B CN 112380723 B CN112380723 B CN 112380723B CN 202011341503 A CN202011341503 A CN 202011341503A CN 112380723 B CN112380723 B CN 112380723B
Authority
CN
China
Prior art keywords
propagation
plane
matrix
fresnel diffraction
field matrix
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202011341503.7A
Other languages
English (en)
Other versions
CN112380723A (zh
Inventor
金欣
范逸辉
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shenzhen International Graduate School of Tsinghua University
Original Assignee
Shenzhen International Graduate School of Tsinghua University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Shenzhen International Graduate School of Tsinghua University filed Critical Shenzhen International Graduate School of Tsinghua University
Priority to CN202011341503.7A priority Critical patent/CN112380723B/zh
Publication of CN112380723A publication Critical patent/CN112380723A/zh
Application granted granted Critical
Publication of CN112380723B publication Critical patent/CN112380723B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
    • G06F17/141—Discrete Fourier transforms
    • G06F17/142—Fast Fourier transforms, e.g. using a Cooley-Tukey type algorithm
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/15—Correlation function computation including computation of convolution operations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Algebra (AREA)
  • Computing Systems (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Discrete Mathematics (AREA)
  • Optical Integrated Circuits (AREA)

Abstract

本发明公开一种菲涅尔衍射分段传播的快速仿真方法：A1、输入M*M的包含光源信息的第一全视场矩阵MA，对其建立大小2M*2M的全零的虚拟视场矩阵，将MA放置在虚拟视场矩阵中心；A2、根据当前平面与下一平面间的距离及平面内的像素大小设计光场脉冲响应函数，使用三次快速傅里叶变换进行菲涅尔衍射计算，将当前平面的虚拟视场矩阵传播到下一平面中，得到新的虚拟视场矩阵；A3、构建2M*2M的掩膜矩阵，中心直径为D的圆形部分为1，其余部分为0，将该掩膜矩阵与新的虚拟视场矩阵点乘得到一新矩阵，作为下一平面的虚拟视场矩阵；A4、将下一平面作为当前平面返回A2，重复A2和A3进行分段传播直至到达目标平面，从中心提取M*M的第二全视场矩阵为最终传播结果。

Description

一种菲涅尔衍射分段传播的快速仿真方法
技术领域
本发明涉及光学与数字图像处理领域，具体涉及一种菲涅尔衍射分段传播的快速仿真方法。
背景技术
光学研究中往往需要进行一些仿真实验，在进行光线传播实验中往往会用到菲涅尔衍射传播的仿真。传统的菲涅尔衍射传播仿真中使用的是卷积运算，存在计算时间较长这一问题。之后有团队提出使用三次快速傅里叶变换代替卷积运算进行加速，但该方法需要满足以下条件：光源的最大线性宽度小于或等于用于传播的矩阵边长的一半，光源大小、波长和传播距离之间满足菲涅尔系数，否则会出现频谱混叠现象。通过合理设计光源与传播面的尺寸，该模型可以完成一次从平面1到平面2的菲涅尔传播，但在进行第二次传播的时候，平面2已经被光源所充斥，无法满足上述条件，故不能进行分段传播。
发明内容
本发明的主要目的在于提出一种菲涅尔衍射分段传播的快速仿真方法，通过使用虚拟视场矩阵和掩膜将多次传播时传播面内光源的最大线性宽度进行合理控制，以解决现有的菲涅尔衍射分段传播仿真中存在的频谱混叠问题。
为达上述目的，本发明采用以下技术方案：
一种菲涅尔衍射分段传播的快速仿真方法，包括如下步骤：
A1、输入一个大小为M*M的包含光源信息的第一全视场矩阵，针对所述第一全视场矩阵建立一个全零的大小为2M*2M的虚拟视场矩阵，将输入的所述第一全视场矩阵放置在所述虚拟视场矩阵的中心；
A2、根据当前传播平面与下一传播平面间的距离以及平面内的像素大小设计菲涅尔衍射的光场脉冲响应函数，使用三次快速傅里叶变换进行菲涅尔衍射计算，将当前传播平面的虚拟视场矩阵传播到下一传播平面中，得到新的虚拟视场矩阵；
A3、构建一个大小为2M*2M的掩膜矩阵，其中，中心直径为D的圆形部分的元素值为1，其余部分的元素值为0，将该掩膜矩阵与所述新的虚拟视场矩阵进行点乘得到一新矩阵，作为所述下一传播平面的虚拟视场矩阵；
A4、将所述下一传播平面作为当前传播平面返回步骤A2，以循环步骤A2和A3进行分段传播，直至传播到达目标平面，并从所述目标平面的中心提取出大小为M*M的第二全视场矩阵，即为最终的传播结果。
在一些实施例中，步骤A2中利用三次快速傅里叶变换计算菲涅尔衍射的条件为：从当前传播平面传播到下一传播平面时，在当前传播平面上需要满足式(1)：
Figure BDA0002798745980000021
其中，ω为光源的最大线性宽度，L为进行傅里叶变换的矩阵数组边长对应的长度。
在一些实施例中，步骤A2中设计的所述光场脉冲响应函数为：
Figure BDA0002798745980000022
其中，(x1,y1)表示源平面内的像素坐标，(x2,y2)表示目标平面内的像素坐标，e为自然常数，i为虚数，k为波数，z为相邻两个传播平面间的传播距离，λ为波长。
在一些实施例中，步骤A2中三次快速傅里叶变换计算菲涅尔衍射包括：第一次变换：对当前传播平面内的光场进行快速傅里叶变换；第二次变换：对菲涅尔衍射的光场脉冲响应函数进行快速傅里叶变换；第三次变换：对上述第一次变换和第二次变换的结果进行点乘得到新的结果，对所述新的结果进行快速反傅里叶变换，得到下一传播平面的光场。
在一些实施例中，步骤A3中所述中心直径D的设计方法为：
Figure BDA0002798745980000023
其中，NF为菲涅尔系数，λ为波长，z为相邻两个传播平面间的传播距离；同时，所述中心直径D还需满足D≤M。
本发明的有益效果在于：传统的分段传播模型由于使用了卷积运算，耗时较长，而本发明中使用三次快速傅里叶变换代替卷积运算来加速菲涅尔衍射的计算过程。另外，为了解决频谱混叠问题，本发明在使用三次快速傅里叶变换时不仅需满足光源的最大线性宽度小于或等于用于传播的矩阵边长一半这个条件，还需要光源大小、波长和传播距离之间满足菲涅耳系数，但由于快速傅里叶变换进行传播的模型仅能进行单次传播，在进行第二次传播时不满足上述条件，为此，本发明通过使用虚拟视场矩阵和掩膜矩阵将分段传播时传播面内光源的最大线性宽度控制在满足上述条件的范围内，从而解决了频谱混叠问题的同时保证了可以实现菲涅尔衍射分段传播的快速仿真，使得本发明可以应用于光路需要进行分段传播的复杂仿真实验中。
附图说明
图1是本发明实施例提出的菲涅尔衍射分段传播的快速仿真方法流程图；
图2是分段传播仿真模型示意图；
图3是用于防止频谱混叠的掩膜示意图。
具体实施方式
下面结合附图和具体的实施方式对本发明作进一步说明。
本发明使用三次快速傅里叶变换来加速菲涅尔衍射的计算，同时为解决菲涅尔衍射分段传播仿真中出现频谱混叠的问题，使用虚拟视场矩阵和掩膜将多次传播时传播面内光源的最大线性宽度进行合理控制，解决了频谱混叠问题。本发明实施例为了解决上述问题提出一种菲涅尔衍射分段传播的快速仿真方法，参考图1，该快速仿真方法包括如下步骤A1～A4：
A1、输入一个大小为M*M的包含光源信息的第一全视场矩阵，针对该第一全视场矩阵建立一个全零的大小为2M*2M的虚拟视场矩阵，并将输入的第一全视场矩阵放置在虚拟视场矩阵的中心。此处矩阵的行列大小M根据光源直径d(单位“米”)来设置，具体可以根据实际的光源定义矩阵内每一格所代表的物理尺寸，使M满足M≥2d即可。
A2、根据当前传播平面与下一传播平面间的距离以及平面内的像素大小设计菲涅尔衍射的光场脉冲响应函数，使用三次快速傅里叶变换进行菲涅尔衍射计算，将当前传播平面的虚拟视场矩阵传播到下一传播平面中，得到新的虚拟视场矩阵。其中，利用三次快速傅里叶变换计算菲涅尔衍射的条件为：从当前传播平面传播到下一传播平面时，在当前传播平面上需要满足式(1)：
Figure BDA0002798745980000041
其中，ω为光源的最大线性宽度，L为进行傅里叶变换的矩阵数组边长对应的长度。步骤A1中将大小为M*M的所述第一全视场矩阵放置在大小为2M*2M的虚拟视场矩阵的中心，可以使得输入的传播面满足式(1)的条件，从而可以使用快速傅里叶变换计算菲涅尔衍射，防止频谱混叠现象发生。
A3、构建一个大小为2M*2M的掩膜矩阵，该掩膜矩阵中，中心直径为D的圆形部分的元素值为1，其余部分的元素值为0，将该掩膜矩阵与所述新的虚拟视场矩阵进行点乘得到一新矩阵，作为所述下一传播平面的虚拟视场矩阵。在一些实施例中，中心直径D的设计方法为：
Figure BDA0002798745980000042
其中，NF为菲涅尔系数，λ为波长，z为相邻两个传播平面间的传播距离；同时，所述中心直径D还需满足式(1)，即利用掩膜后当光源的最大线性宽度ω＝D、进行傅里叶变换的矩阵数组边长L＝2M时，D满足
Figure BDA0002798745980000043
即D≤M。
叠加掩膜防止频谱混叠的示意图如图3所示，圆形部分表示叠加掩膜后保留的部分。使用掩膜的目的是通过生成上述掩膜并点乘新的虚拟视场矩阵，可以使得处理后的新矩阵满足式(1)和式(2)的要求。
A4、将所述下一传播平面作为当前传播平面返回步骤A2，以循环步骤A2和A3进行分段传播，直至传播到达目标平面，并从所述目标平面的中心提取出大小为M*M的第二全视场矩阵，即为最终的传播结果。
比如，参考图2，平面1表示源平面，平面n表示目标平面，则传播到目标平面n需要进行n-1次传播，即分n-1段传播，源平面即为第一次传播时的当前传播平面，平面2即为第一次传播时的下一传播平面。传播过程为：
第一次传播：输入大小为M*M的包含光源信息(光源O)的全视场矩阵MA，并针对全视场矩阵MA建立一个全零的大小为2M*2M的虚拟视场矩阵MB1，将输入的全视场矩阵MA放置在虚拟视场矩阵MB1的中心。根据平面1与平面2间的传播距离z1以及平面内的像素大小设计此时菲涅尔衍射的光场脉冲响应函数，使用三次快速傅里叶变换进行菲涅尔衍射计算，将平面1的虚拟视场矩阵MB1传播到平面2中，得到新的虚拟视场矩阵MB1’；判断是否到达目标平面，如果是，则直接从新的虚拟视场矩阵MB1’中心提取出一个M*M全视场矩阵作为最终传播结果；如果否，则：构建一个大小为2M*2M的掩膜矩阵MASK，将掩膜矩阵MASK与新的虚拟视场矩阵MB1’进行点乘得到一新矩阵MB2，P即为叠加掩膜后保留的部分，MB2作为平面2的虚拟视场矩阵，进入第二次传播；
第二次传播：平面2为当前传播平面，平面3为下一传播平面，根据平面2与平面3间的传播距离z2以及平面内的像素大小设计此时菲涅尔衍射的光场脉冲响应函数，使用三次快速傅里叶变换进行菲涅尔衍射计算，将平面2的虚拟视场矩阵MB2传播到平面3中，得到新的虚拟视场矩阵MB2’；同样地，判断是否到达目标平面，如果是，则直接从新的虚拟视场矩阵MB2’中心提取出一个M*M全视场矩阵作为最终传播结果；如果否，则：将掩膜矩阵MASK与新的虚拟视场矩阵MB2’进行点乘得到一新矩阵MB3，作为平面3的虚拟视场矩阵，继续进行传播；
以此类推实现分段传播。
在一些实施例中，步骤A2根据当前传播平面与下一传播平面间的传播距离、各传播平面内的像素大小，可以设计出菲涅尔衍射的光场脉冲响应函数，即：
Figure BDA0002798745980000051
其中，(x1,y1)表示源平面内的坐标，(x2,y2)表示目标平面内的坐标，e为自然常数，i为虚数，k为波数，z为相邻两个传播平面间的传播距离，λ为波长。
在一些实施例中，三次快速傅里叶变换计算菲涅尔衍射包括：对当前传播平面内的光场进行快速傅里叶变换；对菲涅尔衍射的光场脉冲响应函数进行快速傅里叶变换；将上述两次变换的结果进行点乘得到新的结果，对所述新的结果进行快速反傅里叶变换，得到下一传播平面的光场。具体为(以第一次传播为例)：
将当前传播平面1的光场E1(x1,y1)和光场的脉冲响应函数h(x2,y2)分别通过傅里叶变换到频域下并进行点乘，再通过傅里叶反变换得到空域下平面2的光场E2(x2,y2)，从而代替将平面1的光场E1(x1,y1)和光场的脉冲响应函数h(x2,y2)进行卷积运算得到平面2的光场E2(x2,y2)这一运算过程：
原先的卷积运算过程如式(4)
E2(x2,y2)＝E1(x1,y1)*h(x2,y2) (4)
式(4)中“*”表示卷积运算；
采用三次快速傅里叶变换替代上述卷积运算，即
E2(x2,y2)＝F-1(F(E1(x1,y1))·F(h(x2,y2))) (5)
其中，F()表示对括号内的对象进行傅里叶变换，F-1()表示对括号内的对象进行傅里叶逆变换。
以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的技术人员来说，在不脱离本发明构思的前提下，还可以做出若干等同替代或明显变型，而且性能或用途相同，都应当视为属于本发明的保护范围。

Claims (5)

1.一种菲涅尔衍射分段传播的快速仿真方法，其特征在于，包括如下步骤：
A1、输入一个大小为M*M的包含光源信息的第一全视场矩阵，针对所述第一全视场矩阵建立一个全零的大小为2M*2M的虚拟视场矩阵，将输入的所述第一全视场矩阵放置在所述虚拟视场矩阵的中心；
A2、根据当前传播平面与下一传播平面间的距离以及平面内的像素大小设计菲涅尔衍射的光场脉冲响应函数，使用三次快速傅里叶变换进行菲涅尔衍射计算，将当前传播平面的虚拟视场矩阵传播到下一传播平面中，得到新的虚拟视场矩阵；
A3、构建一个大小为2M*2M的掩膜矩阵，其中，中心直径为D的圆形部分的元素值为1，其余部分的元素值为0，将该掩膜矩阵与所述新的虚拟视场矩阵进行点乘得到一新矩阵，作为所述下一传播平面的虚拟视场矩阵；
A4、将所述下一传播平面作为当前传播平面返回步骤A2，以循环步骤A2和A3进行分段传播，直至传播到达目标平面，并从所述目标平面的中心提取出大小为M*M的第二全视场矩阵，即为最终的传播结果。
2.如权利要求1所述的菲涅尔衍射分段传播的快速仿真方法，其特征在于，步骤A2中利用三次快速傅里叶变换计算菲涅尔衍射的条件为：
从当前传播平面传播到下一传播平面时，在当前传播平面上需要满足式(1)：
Figure FDA0002798745970000011
其中，ω为光源的最大线性宽度，L为进行傅里叶变换的矩阵数组边长对应的长度。
3.如权利要求1所述的菲涅尔衍射分段传播的快速仿真方法，其特征在于，步骤A2中设计的所述光场脉冲响应函数为：
Figure FDA0002798745970000012
其中，(x1,y1)表示源平面内的像素坐标，(x2,y2)表示目标平面内的像素坐标，e为自然常数，i为虚数，k为波数，z为相邻两个传播平面间的传播距离，λ为波长。
4.如权利要求1所述的菲涅尔衍射分段传播的快速仿真方法，其特征在于，步骤A2中三次快速傅里叶变换计算菲涅尔衍射包括：
第一次变换：对当前传播平面内的光场进行快速傅里叶变换；
第二次变换：对菲涅尔衍射的光场脉冲响应函数进行快速傅里叶变换；
第三次变换：对上述第一次变换和第二次变换的结果进行点乘得到新的结果，对所述新的结果进行快速反傅里叶变换，得到下一传播平面的光场。
5.如权利要求1所述的菲涅尔衍射分段传播的快速仿真方法，其特征在于，步骤A3中所述中心直径D的设计方法为：
Figure FDA0002798745970000021
其中，NF为菲涅尔系数，λ为波长，z为相邻两个传播平面间的传播距离；
同时，所述中心直径D还需满足D≤M。
CN202011341503.7A 2020-11-25 2020-11-25 一种菲涅尔衍射分段传播的快速仿真方法 Active CN112380723B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011341503.7A CN112380723B (zh) 2020-11-25 2020-11-25 一种菲涅尔衍射分段传播的快速仿真方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011341503.7A CN112380723B (zh) 2020-11-25 2020-11-25 一种菲涅尔衍射分段传播的快速仿真方法

Publications (2)

Publication Number Publication Date
CN112380723A CN112380723A (zh) 2021-02-19
CN112380723B true CN112380723B (zh) 2022-09-16

Family

ID=74587809

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011341503.7A Active CN112380723B (zh) 2020-11-25 2020-11-25 一种菲涅尔衍射分段传播的快速仿真方法

Country Status (1)

Country Link
CN (1) CN112380723B (zh)

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10372845B1 (en) * 2015-10-30 2019-08-06 Nova Measuring Instruments Ltd. Scatterometry system and method
FR3080469A1 (fr) * 2018-04-23 2019-10-25 B<>Com Procede de traitement d'un hologramme, dispositif, systeme de restitution holographique et programme d'ordinateur associes

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10372845B1 (en) * 2015-10-30 2019-08-06 Nova Measuring Instruments Ltd. Scatterometry system and method
FR3080469A1 (fr) * 2018-04-23 2019-10-25 B<>Com Procede de traitement d'un hologramme, dispositif, systeme de restitution holographique et programme d'ordinateur associes

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于微焦斑源X射线传播的相衬成像模拟;刘鑫等;《深圳大学学报(理工版)》;20070731(第03期);全文 *

Also Published As

Publication number Publication date
CN112380723A (zh) 2021-02-19

Similar Documents

Publication Publication Date Title
KR102672004B1 (ko) 저 정밀도 뉴럴 네트워크 학습을 위한 방법 및 장치
AU2022200600B2 (en) Superpixel methods for convolutional neural networks
TWI638272B (zh) 用於對類神經網路執行類神經網路計算之系統與方法及相關正規化電路
CN112673383A (zh) 神经网络核中动态精度的数据表示
JP2020532777A (ja) ディープニューラルネットワークの実行方法、実行装置、学習方法、学習装置及びプログラム
US10757446B2 (en) Data compression apparatus and data compression method
US20200293857A1 (en) Cnn processing device, cnn processing method, and program
CN110782393A (zh) 一种基于可逆网络的图像分辨率压缩及重建方法
CN112380723B (zh) 一种菲涅尔衍射分段传播的快速仿真方法
CN111951171A (zh) Hdr图像生成方法、装置、可读存储介质及终端设备
Yang et al. Fractional‐order tensor regularisation for image inpainting
CN113222856A (zh) 一种逆半色调图像处理方法、终端设备及可读存储介质
CN110245706B (zh) 一种针对嵌入式应用的轻量化目标检测方法
US8094956B2 (en) Method and device for down-sampling a DCT image in the DCT domain
CN110084759B (zh) 一种图像填补方法、终端设备及存储介质
Puchala et al. Numerical accuracy of integral images computation algorithms
CN114549314A (zh) 一种提高图像分辨率的方法
CN112949841A (zh) 一种基于Attention的CNN神经网络的训练方法
CN115861099B (zh) 一种引入物理成像先验知识约束的卫星云图图像复原方法
Hunt Minimizing the computation time for using the technique of sectioning for digital filtering of pictures
CN105976317A (zh) 一种图像空间退化模拟方法及系统
Chen et al. Locally edge-adapted distance for image interpolation based on genetic fuzzy system
US8232993B2 (en) Systems and methods for graphical rendering
CN116503252A (zh) 生成图像超分数据集的方法、图像超分模型及训练方法
US20050213851A1 (en) Scaling device and method for scaling a digital picture

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant