光电容积脉搏波法（PPG）是一种低成本，无创的光学技术，可在皮肤表面进行生理测量。其极广泛的应用之一是商用智能手表和运动手环中包含的可穿戴心率传感器，它在日常环境下可提供舒适和连续的脉搏监测。本文演示了如何在 Zemax OpticStudio 中对人体皮肤建模以进行生理测量，并说明了使用 ZOS-API 对基于 PPG 的心率传感器进行的时间相关模拟。
PPG 器件由红外或可见光波长范围内的发光二极管 (LED) 和光电探测器组成。它们提供了一种简单的光学技术来检测组织中的血容量变化，因为血液比周围的组织对光具有更强烈地吸收和散射效应。因此，血液的脉动将导致检测器信号发生相反的相位变化。本文介绍如何在 OpticStudio 中模拟人体皮肤组织模型，并演示如何使用 ZOS-API 应用程序模拟 PPG 设备随时间推移的测量信号。
PPG 传感器可设计为反射或透射模式。由于光的穿透深度取决于其波长，因此绿色和黄色 LED 光线极适合在浅表血流中进行测量，并且通常以反射模式使用。另一方面，红外和近红外波长更适合测量深层组织血流，可用于透射模式。在次案例中，我们展示了一个反射 PPG 设备。
我们的目标是根据相关文献中发表的数据开发一个逼真的皮肤模型。因此，我们打算应用某种波长，通常设置为对应皮肤和血液的光学参数在文献中广泛可见的波长，并且也接近商业设备中极常用的波长。因此，我们建模选择了 575nm 的波长，并使用 QSMF-C160 LED (Avago Technologies) 作为光源。此 LED 的模型可以直接从 Radiant Source Model 数据库下载，并且可以通过从 Radiant Source Model 文件生成的光线来创建光源文件。此功能在 OpticStudio 的旗舰版中可用，知识库文章 如何使用 RSMX 光源模型生成光线集（https://support.zemax.com/hc/zh-cn/articles/1500005487321） 中详细讨论了该过程。
为了模拟人体组织介质中的光传输，我们创建了分层皮肤模型，该模型考虑了表皮、真皮和皮下脂肪的组织结构。由于此例的主要目标是模拟基于 PPG 的心率传感器，其中关键点是测量由血液脉动引起的变化，因此我们专注于准确地建模可以观察到这种脉动的皮肤层。因此，我们分别对不同血含量值的真皮亚层进行建模，即真皮乳头层、毛细血管下皮层、上层血网真皮层、网状真皮层和深部血网真皮层。另一方面，由于表皮中没有血液含量，为了保持模型简单我们只使用一个厚表皮层，它包含所有角质层、颗粒层、棘层和基底层。然后，与大多数已发布的皮肤模型类似，我们也用一层结构表示皮下脂肪。
上述所有皮肤层在 OpticStudio 中都建模为矩形物体，每层的厚度值基于文献数据，通过使得侧面没有漏光的方式来确定横截面的尺寸。通过使用上一个层作为参考对象并对上一个层的 Z 长度单元格中的 Z 位置值应用 Pickup 求解来放置后续层结构。该解决方案确保了各层结构紧贴在一起，之间没有任何间隙。
自定义组织层结构
由于本案例研究仅依赖于文献中发表的数据，因此我们在整个建模过程中并未进行任何新的测量。尽管模型参数基于已发布的数据，但需要注意的是，人体皮肤的光学参数在不同人群中可能存在明显差异。因此，特定主题可能需要使用不同的参数。所以，如果您的特定应用可以获得更准确的数据，请制定相应的组织结构模型。
详细表示皮肤中的所有每条血管将需要添加数百个具有复杂空间排列的物体，并且会降低模型的通用性，因此这种类型的建模在文献中并不进行考虑，所以我们也没有应用这种建模方式。相反，我们通过计算血液和周围组织结构的光学参数的加权平均值来考虑不同皮肤层的血液含量。
因此，我们基于以下原始数据，使用 OpticStudio 中的模型材料求解对表皮层的材料进行了建模：
皮肤层结构 | 575 nm 波长 折射率 | 血液含量 [%] |
活表皮 | 1.45 | 0 |
真皮层 | 1.40 | 2 |
上层血网真皮层 | 1.40 | 5 |
网状真皮层 | 1.40 | 1 |
深层血网真皮层 | 1.40 | 5 |
皮下脂肪 | 1.44 | 5 |
血液 | 1.35 | n/a |
人体组织中的体散射
Henyey-Greenstein 分布函数可以准确描述混浊介质（例如生物组织）中小颗粒的光线散射。Henyey-Greenstein 模型只有一个自由参数，即各向异性因子 g。该参数的域区间为 [-1, 1]，其中 g=-1 对应反向散射，g=0 表示各向同性散射，g=1 表示正向散射。散射光的角度分布定义为：
在 OpticStudio 的非序列模式下，Henyey-Greenstein 体散射模型以 DLL (Henyey-Greenstein-bulk.DLL) 的形式提供，包含在 OpticStudio 安装文件中。知识库文章 使用 Henyey-Greenstein 分布对体散射进行建模（https://support.zemax.com/hc/zh-cn/articles/1500005577062） 对模型进行了进一步的解释以及对 DLL 进行了全面分析。此外，有关 OpticStudio 中表面和体散射模型的详细讨论可以在此知识库文章中找到：OpticStudio 中有哪些可用散射模型？（https://support.zemax.com/hc/zh-cn/articles/1500005575662）
在文章附件的多层皮肤模型中，每一层的散射参数都是根据文献中给出的实际值设置的。虽然 Henyey-Greenstein 散射 DLL 的输入参数是平均路径、传输分数和各向异性参数 g，但在文献中通常散射和吸收系数分别定义为 μs 与 μa，与各向异性因子一起展示。因此，我们使用以下公式来计算模型的输入参数：
类似地，与之前的折射率、散射和吸收系数以及不同皮肤层的各向异性因子一样，为血液和组织其余部分相应值的加权平均值进行计算，以下为对应于 575 nm 应用波长的原始数据：
皮肤层结构 | μa [1/mm] | μs [1/mm] | g [-] |
活皮层 | 1 | 20 | 0.79 |
真皮层 | 0.28 | 21.5 | 0.79 |
上层血网真皮层 | 0.28 | 21.5 | 0.79 |
网状真皮层 | 0.28 | 21.5 | 0.79 |
深层血网真皮层 | 0.28 | 21.5 | 0.79 |
皮下脂肪 | 0.081 | 1.396 | 0.75 |
血液 | 32.4 | 50 | 0.98 |
多层皮肤模型和组织中的光线传输显示在下面的 3D 布局图中。为了说明每个皮层的散射情况，光线在图上按光线分段着色。
为了提供除显示图之外的数值结果，我们在设计中添加了三个矩形探测器。它们与皮肤表面之间被一个薄薄的空气间隙隔开（出于非序列模式建模嵌套规则的考虑）。两个探测器具有与皮肤层相同的截面尺寸，一个朝向光源，另一个朝向皮肤模型，以分别测量所有入射光和背向散射光供后续参考。第三个探测器是一个同样面向皮肤的小探测器 (2mm x 2mm)，代表了 PPG 设备中的典型光电探测器。
上述设计可以在文章附件 (skinModel.zar) 中找到，当测量/模拟与时间依赖性不相关时，它可以用作现成的皮肤模型。另一方面，对时间相关效应建模的方法——例如，在心率传感器的应用情况下——将在下一节中讨论。
为了模拟心率监测，我们可以使用 ZOS-API 来模拟组织中的脉动血流。我们通过调整皮肤层的血液含量来模拟心脏周期的不同阶段，然后我们检查探测到的背散射光作为时间步长的函数，将各层血液含量的时间变化通过倍增因子考虑在内，假设血容量在每一层中按比例且同时变化。在这个例子中，我们使用 Python API（通过 .NET 连接到 OpticStudio）来修改模型参数，使用微调设置运行光线追迹，然后分析和绘制结果。应用到的 Python 程序脚本可以从文章附件 (PPGsimulation.py) 下载。
通过 API 修改人体组织参数
根据文献，当心脏在收缩期将血液泵入血管时，皮肤层的相对血液含量会翻倍。我们使用经验函数来表征这种脉动，下图显示了 5 个心脏周期的结果（10步长/周期）。
基于此，我们首先计算了每一层的血液含量，然后我们相应地更新了折射率、平均自由程、透射率和Henyey-Greenstein 散射分布的各向异性因子 g。对应的Python代码如下：
layer = TheNCE.GetObjectAt(layerNum)
solver = layer.MaterialCell.CreateSolveType(ZOSAPI.Editors.SolveType.MaterialModel)
solver._S_MaterialModel.IndexNd = n
layer.MaterialCell.SetSolveData(solver)
volPhysData = layer.VolumePhysicsData
volPhysData.ModelSettings._S_DLLDefinedScattering.MeanPath = meanPath
volPhysData.ModelSettings._S_DLLDefinedScattering.SetParameterValue(0, transmission)
volPhysData.ModelSettings._S_DLLDefinedScattering.SetParameterValue(1, g)
结尾处，我们对每个时间步长运行光线追迹，并从对应 PPG 的小尺寸探测器中提取结果，即光线命中总数和吸收功率。为了提高光线追迹速度，我们应用了简单光线分裂 (Simple Ray Splitting)，其细节在知识库文章 什么是简单光线分裂？（https://support.zemax.com/hc/zh-cn/articles/1500005575742）中进行了讨论。光线追迹和数据提取可以使用以下 Python 命令完成：
NSCRayTrace = TheSystem.Tools.OpenNSCRayTrace()
NSCRayTrace.ClearDetectors(0)
NSCRayTrace.SplitNSCRays = True
NSCRayTrace.ScatterNSCRays = True
NSCRayTrace.UsePolarization = True
NSCRayTrace.IgnoreErrors = True
NSCRayTrace.SaveRays = False
NSCRayTrace.Run()
NSCRayTrace.WaitForCompletion()
NSCRayTrace.Close()
hits = -3 # pixel =-3 for total hits
power = 0 # pixel = 0 for total power
Data = 0
hits_bool_return, total_hits = TheNCE.GetDetectorData(detectorNum, hits, Data, 0)
power_bool_return, total_power = TheNCE.GetDetectorData(detectorNum, power, Data, 0)
由于血液的吸收和散射系数比组织的其余部分大得多，因此血液含量越高，背散射光能量越低，检测器的测量功率也越低，反之亦然。根据我们的模拟，在每个时间步长使用 1W 的光源和 105 条分析光线，建模的血含量脉动导致检测器信号发生 10~15% 的变化，如下图所示。
这些结果与文献中的实验观察结果非常吻合。在商用心率传感器中，作为后处理的一部分，通常会对该数据应用简单的信号处理算法。第 一步是对噪声信号进行平滑处理，然后计算每分钟超过某个预定义阈值的峰值数来计算心率。
参考文献
1. Maeda, et al. Monte Carlo simulation of spectral reflectance using a multilayered skin tissue model. Optical Review (2010)