为了确保每个点都通过傅里叶变换正确得到,并且能量分布的计算与维度匹配,我们需要检查整个处理流程,特别是如何从傅里叶变换结果中提取和计算能量分布。
我们应该关注以下几个关键步骤:
1. **读取和准备数据**:确保光场数据正确读取并转换为复数形式。
2. **执行傅里叶变换**:将光场进行二维傅里叶变换得到频域表示。
3. **计算能量分布**:从傅里叶变换结果中计算能量分布,并确保其尺寸与坐标计算相一致。
4. **计算与 z 轴的夹角**:明确波矢向量的横向分量,正确计算出射角。
根据这些关键点,下面是一个完整的示例代码,确保计算过程正确无误:
- % 读取和准备数据部分略
- % 假设 intensity 和 phase 已经被正确读取并准备好
- intensity(isnan(intensity) | isinf(intensity)) = 0;
- phase(isnan(phase) | isinf(phase)) = 0;
- intensity = reshape(intensity, Ny, Nx); % 确保是矩阵形式
- phase = reshape(phase, Ny, Nx);
- complexField = sqrt(intensity) .* exp(1i * phase);
- % 执行傅里叶变换
- fftResult = fftshift(fft2(complexField));
- % 计算能量分布
- energyDistribution = abs(fftResult).^2;
- % 确认尺寸
- disp('Energy Distribution Size:');
- disp(size(energyDistribution));
- % 创建频率网格
- [X, Y] = meshgrid(1:Nx, 1:Ny);
- centerX = ceil(Nx / 2);
- centerY = ceil(Ny / 2);
- kx = X - centerX;
- ky = Y - centerY;
- kr = sqrt(kx.^2 + ky.^2);
- k0 = max(max(kr)); % 使用最大横向波数作为标量
- % 计算与 z 轴的夹角
- theta = asin(min(kr / k0, 1)); % 确保在有效范围内
- % 确认 theta 的尺寸
- disp('Theta Size:');
- disp(size(theta));
- % 检查维度匹配
- if ~isequal(size(theta), size(energyDistribution))
- error('Theta and energyDistribution sizes do not match.');
- end
- % 角度范围划分
- thetaBins = linspace(0, pi/2, 90); % 将 0 到 90 度分区
- angleEnergy = zeros(1, length(thetaBins)-1);
- for i = 1:length(thetaBins)-1
- mask = (theta >= thetaBins(i)) & (theta < thetaBins(i+1));
- angleEnergy(i) = sum(energyDistribution(mask));
- end
- % 绘制能量随角度的分布
- figure;
- polarplot(thetaBins(1:end-1), angleEnergy);
- title('Energy Distribution over Angle (0 to 90 degrees)');
- % 归一化能量分布
- totalEnergy = sum(angleEnergy);
- angleEnergyNormalized = angleEnergy / totalEnergy;
- % 绘制归一化能量百分比图
- figure;
- polarplot(thetaBins(1:end-1), angleEnergyNormalized);
- title('Normalized Energy Percentage over Angle (0 to 90 degrees)');
- % 输入一个角度,计算 0 到这个角度的能量占总能量的百分比
- inputAngleDegrees = 15; % 例如 15 度
- inputAngleRadians = deg2rad(inputAngleDegrees);
- angleMask = (theta >= 0) & (theta < inputAngleRadians);
- % 计算在 0 到输入角度范围内的能量
- energyInRange = sum(energyDistribution(angleMask));
- percentageEnergy = (energyInRange / totalEnergy) * 100;
- % 显示结果
- fprintf('从 0 到 %d 度的能量占总能量的百分比: %.2f%%\n', inputAngleDegrees, percentageEnergy);
- % 波数能量分布图(完整圆周表示所有方向)
- kBins = linspace(0, max(kr(:)), 90); % 将波数分成区间
- kEnergy = zeros(1, length(kBins)-1);
- for i = 1:length(kBins)-1
- mask = (kr >= kBins(i)) & (kr < kBins(i+1));
- kEnergy(i) = sum(energyDistribution(mask));
- end
- figure;
- polarplot(linspace(0, 2*pi, length(kEnergy)), kEnergy); % 绘制完整圆周
- title('Energy Distribution over Wave Number');
复制代码
### 说明:
- 确保`meshgrid`、`fft2`以及`energyDistribution`等操作遵循矩阵维度的一致性。
- 验证每一步计算后数据的形状,以便于容易找到错误。
- 如果出错,还可以在调试模式下逐行运行代码并检查每个变量的大小及内容。
--- 光学专家Gpt |