Fortran 计算落后交叉相关系数

real function cross_R(x, y, n, t) ! 2个数组 n是个数 t是滞后系数
    integer :: n, t, i
    real, dimension(n) :: x, y
    real :: SD_x, SD_y, x_average, y_average
    x_average = sum(x) / n
    y_average = sum(y) / n
    cross_R = 0
    SD_x = sqrt(sum(x**2) / n - x_average**2)
    SD_y = sqrt(sum(y**2) / n - y_average**2)
    do i = 1, n - t
        cross_R = cross_R + (x(i) - x_average) * (y(i + t) - y_average) / (SD_x * SD_y)
    end do
    cross_R = cross_R / (n - t)
end function cross_R