我正在编写一个函数以从给定的时间序列中获取随机时间序列$\frac{1}{f}$法律。通过在函数中引入随机相位来获得随机化。

我在模拟时间序列中遇到问题，值系统地上升。我想这与函数本身的定义有关。

编写函数报告如下：

def Pk(k, alpha=-11.0 / 3, fknee=1):
    """Simple power law formula with ASD"""
    return (k / fknee) ** alpha


def gen_pkfield(nsample=100, alpha=-11.0 / 3, fknee=1, fsample=2.5):
    
    ufreq = np.fft.fftfreq(n=nsample+1, d=1./fsample)   # the Discrete Fourier Transform sample frequencies (with negative values)
    kfreq = np.sqrt(ufreq ** 2) # to create the positive one
    

    with np.errstate(divide="ignore"):
        asd_per = Pk(kfreq, alpha=alpha, fknee=fknee)   # there is a factor sqrt(2) for the ASD


    pha = np.random.uniform(low=-np.pi, high=np.pi, size=(nsample+1) )
    
    asd_random = asd_per * (np.cos(pha) + 1j * np.sin(pha))     # frequency domain signal with random phase
        
    psd_per = asd_random[1:]**2.    # psd periodogram with random phase 
    psd_norm =  psd_per  * (fsample*nsample/2.) 
    asd_norm = np.sqrt(psd_norm ) 
    
    plt.figure()
    plt.plot(asd_norm)
    plt.ylabel(r'Amplitude spectral density [Hz/$\sqrt{Hz}$] ')
    plt.xlabel('Frequency [Hz]')
    plt.tight_layout()  
    
    ifft = np.real( np.fft.ifft( asd_norm )  )     # reverse the normalization for ifft coherence   
    
    return ifft 

我介绍的各种规范化是为了numpy.ifft与scipy.periodogram（我从中获取$\alpha$和$knee$值）。

这里获得的 ASD（asd_norm在函数中命名）。这是通过$\alpha=0.9$和$knee=0.148$

这里是模拟时间序列的结果：

如您所见，时间序列值在边界处上升。我想知道为什么。

新的

我更正了函数内部的内容。我混合了PSD和ASD。这里是新定义：

def gen_pkfield(nsample=100, alpha=-11.0 / 3, fknee=1, fsample=2.5): # sampling per block
    
    ufreq = np.fft.fftfreq(n=nsample, d=1./fsample)     # the Discrete Fourier Transform sample frequencies (with negative values)
    kfreq = np.sqrt(ufreq[:] ** 2 + ufreq ** 2)
    with np.errstate(divide="ignore"):
        psd_per = Pk(kfreq, alpha=alpha, fknee=fknee)   # there is a factor sqrt(2) for the ASD
            
    asd_per = psd_per * np.sqrt( fsample*nsample/2. )   # psd periodogram with random phase
    
            
    pha = np.random.uniform(low=-np.pi, high=np.pi, size=(nsample) )
    asd_norm = asd_per * (np.cos(pha) + 1j * np.sin(pha))   # frequency domain signal with random phase

    
    plt.figure()
    plt.plot(kfreq,asd_norm,'k-')
    plt.ylabel(r'Amplitude spectral density [Hz/$\sqrt{Hz}$] ')
    plt.xlabel('Frequency [Hz]')
    plt.tight_layout()  
    
    ifft = np.real( np.fft.ifft( asd_norm )  )     # reverse the normalization for ifft coherence   
    
    
    return ifft 

现在的主要问题是我获得了公正的inf价值观。

===============

解决了

通过删除直流分量（频率 = 0 Hz）

def gen_pkfield(nsample=100, alpha=-11.0 / 3, fknee=1, fsample=2.5): # sampling per block
    
    ufreq = np.fft.fftfreq(n=nsample+1, d=1./fsample)[1:]   # the Discrete Fourier Transform sample frequencies (with negative values), [1:] to solve the DC
    kfreq = np.sqrt(ufreq[:] ** 2 + ufreq ** 2)
        
    with np.errstate(divide="ignore"):
        psd_per = Pk(kfreq, alpha=alpha, fknee=fknee)   # there is a factor sqrt(2) for the ASD
            
    asd_per = np.sqrt( psd_per * fsample*nsample/2. )   # psd periodogram with random phase
            
    pha = np.random.uniform(low=-np.pi, high=np.pi, size=(nsample) )
    asd_norm = asd_per * (np.cos(pha) + 1j * np.sin(pha))   # frequency domain signal with random phase 

    ifft = np.real( np.fft.ifft( asd_norm )  )     # reverse the normalization for ifft coherence   
    
    return ifft