obspy.signal.cross_correlation.correlate

    """
    Cross-correlation of two signals up to a specified maximal shift.

    This function only allows 'naive' normalization with the overall
    standard deviations. This is a reasonable approximation for signals of
    similar length and a relatively small shift parameter
    (e.g. noise cross-correlation).
    If you are interested in the full cross-correlation function better use
    :func:`~obspy.signal.cross_correlation.correlate_template` which also
    provides correct normalization.

    :type a: :class:`~numpy.ndarray`, :class:`~obspy.core.trace.Trace`
    :param a: first signal
    :type b: :class:`~numpy.ndarray`, :class:`~obspy.core.trace.Trace`
    :param b: second signal to correlate with first signal
    :param int shift: Number of samples to shift for cross correlation.
        The cross-correlation will consist of ``2*shift+1`` or
        ``2*shift`` samples. The sample with zero shift will be in the middle.
    :param bool demean: Demean data beforehand.
    :param normalize: Method for normalization of cross-correlation.
        One of ``'naive'`` or ``None``
        (``True`` and ``False`` are supported for backwards compatibility).
        ``'naive'`` normalizes by the overall standard deviation.
        ``None`` does not normalize.
    :param str method: Method to use to calculate the correlation.
         ``'direct'``: The correlation is determined directly from sums,
         the definition of correlation.
         ``'fft'`` The Fast Fourier Transform is used to perform the
         correlation more quickly.
         ``'auto'`` Automatically chooses direct or Fourier method based on an
         estimate of which is faster. (Only availlable for SciPy versions >=
         0.19. For older Scipy version method defaults to ``'fft'``.)

    :return: cross-correlation function.

    To calculate shift and value of the maximum of the returned
    cross-correlation function use
    :func:`~obspy.signal.cross_correlation.xcorr_max`.

    .. note::

        For most input parameters cross-correlation using the FFT is much
        faster.
        Only for small values of ``shift`` (approximately less than 100)
        direct time domain cross-correlation migth save some time.

    .. note::

        If the signals have different length, they will be aligned around
        their middle. The sample with zero shift in the cross-correlation
        function corresponds to this correlation:

        ::

            --aaaa--
            bbbbbbbb

        For odd ``len(a)-len(b)`` the cross-correlation function will
        consist of only ``2*shift`` samples because a shift of 0
        corresponds to the middle between two samples.
    """

###helper

def _pad_zeros(a, num, num2=None):
    """Pad num zeros at both sides of array a"""
    if num2 is None:
        num2 = num
    hstack = [np.zeros(num, dtype=a.dtype), a, np.zeros(num2, dtype=a.dtype)]
    return np.hstack(hstack)

##################################################################################

def _xcorr_padzeros(a, b, shift, method):
    """
    Cross-correlation using SciPy with mode='valid' and precedent zero padding.
    """
    if shift is None:
        shift = (len(a) + len(b) - 1) // 2
    dif = len(a) - len(b) - 2 * shift
    if dif > 0:
        b = _pad_zeros(b, dif // 2)
    else:
        a = _pad_zeros(a, -dif // 2)
    return scipy.signal.correlate(a, b, mode='valid', method=method)
    
##################################################################################

def _xcorr_slice(a, b, shift, method):
    """
    Cross-correlation using SciPy with mode='full' and subsequent slicing.
    """
    mid = (len(a) + len(b) - 1) // 2
    if shift is None:
        shift = mid
    if shift > mid:
        # Such a large shift is not possible without zero padding
        return _xcorr_padzeros(a, b, shift, method)
    cc = scipy.signal.correlate(a, b, mode='full', method=method)
    return cc[mid - shift:mid + shift + len(cc) % 2]
    
##################################################################################
##################################################################################

### main

def correlate(a, b, shift, demean=True, normalize='naive', method='auto'):
  
    if normalize is False:
        normalize = None
    if normalize is True:
        normalize = 'naive'
    # if we get Trace objects, use their data arrays
    if isinstance(a, Trace):
        a = a.data
    if isinstance(b, Trace):
        b = b.data
    a = np.asarray(a)
    b = np.asarray(b)
    if demean:
        a = a - np.mean(a)
        b = b - np.mean(b)
        
    # choose the usually faster xcorr function for each method
    
    _xcorr = _xcorr_padzeros if method == 'direct' else _xcorr_slice
    cc = _xcorr(a, b, shift, method)
    
    if normalize == 'naive':
        norm = (np.sum(a ** 2) * np.sum(b ** 2)) ** 0.5
        if norm <= np.finfo(float).eps:
            # norm is zero
            # => cross-correlation function will have only zeros
            cc[:] = 0
        elif cc.dtype == float:
            cc /= norm
        else:
            cc = cc / norm
    elif normalize is not None:
        raise ValueError("normalize has to be one of (None, 'naive'))")
    return cc

 

obpsy correlation은 scipy이 correlate 함수를 쓴다.

 

참고 

2026.09.09 - [Python/Scipy] - scipy.signal.correlate

scipy.signal.correlate

 

obspy correlate의 작동방식을 보면,

 

먼저  전체적인 파이프라인은 아래와 같음

 

 

1. 입력 a,b
2. Trace면 .data 추출 및  demean
3. method에 따라 계산 방식 선택 (direct, fft, auto)
4. 원하는 ±shift 구간의 cross-correlation 계산 및 normalization
5. return cc

method	obspy 함수	scipy mode
'direct'	_xcorr_padzeros()	'valid'
'fft'	_xcorr_slice()	'full'
'auto'	_xcorr_slice()	'full'

 

_xcorr_slice()

full cross correlation을 계산하고 slice(사용자 지정 shift 만큼)

def _xcorr_slice(a, b, shift, method):

    mid = (len(a) + len(b) - 1) // 2

    if shift is None:
        shift = mid

    if shift > mid:
        return _xcorr_padzeros(a, b, shift, method)

    cc = scipy.signal.correlate(
        a,
        b,
        mode='full',
        method=method
    )

    return cc[
        mid - shift:
        mid + shift + len(cc) % 2
    ]

 

scipy.correlate을 이용해서 full correlation(lag에 대한 correlation 계산)을 계산한 후 그 간운데 ±shift만 잘라냄.

예시,

len(a) = 1001
len(b) = 1001

# full correlation length : 1001+1001−1=2001


mid = 2001 // 2
    = 1000
    
# full CC index 

index:
0 -------------------------------- 1000 -------------------------------- 2000

lag:
-1000 ----------------------------- 0 -------------------------------- +1000
                                      ↑
                                     mid
                                     
shift = 100 이라고 지정해주면

전체 full correlation

-1000 ================================================ +1000
                         |-----------|
                           ObsPy 반환
                          -100 ~ +100 sample
                          
                          
## ouput length : 2xshift + 1

참고 lag=0 은 두 배열이 첫항부터 정렬되었을때를 말함. 아래 참고

lag = -500

b: b0 ... b500
                ↓
a:              a0 a1 ... a1000

딱 1 sample 겹침


a: a0 a1 a2 ... a500 ... a1000
b: b0 b1 b2 ... b500
   ↑  ↑  ↑       ↑

a[0]과 b[0]이 맞음

lag = +500

a: a0 ... a500 ... a1000
             ↑
b:           b0 ... b500

이때 b가 a 안쪽 뒤쪽에 위치

lag = +1000

a: a0 ....................... a1000
                              ↑
b:                            b0 ... b500

딱 1 sample 겹침

_xcorr_padzeros()

valid mode의 결과 길이가 정확하게 원하는 2*shift + 1 정도가 되도록 두 신호 중 하나를 zero padding한다.

def _xcorr_padzeros(a, b, shift, method):

    if shift is None:
        shift = (len(a) + len(b) - 1) // 2

    dif = len(a) - len(b) - 2 * shift

    if dif > 0:
        b = _pad_zeros(b, dif // 2)
    else:
        a = _pad_zeros(a, -dif // 2)

    return scipy.signal.correlate(
        a,
        b,
        mode='valid',
        method=method
    )

 

normalize='naive'

norm = (
    np.sum(a ** 2) *
    np.sum(b ** 2)
) ** 0.5


cc /= norm

 

 

 

따라서 아래와 같은 pipeline으로 진행

 

cc = correlate(
    a,
    b,
    shift=100,
    demean=True,
    normalize='naive',
    method='auto'
)

correlate()
   │
   ├─ a,b numpy array 변환
   │
   ├─ mean 제거
   │
   │    a -= mean(a)
   │    b -= mean(b)
   │
   ├─ method != 'direct'
   │
   └─ _xcorr_slice()
          │
          ├─ signal.correlate(
          │      a,
          │      b,
          │      mode='full',
          │      method='auto'
          │   )
          │
          │        ↓
          │   SciPy가 direct / fft 선택
          │
          └─ ±100 samples slice

         ↓

normalize

         ↓

201개 CC 반환