0

这是我计算菲涅耳衍射的 MATLAB 代码。我想在 C++ 中实现相同的代码。

clc;clear all;close all;

N=512;
M=512;
lambda=632e-9;
X = 12.1e-6;
k=2*pi/lambda;
z0 = (N*X^2)/lambda;
a = (exp(j*k*z0)/(j*lambda))
hX=(exp(j*k*z0)/(j*lambda))*exp(j*pi/(z0*lambda)*X^2*[0:N-1].^2);
hY=(exp(j*k*z0)/(j*lambda))*exp(j*pi/(z0*lambda)*X^2*[0:M-1].^2);
h = hX.'*hY;
figure; imshow(real(h), []);

当我试图在 C++ 中实现相同的代码时,我尝试了这个:

int main (){
std::complex<double> com_two(0,1);
double mycomplex = 1;
double pi = 3.1415926535897;
double N = 512;
double M = 512;
double lambda = 632e-9;
double X = 12.1e-6;
double k = (2*pi)/lambda;
double z0=(N*pow(X,2))/lambda;
//std::complex<double> hx; // Both definitions don't work
//double hy;
/*for (int j=0; j < N; j++){
    hy = (exp(mycomplex*k*z0) / (mycomplex*lambda))*exp(mycomplex*pi/(z0*lambda)*pow(X,2)*pow(j,2));
    cout << hy <<endl;
}*/
system("pause");
return 0;

}

但问题是,在执行 hx 和 hy 值的计算时会返回复数。

那么我应该如何在 MATLAB 代码中定义像 i,j 这样的“mycomplex”。我还发现了一些作业,例如std::complex<double> complex;

但我想这行不通。在我看来,我必须将 hx 和 hy 值作为复数。但我找不到合适的解决方案。

如果有人对此提供帮助,我将非常高兴。

4

2 回答 2

1

问题不是复数本身,而是您执行的计算。首先,您在 Matlab alg(至少您展示的版本)和 C++ 中做不同的事情。这是您应该如何将此代码移植到 C++ 的方法。

#include <cstdlib>
#include <iostream>
#include <complex>
using namespace std;
/*
 * 
 */
int main(int argc, char** argv) {
    double pi = 3.1415926535897;
    double N = 512;
    double M = 512;
    double lambda = 632e-9;
    double X = 12.1e-6;
    double k = (2 * pi) / lambda;
    double z0 = (N * pow(X, 2)) / lambda;
    std::complex<double> hx, hy; // Both definitions don't work
    //double hy;
    int j = 1;
    //for (int j=0; j < N; j++){
        hy = (exp(j*k*z0) / (j*lambda))*exp(j*pi/(z0*lambda)*pow(X,2)*pow(j,2));
        cout << hy <<endl;
    //}
    system("pause");
    return 0;
}

现在问题出在你lambda = 632e-9的接近于零,然后 k 非常大k = (2 * pi) / lambda;,因为它是接近 0 的倒数,反过来

exp(j*k*z0)

很大,而且这个

(exp(j*k*z0) / (j*lambda))

更。那么这个

(exp(j*k*z0) / (j*lambda))*exp(j*pi/(z0*lambda) 

是巨大的,剩下的因素pow(X,2)*pow(j,2))是有意义的。那么那么

cout << hy <<endl; 

(inf,0)即使 lambda = 0.0001也会打印;(632e-9 怎么样!)。

但是,如果您使用足够小的j值来吸收kexp 中的巨大值,您仍然可以打印它(j也在分母中,但de l'Hopital 规则将确保这会更小)。

从您的评论(不是您的代码)中可以看出,您希望j成为复数 z=1i。那么解决方法就是改变

double j =1;

std::complex<double> j(0,1); 

解决方案:

int main(int argc, char** argv) {
    double pi = 3.1415926535897;
    double N = 512;
    double M = 512;
    double lambda = 0.00001;//632e-9;
    double X = 12.1e-6;
    double k = (2 * pi) / lambda;
    double z0 = (N * pow(X, 2)) / lambda;
    std::complex<double> hx, hy;
    std::complex<double> j(0,1);
        hy = (exp(j*k*z0) / (j*lambda))*exp(j*pi/(z0*lambda)*pow(X,2)*pow(j,2));
        cout << hy <<endl;
    return 0;
}
于 2013-11-01T16:00:13.230 回答
0

我在您的代码中更改了几行,它给出的结果与 1#INF 不同

1)mycomplexdoubletostd::complex<double>和 assign(0,1)使其变为虚构 1. 实际上com_two具有您需要的相同类型和值mycomplex并且您没有使用它。

2) 删除评论。hy将类型从更改doublestd::complex<double>。删除未使用hx的声明。

pow(j,2)3)由于歧义,我的编译器发现错误,所以我也修复了它。

#include <iostream>
#include <complex>

using namespace std;

int main (){
  complex<double> mycomplex(0,1);
  double pi = 3.1415926535897;
  double N = 512;
  double M = 512;
  double lambda = 632e-9;
  double X = 12.1e-6;
  double k = (2*pi)/lambda;
  double z0=(N*pow(X,2))/lambda;
  complex<double> hy;
  for (int j=0; j < N; j++){
      hy = (exp(mycomplex*k*z0) / (mycomplex*lambda))*exp(mycomplex*pi/(z0*lambda)*pow(X,2)*pow((double)j,2));
      cout << hy <<endl;
  }
  system("pause");
  return 0;
}
于 2013-11-01T15:54:23.427 回答