-
Notifications
You must be signed in to change notification settings - Fork 17
/
harminv.1
166 lines (155 loc) · 6.64 KB
/
harminv.1
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
.\" Copyright (c) 2000 Massachusetts Institute of Technology
.\"
.\" This program is free software; you can redistribute it and/or modify
.\" it under the terms of the GNU General Public License as published by
.\" the Free Software Foundation; either version 2 of the License, or
.\" (at your option) any later version.
.\"
.\" This program is distributed in the hope that it will be useful,
.\" but WITHOUT ANY WARRANTY; without even the implied warranty of
.\" MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
.\" GNU General Public License for more details.
.\"
.\" You should have received a copy of the GNU General Public License
.\" along with this program; if not, write to the Free Software
.\" Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
.\"
.TH HARMINV 1 "June 4, 2004" "harminv" "harminv"
.SH NAME
harminv \- extract mode frequencies from time-series data
.SH SYNOPSIS
.B harminv
[\fIOPTION\fR]... [\fIfreq-min\fR-\fIfreq-max\fR]...
.SH DESCRIPTION
.PP
." Add any additional description here
\fBharminv\fR is a program designed to solve the problem of "harmonic
inversion": given a time series consisting of a sum of sinusoids
("modes"), extract their frequencies and amplitudes. It can also
handle the case of exponentially-decaying sinusoids, in which case it
extracts their decay rates as well.
\fBharminv\fR is often able to achieve much greater accuracy and
robustness than Fourier-transform methods, essentially because it
assumes a specific form for the input.
It uses a low-storage "filter-diagonalization method" (FDM), as
described in V. A. Mandelshtam and H. S. Taylor, "Harmonic inversion
of time signals," \fIJ. Chem. Phys.\fR \fB107\fR, 6756 (1997). See
also erratum, \fIibid\fR \fB109\fR, 4128 (1998).
.SH INPUT
\fBharminv\fR reads in a sequence of whitespace-separated real or
complex numbers from standard input, as well as command-line arguments
indicating one or more frequency ranges to search, and outputs the
modes that it extracts from the data. (It preferentially finds modes
in the frequency range you specify, but may sometimes find additional
modes outside of that range.) The data should correspond to
equally-spaced time intervals, but there is no constraint on the
number of points.
Complex numbers in the input should be expressed in the format
\fIRE\fR+\fIIM\fRi (no whitespace). Otherwise, whitespace is ignored.
Also, comments beginning with "#" and extending to the end of the line
are ignored.
A typical invocation is something like
.IP "" 4
harminv -t 0.02 1-5 < input.dat
.PP
which reads a sequence of samples, spaced at 0.02 time intervals (in
ms, say, corresponding to 50 kHz), and searches for modes in the
frequency range 1-5 kHz. (See below on units.)
.SH OUTPUT
\fBharminv\fR writes four comma-delimited columns to standard output, one
line for each mode: frequency, decay constant, Q, amplitude, phase,
and error. Each mode corresponds to a function of the form:
\fIamplitude\fR * exp[-i (2 pi \fIfrequency\fR t - \fIphase\fR) - \fIdecay\fR t]
Here, i is sqrt(-1), t is the time (see below for units), and the
other parameters in the output columns are:
.TP
.B frequency
The frequency of the mode. If you don't recognize that from the
expression above, you should recall Euler's formula: exp(i x) = cos(x)
+ i sin(x). Note that for complex data, there is a distinction between
positive and negative frequencies.
.TP
.B decay constant
The exponential decay constant, indicated by
.I decay
in the above formula. The inverse of this is often called the
"lifetime" of the mode. The "half-life" is ln(2)/\fIdecay\fR.
.TP
.B Q
A conventional, dimensionless expression of the decay lifetime: Q = pi
|\fIfrequency\fI| / \fIdecay\fR. Q is the number of periods for the
"energy" in the mode (the squared amplitude) to decay by exp(-2 pi).
Equivalently, if you look at the power spectrum (|Fourier
transform|^2), 1/Q is the fractional width of the peak at half maximum.
.TP
.B amplitude
The (real, positive) amplitude of the sinusoids. The amplitude (and
phase) information generally seems to be less accurate than the
frequency and decay constant.
.TP
.B phase
The phase shift (in radians) of the sinusoids, as given by the formula
above.
.TP
.B error
An estimate of the degree of confidence in the mode information
(larger error means less confidence). A precise quantitative
interpretation is difficult, unfortunately. However, modes with
errors much larger than the modes with the smallest errors may be
spurious (e.g. because your data don't exactly fit the assumed model).
.SH UNITS
The frequency (and decay) values, both input and output, are specified
in units of 1/time, where the units of time are determined by the
sampling interval \fIdt\fR (the time between consecutive inputs).
\fIdt\fR is by default 1, unless you specify it with the
.B -t
.I dt
option.
In other words, pick some units (e.g. ms in the example above) and use
them to express the time step. Then, be consistent and use the
inverse of those units (e.g. kHz = 1/ms) for frequency.
Note that the frequency is the usual 1/period definition; it is not
the angular frequency.
.SH OPTIONS
.TP
.B -h
Display help on the command-line options and usage.
.TP
.B -V
Print the version number and copyright info for \fBharminv\fR.
.TP
.B -v
Enable verbose output, printed to standard output as comment lines
(starting with a "#" character). Also, any "#" comments in the input
are echoed to the output.
.TP
.B -T
Specify period-ranges instead of frequency-ranges on the command line
(in units of time corresponding to those specified by \fB-t\fR). The
output is still frequency and not period, however.
.TP
\fB\-t\fR \fIdt\fR
Specify the sampling interval \fIdt\fR; this determines the units of
time used throughout the input and output. Defaults to 1.0.
.TP
\fB\-f\fR \fInf\fR
Specify the "spectral density," essentially the number of frequencies
to look for modes nearby within the frequency range. This should be
at least as great as the number of modes you expect to find, and
usually you should just set it to a much greater number. \fInf\fR
must be at least 2, and defaults to 100.
Don't set it to something too large, however: the computation time
scales as O(N * \fInf\fR + \fInf\fR^3), where N is the number of data
points. Also, there may be numerical problems if \fInf\fR is too
large.
.TP
\fB\-s\fR \fIsort\fR
Specify how the outputs are sorted, where \fIsort\fR is one of
frequency/error/decay/amplitude. (Only the first character of
\fIsort\fR matters.) All sorts are in ascending order. The default
is to sort by frequency.
.SH BUGS
Send bug reports to S. G. Johnson, stevenj@alum.mit.edu.
.SH AUTHORS
Written by Steven G. Johnson. Copyright (c) 2004 by the Massachusetts
Institute of Technology.