C++(ROOT)数据处理之实现傅里叶变化及滤波分析

SIPM

1. 背景介绍
   前放低温下出现振荡低频噪声，致使信号分析出现偏差。如积分电荷或计算过阈脉冲个数时，振荡噪声使得p.e电荷谱前后带有小尾巴。使用傅里叶变化滤掉低频噪声可简化分析难度，提升结果准确性。
   傅里叶变换无须过多介绍，一言概之： 时空转为频谱。将频谱中不需要的频率分量去掉再反傅里叶变化即可得到所需的信号。详细推导可见kkk

2. 实现方法
   2.1 傅里叶变换：
   使用ROOT（基于C++）实现，主要代码如下：
   ```
   //A function to sample
   TF1 *fsin = new TF1("fsin", "sin(10*x)+sin(5*x)+sin(2*x)", 0, 2*TMath::Pi());
   fsin->Draw();

   Int_t n=250;
   TH1D *hsin = new TH1D("hsin", "hsin", n+1, 0, 2*TMath::Pi());
   Double_t x;

   //Fill the histogram with function values
   for (Int_t i=0; i<=n; i++){
   x = (Double_t(i)/n)*(2*TMath::Pi());
   hsin->SetBinContent(i+1, fsin->Eval(x));
   }
   //构建一个sin(10*x)+sin(5*x)+sin(2*x)的信号数据，包含三个频率
   //-----------------------------------------------------------------
   TH1 *hm =0;
   TVirtualFFT::SetTransform(0);
   hm = hsin->FFT(hm, "MAG");
   TH1 *hp = 0;
   hp = hsin->FFT(hp, "PH");
   //FFT即为快速傅里叶变换，变换后的mangnitude（频率相关）存为hm，相位相关信息存为hp
   //-----------------------------------------------------------

需注意对于hm：幅度分量（纵轴）需归一后方为真实值。对于直流分量，归一常数为1/N，N为直方图Bin数。对于交流分量，归一常数为2/N。
频率分量（横轴），考虑到上下限为0~$2*\pi$，其横坐标f真实值即为$2*\pi/f$,恰好对应10,5,2
2.2 滤波

Double_t *re_full = new Double_t[n];
Double_t *im_full = new Double_t[n];
fft->GetPointsComplex(re_full,im_full);
for(int i=3; i< 251; i++){
re_full[i] = 0;
// cout<<"-----"< im_full[i] = 0;
}

2.3 傅里叶逆变换

//Now let's make a backward transform:
TVirtualFFT *fft_back = TVirtualFFT::FFT(1, &n, "C2R M K");
fft_back->SetPointsComplex(re_full,im_full);
fft_back->Transform();
TH1 *hb = 0;
//Let's look at the output
hb = TH1::TransformHisto(fft_back,hb,"Re");
hb->SetTitle("The backward transform result");
hb->Draw();
//NOTE: here you get at the x-axes number of bins and not real values
//(in this case 25 bins has to be rescaled to a range between 0 and 4*Pi;
//also here the y-axes has to be rescaled (factor 1/bins)
``
3. 效果
4. summary