Cross-Correlation Differential Time ESTimator (ccdtest)
ccdtest v1.0.1 (January 11, 2026) is permanently archived on Zenodo with DOI: 10.5281/zenodo.10546486
This program measures the relative times between similar signals on a single station (or multiple closely spaced stations - i.e. an array) for two different events. We need an input file of the following format:
WS 3.5 0.2
WS 3.0 0.2
WS 3.0 0.2
WS 2.5 0.25
FS 0.8 2.2 4 2
FS 1.0 2.5 2 2
FS 1.2 2.8 4 2
FS 1.4 3.5 4 2
FS 1.8 4.0 4 2
FS 2.2 4.5 4 2
WF DPRK5_DPRK6_ILAR_example/DPRK6_IL01_IM_--_SHZ_.sac DPRK5_DPRK6_ILAR_example/DPRK5_IL01_IM_--_SHZ_.sac
Each line beginning WS is a window specification and is followed by two numbers:
(a) the length of the window in seconds, and
(b) the spacing between multiple occurrences of the window, in seconds.
Each line beginning FS is a filter specification and is followed by:
(a) low frequency (Hz)
(b) high frequency (Hz)
(c) order (integer)
(d) number of passes (see XAPIIR documentation for all these parameters)
Each line beginning WF is followed by the names of two SAC files.
You then need to run the ccdtest program using a little script such as:
#!/bin/sh
ev1=DPRK6
ev2=DPRK5
sta=IL01
pha=P
dtxtem=120.0
dtxtar=120.0
dtmpwn=2.0
dtarwn=3.5
icc=2
ccdtest $ev1 $ev2 $sta $pha $dtxtem $dtxtar $dtmpwn $dtarwn $icc < ccdtest.input
The parameter icc can currently take 4 values: 1, 2, 3, or 4.
icc = 1 calculates the fully normalized correlation coefficient C.
icc = 2 calculates C|C|.
icc = 3 calculates the phase correlation coefficient PCC (based on Schimmel, 1999, https://doi.org/10.1785/BSSA0890051366)
icc = 4 calculates PCC|PCC|.
Note that you need the fftpack library to be able to run icc = 3 or icc = 4.
If you know that you will not need this then you can copy the file no_pcc_ccdtest.f to ccdtest.f.
This should generate the output as follows with the example given:
DPRK6 DPRK5 2017-09-03T03:39:05.6499 2016-09-09T00:39:05.2087 IL01 P 0.6726 -31028400.4412
(only the final line of the output is shown.)
Modification 2024/07/04
Slight change in the output at intermediate levels - the final line is unchanged.
Instead of just reporting the grid node with the highest CC value we interpolate the true maximum time like we do for the stack and we write out this line for each correlation calculation.
The following terms are appended on the end of each CC-time report:
" -0.0008 calc 456 4 6 1 16 1.7500 4.2400"
dt (time relative to the total stack-determined time)
NCALC (456)
IWS (4 = id of window length specification)
IFS (6 = id of frequency specification)
IWP (1 = id of waveform pair)
ICORR (16 = id of window position specification)
Dtemp_start ( 1.7500 = time of template start relative to the reference template start)
DTtemp_end (4.2400 = time of template end relative to the reference template start)
We also calculate stacks for each combination of IFS (frequency band), IWS (window specification), and ICORR (window placement) over all channels present - where only the maximum closest to the maximum of the global stack is sought. These add the following terms to the end of the corresponding CC-time reports:
" -0.0084 local 4 6 16 1 1.7500 4.2400"
dt (time relative to the total stack-determined time)
IWS (4 = id of window length specification)
IFS (6 = id of frequency specification)
ICORR (16 = id of window position specification)
NCFSEL (number of channels collected)
Dtemp_start ( 1.7500 = time of template start relative to the reference template start)
DTtemp_end (4.2400 = time of template end relative to the reference template start)
This program requires SAC and uses the XAPIIR library (Harris, 1990).
This whole library is provided in the single source file XAPIIR.f.
(The XAPIIR routines are entirely third party software.)
It also requires the LAPACK and BLAS libraries - although all the routines needed are included in the linear algebra part of https://github.com/stevenjgibbons/LEOPACK-2022-revision
Modification 2025/12/30
Writes out a file in SAC format outputcc.sac containing the stack CC trace over the short time-window.
The starting epoch time of this trace is target time - template reference time.
References:
Harris, David. XAPiir: A recursive digital filtering package, report, September 21, 1990; California.
(https://digital.library.unt.edu/ark:/67531/metadc1203741/m1/1/: accessed November 5, 2023),
University of North Texas Libraries, UNT Digital Library,
https://digital.library.unt.edu; crediting UNT Libraries Government Documents Department.
If you publish using this code, I kindly ask you to cite the following paper:
Relative location of seismic events using broadband frequency-wavenumber analysis: Application to the North Korean Nuclear Test Site
Steven J. Gibbons, RAS Techniques and Instruments, rzag023, https://doi.org/10.1093/rasti/rzag023
