line |
stmt |
bran |
cond |
sub |
pod |
time |
code |
1
|
|
|
|
|
|
|
# This file generated by mkwindows.pl |
2
|
|
|
|
|
|
|
package PDL::DSP::Windows; |
3
|
|
|
|
|
|
|
$PDL::DSP::Windows::VERSION = '0.007'; |
4
|
2
|
|
|
2
|
|
228614
|
use base 'Exporter'; |
|
2
|
|
|
|
|
4
|
|
|
2
|
|
|
|
|
193
|
|
5
|
2
|
|
|
2
|
|
10
|
use strict; use warnings; |
|
2
|
|
|
2
|
|
4
|
|
|
2
|
|
|
|
|
43
|
|
|
2
|
|
|
|
|
10
|
|
|
2
|
|
|
|
|
6
|
|
|
2
|
|
|
|
|
57
|
|
6
|
2
|
|
|
2
|
|
706
|
use PDL::LiteF; |
|
2
|
|
|
|
|
877
|
|
|
2
|
|
|
|
|
11
|
|
7
|
2
|
|
|
2
|
|
187755
|
use PDL::FFT; |
|
2
|
|
|
|
|
16035
|
|
|
2
|
|
|
|
|
15
|
|
8
|
2
|
|
|
2
|
|
448
|
use PDL::Math qw( acos cosh acosh ); |
|
2
|
|
|
|
|
4
|
|
|
2
|
|
|
|
|
10
|
|
9
|
2
|
|
|
2
|
|
221
|
use PDL::Core qw( topdl ); |
|
2
|
|
|
|
|
4
|
|
|
2
|
|
|
|
|
9
|
|
10
|
2
|
|
|
2
|
|
154
|
use PDL::MatrixOps qw( eigens_sym ); |
|
2
|
|
|
|
|
4
|
|
|
2
|
|
|
|
|
18
|
|
11
|
2
|
|
|
2
|
|
153
|
use PDL::Options qw( iparse ifhref ); |
|
2
|
|
|
|
|
3
|
|
|
2
|
|
|
|
|
279
|
|
12
|
|
|
|
|
|
|
|
13
|
|
|
|
|
|
|
eval { require PDL::LinearAlgebra::Special }; |
14
|
|
|
|
|
|
|
my $HAVE_LinearAlgebra = 1 if !$@; |
15
|
|
|
|
|
|
|
|
16
|
|
|
|
|
|
|
eval { require PDL::GSLSF::BESSEL; }; |
17
|
|
|
|
|
|
|
my $HAVE_BESSEL = 1 if !$@; |
18
|
|
|
|
|
|
|
|
19
|
|
|
|
|
|
|
eval { require PDL::Graphics::Gnuplot; }; |
20
|
|
|
|
|
|
|
my $HAVE_GNUPLOT = 1 if !$@; |
21
|
|
|
|
|
|
|
|
22
|
|
|
|
|
|
|
#eval { require PDL::Graphics::PLplot; }; |
23
|
|
|
|
|
|
|
#my $HAVE_PLPLOT = 1 if !$@; |
24
|
|
|
|
|
|
|
|
25
|
2
|
|
|
2
|
|
1713
|
use PDL::Constants qw(PI); |
|
2
|
|
|
|
|
23601
|
|
|
2
|
|
|
|
|
168
|
|
26
|
2
|
|
|
2
|
|
18
|
use constant TPI => 2 * PI; |
|
2
|
|
|
|
|
2
|
|
|
2
|
|
|
|
|
25187
|
|
27
|
|
|
|
|
|
|
|
28
|
|
|
|
|
|
|
our @ISA = qw(Exporter); |
29
|
|
|
|
|
|
|
|
30
|
|
|
|
|
|
|
our @EXPORT_OK = qw( window list_windows chebpoly cos_mult_to_pow cos_pow_to_mult |
31
|
|
|
|
|
|
|
bartlett bartlett_hann bartlett_hann_per bartlett_per blackman blackman_bnh blackman_bnh_per blackman_ex blackman_ex_per blackman_gen blackman_gen3 blackman_gen3_per blackman_gen4 blackman_gen4_per blackman_gen5 blackman_gen5_per blackman_gen_per blackman_harris blackman_harris4 blackman_harris4_per blackman_harris_per blackman_nuttall blackman_nuttall_per blackman_per bohman bohman_per cauchy cauchy_per chebyshev cos_alpha cos_alpha_per cosine cosine_per dpss dpss_per exponential exponential_per flattop flattop_per gaussian gaussian_per hamming hamming_ex hamming_ex_per hamming_gen hamming_gen_per hamming_per hann hann_matlab hann_per hann_poisson hann_poisson_per kaiser kaiser_per lanczos lanczos_per nuttall nuttall1 nuttall1_per nuttall_per parzen parzen_octave parzen_per poisson poisson_per rectangular rectangular_per triangular triangular_per tukey tukey_per welch welch_per); |
32
|
|
|
|
|
|
|
|
33
|
|
|
|
|
|
|
$PDL::onlinedoc->scan(__FILE__) if $PDL::onlinedoc; |
34
|
|
|
|
|
|
|
|
35
|
|
|
|
|
|
|
our %winsubs; |
36
|
|
|
|
|
|
|
our %winpersubs; |
37
|
|
|
|
|
|
|
our %window_definitions; |
38
|
|
|
|
|
|
|
|
39
|
|
|
|
|
|
|
=encoding utf8 |
40
|
|
|
|
|
|
|
|
41
|
|
|
|
|
|
|
=head1 NAME |
42
|
|
|
|
|
|
|
|
43
|
|
|
|
|
|
|
PDL::DSP::Windows - Window functions for signal processing |
44
|
|
|
|
|
|
|
|
45
|
|
|
|
|
|
|
=head1 SYNOPSIS |
46
|
|
|
|
|
|
|
|
47
|
|
|
|
|
|
|
use PDL; |
48
|
|
|
|
|
|
|
use PDL::DSP::Windows('window'); |
49
|
|
|
|
|
|
|
my $samples = window( 10, 'tukey', { params => .5 }); |
50
|
|
|
|
|
|
|
|
51
|
|
|
|
|
|
|
use PDL; |
52
|
|
|
|
|
|
|
use PDL::DSP::Windows; |
53
|
|
|
|
|
|
|
my $win = new PDL::DSP::Windows(10, 'tukey', { params => .5 }); |
54
|
|
|
|
|
|
|
print $win->coherent_gain , "\n"; |
55
|
|
|
|
|
|
|
$win->plot; |
56
|
|
|
|
|
|
|
|
57
|
|
|
|
|
|
|
=head1 DESCRIPTION |
58
|
|
|
|
|
|
|
|
59
|
|
|
|
|
|
|
This module provides symmetric and periodic (DFT-symmetric) |
60
|
|
|
|
|
|
|
window functions for use in filtering and spectral analysis. |
61
|
|
|
|
|
|
|
It provides a high-level access subroutine |
62
|
|
|
|
|
|
|
L. This functional interface is sufficient for getting the window |
63
|
|
|
|
|
|
|
samples. For analysis and plotting, etc. an object oriented |
64
|
|
|
|
|
|
|
interface is provided. The functional subroutines must be either explicitly exported, or |
65
|
|
|
|
|
|
|
fully qualified. In this document, the word I refers only to the |
66
|
|
|
|
|
|
|
mathematical window functions, while the word I is used to describe |
67
|
|
|
|
|
|
|
code. |
68
|
|
|
|
|
|
|
|
69
|
|
|
|
|
|
|
Window functions are also known as apodization |
70
|
|
|
|
|
|
|
functions or tapering functions. In this module, each of these |
71
|
|
|
|
|
|
|
functions maps a sequence of C<$N> integers to values called |
72
|
|
|
|
|
|
|
a B. (To confuse matters, the word I also has |
73
|
|
|
|
|
|
|
other meanings when describing window functions.) |
74
|
|
|
|
|
|
|
The functions are often named for authors of journal articles. |
75
|
|
|
|
|
|
|
Be aware that across the literature and software, |
76
|
|
|
|
|
|
|
some functions referred to by several different names, and some names |
77
|
|
|
|
|
|
|
refer to several different functions. As a result, the choice |
78
|
|
|
|
|
|
|
of window names is somewhat arbitrary. |
79
|
|
|
|
|
|
|
|
80
|
|
|
|
|
|
|
The L window function requires |
81
|
|
|
|
|
|
|
L. The L window function requires |
82
|
|
|
|
|
|
|
L. But the remaining window functions may |
83
|
|
|
|
|
|
|
be used if these modules are not installed. |
84
|
|
|
|
|
|
|
|
85
|
|
|
|
|
|
|
The most common and easiest usage of this module is indirect, via some |
86
|
|
|
|
|
|
|
higher-level filtering interface, such as L. |
87
|
|
|
|
|
|
|
The next easiest usage is to return a pdl of real-space samples with the subroutine L. |
88
|
|
|
|
|
|
|
Finally, for analyzing window functions, object methods, such as L, |
89
|
|
|
|
|
|
|
L, L are provided. |
90
|
|
|
|
|
|
|
|
91
|
|
|
|
|
|
|
In the following, first the functional interface (non-object oriented) is described in |
92
|
|
|
|
|
|
|
L"FUNCTIONAL INTERFACE">. Next, the object methods are described in L. |
93
|
|
|
|
|
|
|
Next the low-level subroutines returning samples for each named window |
94
|
|
|
|
|
|
|
are described in L"WINDOW FUNCTIONS">. Finally, |
95
|
|
|
|
|
|
|
some support routines that may be of interest are described in |
96
|
|
|
|
|
|
|
L"AUXILIARY SUBROUTINES">. |
97
|
|
|
|
|
|
|
|
98
|
|
|
|
|
|
|
=head1 FUNCTIONAL INTERFACE |
99
|
|
|
|
|
|
|
|
100
|
|
|
|
|
|
|
=head2 window |
101
|
|
|
|
|
|
|
|
102
|
|
|
|
|
|
|
$win = window({OPTIONS}); |
103
|
|
|
|
|
|
|
$win = window($N,{OPTIONS}); |
104
|
|
|
|
|
|
|
$win = window($N,$name,{OPTIONS}); |
105
|
|
|
|
|
|
|
$win = window($N,$name,$params,{OPTIONS}); |
106
|
|
|
|
|
|
|
$win = window($N,$name,$params,$periodic); |
107
|
|
|
|
|
|
|
|
108
|
|
|
|
|
|
|
Returns an C<$N> point window of type C<$name>. |
109
|
|
|
|
|
|
|
The arguments may be passed positionally in the order |
110
|
|
|
|
|
|
|
C<$N,$name,$params,$periodic>, or they may be passed by |
111
|
|
|
|
|
|
|
name in the hash C. |
112
|
|
|
|
|
|
|
|
113
|
|
|
|
|
|
|
=head3 EXAMPLES |
114
|
|
|
|
|
|
|
|
115
|
|
|
|
|
|
|
# Each of the following return a 100 point symmetric hamming window. |
116
|
|
|
|
|
|
|
|
117
|
|
|
|
|
|
|
$win = window(100); |
118
|
|
|
|
|
|
|
$win = window(100, 'hamming'); |
119
|
|
|
|
|
|
|
$win = window(100, { name => 'hamming' ); |
120
|
|
|
|
|
|
|
$win = window({ N=> 100, name => 'hamming' ); |
121
|
|
|
|
|
|
|
|
122
|
|
|
|
|
|
|
# Each of the following returns a 100 point symmetric hann window. |
123
|
|
|
|
|
|
|
|
124
|
|
|
|
|
|
|
$win = window(100, 'hann'); |
125
|
|
|
|
|
|
|
$win = window(100, { name => 'hann' ); |
126
|
|
|
|
|
|
|
|
127
|
|
|
|
|
|
|
# Returns a 100 point periodic hann window. |
128
|
|
|
|
|
|
|
|
129
|
|
|
|
|
|
|
$win = window(100, 'hann', { periodic => 1 } ); |
130
|
|
|
|
|
|
|
|
131
|
|
|
|
|
|
|
# Returns a 100 point symmetric Kaiser window with alpha=2. |
132
|
|
|
|
|
|
|
|
133
|
|
|
|
|
|
|
$win = window(100, 'kaiser', { params => 2 }); |
134
|
|
|
|
|
|
|
|
135
|
|
|
|
|
|
|
=head3 OPTIONS |
136
|
|
|
|
|
|
|
|
137
|
|
|
|
|
|
|
The options follow default PDL::Options rules-- They may be abbreviated, |
138
|
|
|
|
|
|
|
and are case-insensitive. |
139
|
|
|
|
|
|
|
|
140
|
|
|
|
|
|
|
=over |
141
|
|
|
|
|
|
|
|
142
|
|
|
|
|
|
|
=item B |
143
|
|
|
|
|
|
|
|
144
|
|
|
|
|
|
|
(string) name of window function. Default: C. |
145
|
|
|
|
|
|
|
This selects one of the window functions listed below. Note |
146
|
|
|
|
|
|
|
that the suffix '_per', for periodic, may be ommitted. It |
147
|
|
|
|
|
|
|
is specified with the option C<< periodic => 1 >> |
148
|
|
|
|
|
|
|
|
149
|
|
|
|
|
|
|
|
150
|
|
|
|
|
|
|
=item B |
151
|
|
|
|
|
|
|
|
152
|
|
|
|
|
|
|
|
153
|
|
|
|
|
|
|
ref to array of parameter or parameters for the window-function |
154
|
|
|
|
|
|
|
subroutine. Only some window-function subroutines take |
155
|
|
|
|
|
|
|
parameters. If the subroutine takes a single parameter, |
156
|
|
|
|
|
|
|
it may be given either as a number, or a list of one |
157
|
|
|
|
|
|
|
number. For example C<3> or C<[3]>. |
158
|
|
|
|
|
|
|
|
159
|
|
|
|
|
|
|
=item B |
160
|
|
|
|
|
|
|
|
161
|
|
|
|
|
|
|
number of points in window function (the same as the order |
162
|
|
|
|
|
|
|
of the filter) No default value. |
163
|
|
|
|
|
|
|
|
164
|
|
|
|
|
|
|
=item B |
165
|
|
|
|
|
|
|
|
166
|
|
|
|
|
|
|
If value is true, return a periodic rather than a symmetric window function. Default: 0 |
167
|
|
|
|
|
|
|
(that is, false. that is, symmetric.) |
168
|
|
|
|
|
|
|
|
169
|
|
|
|
|
|
|
=back |
170
|
|
|
|
|
|
|
|
171
|
|
|
|
|
|
|
=cut |
172
|
|
|
|
|
|
|
|
173
|
|
|
|
|
|
|
sub window { |
174
|
143
|
|
|
143
|
1
|
62112
|
my $win = new PDL::DSP::Windows(@_); |
175
|
143
|
|
|
|
|
320
|
$win->samples(); |
176
|
|
|
|
|
|
|
} |
177
|
|
|
|
|
|
|
|
178
|
|
|
|
|
|
|
=head2 list_windows |
179
|
|
|
|
|
|
|
|
180
|
|
|
|
|
|
|
list_windows |
181
|
|
|
|
|
|
|
list_windows STR |
182
|
|
|
|
|
|
|
|
183
|
|
|
|
|
|
|
C prints the names all of the available windows. |
184
|
|
|
|
|
|
|
C prints only the names of windows matching |
185
|
|
|
|
|
|
|
the string C. |
186
|
|
|
|
|
|
|
|
187
|
|
|
|
|
|
|
=cut |
188
|
|
|
|
|
|
|
|
189
|
|
|
|
|
|
|
sub list_windows { |
190
|
0
|
|
|
0
|
1
|
0
|
my ($expr) = @_; |
191
|
0
|
|
|
|
|
0
|
my @match; |
192
|
0
|
0
|
|
|
|
0
|
if ($expr) { |
193
|
0
|
|
|
|
|
0
|
my @alias; |
194
|
0
|
|
|
|
|
0
|
foreach (sort keys %winsubs) { |
195
|
0
|
0
|
|
|
|
0
|
push(@match,$_) , next if /$expr/i; |
196
|
0
|
0
|
|
|
|
0
|
push(@match, $_ . ' (alias ' . $alias[0] . ')') if @alias = grep(/$expr/i,@{$window_definitions{$_}->{alias}}); |
|
0
|
|
|
|
|
0
|
|
197
|
|
|
|
|
|
|
} |
198
|
|
|
|
|
|
|
} |
199
|
|
|
|
|
|
|
else { |
200
|
0
|
|
|
|
|
0
|
@match = sort keys %winsubs; |
201
|
|
|
|
|
|
|
} |
202
|
0
|
|
|
|
|
0
|
print join(', ',@match),"\n"; |
203
|
|
|
|
|
|
|
} |
204
|
|
|
|
|
|
|
|
205
|
|
|
|
|
|
|
|
206
|
|
|
|
|
|
|
=head1 METHODS |
207
|
|
|
|
|
|
|
|
208
|
|
|
|
|
|
|
=head2 new |
209
|
|
|
|
|
|
|
|
210
|
|
|
|
|
|
|
=for usage |
211
|
|
|
|
|
|
|
|
212
|
|
|
|
|
|
|
my $win = new PDL::DSP::Windows(ARGS); |
213
|
|
|
|
|
|
|
|
214
|
|
|
|
|
|
|
=for ref |
215
|
|
|
|
|
|
|
|
216
|
|
|
|
|
|
|
Create an instance of a Windows object. If C are given, the instance |
217
|
|
|
|
|
|
|
is initialized. C are interpreted in exactly the |
218
|
|
|
|
|
|
|
same way as arguments the subroutine L. |
219
|
|
|
|
|
|
|
|
220
|
|
|
|
|
|
|
=for example |
221
|
|
|
|
|
|
|
|
222
|
|
|
|
|
|
|
For example: |
223
|
|
|
|
|
|
|
|
224
|
|
|
|
|
|
|
my $win1 = new PDL::DSP::Windows(8,'hann'); |
225
|
|
|
|
|
|
|
my $win2 = new PDL::DSP::Windows( { N => 8, name => 'hann' } ); |
226
|
|
|
|
|
|
|
|
227
|
|
|
|
|
|
|
=cut |
228
|
|
|
|
|
|
|
|
229
|
|
|
|
|
|
|
sub new { |
230
|
146
|
|
|
146
|
1
|
3409
|
my $proto = shift; |
231
|
146
|
|
33
|
|
|
626
|
my $class = ref($proto) || $proto; |
232
|
146
|
|
|
|
|
251
|
my $self = {}; |
233
|
146
|
|
|
|
|
268
|
bless ($self, $class); |
234
|
146
|
100
|
|
|
|
539
|
$self->init(@_) if (@_); |
235
|
146
|
|
|
|
|
271
|
return $self; |
236
|
|
|
|
|
|
|
} |
237
|
|
|
|
|
|
|
|
238
|
|
|
|
|
|
|
=head2 init |
239
|
|
|
|
|
|
|
|
240
|
|
|
|
|
|
|
=for usage |
241
|
|
|
|
|
|
|
|
242
|
|
|
|
|
|
|
$win->init(ARGS); |
243
|
|
|
|
|
|
|
|
244
|
|
|
|
|
|
|
=for ref |
245
|
|
|
|
|
|
|
|
246
|
|
|
|
|
|
|
Initialize (or reinitialize) a Windows object. ARGS are interpreted in exactly the |
247
|
|
|
|
|
|
|
same way as arguments the subroutine L. |
248
|
|
|
|
|
|
|
|
249
|
|
|
|
|
|
|
=for example |
250
|
|
|
|
|
|
|
|
251
|
|
|
|
|
|
|
For example: |
252
|
|
|
|
|
|
|
|
253
|
|
|
|
|
|
|
my $win = new PDL::DSP::Windows(8,'hann'); |
254
|
|
|
|
|
|
|
$win->init(10,'hamming'); |
255
|
|
|
|
|
|
|
|
256
|
|
|
|
|
|
|
=cut |
257
|
|
|
|
|
|
|
|
258
|
|
|
|
|
|
|
sub init { |
259
|
159
|
|
|
159
|
1
|
8897
|
my $self = shift; |
260
|
|
|
|
|
|
|
|
261
|
159
|
|
|
|
|
900
|
my $opt = new PDL::Options( |
262
|
|
|
|
|
|
|
{ |
263
|
|
|
|
|
|
|
name => 'hamming', |
264
|
|
|
|
|
|
|
periodic => 0, # symmetric or periodic |
265
|
|
|
|
|
|
|
N => undef, # order |
266
|
|
|
|
|
|
|
params => undef, |
267
|
|
|
|
|
|
|
}); |
268
|
159
|
|
|
|
|
7688
|
my ($N,$name,$params,$periodic); |
269
|
159
|
100
|
|
|
|
419
|
$N = shift unless ref ($_[0]); |
270
|
159
|
100
|
|
|
|
389
|
$name = shift unless ref ($_[0]); |
271
|
159
|
100
|
|
|
|
357
|
$params = shift unless ref ($_[0]) eq 'HASH'; |
272
|
159
|
100
|
|
|
|
310
|
$periodic = shift unless ref ($_[0]); |
273
|
159
|
100
|
|
|
|
342
|
my $iopts = @_ ? shift : {}; |
274
|
159
|
|
|
|
|
470
|
my $opts = $opt->options($iopts); |
275
|
159
|
100
|
|
|
|
22864
|
$name = $opts->{name} unless $name; |
276
|
159
|
|
|
|
|
302
|
$name =~ s/_per$//; |
277
|
159
|
100
|
|
|
|
341
|
$N = $opts->{N} unless $N; |
278
|
159
|
100
|
|
|
|
454
|
$params = $opts->{params} unless defined $params; |
279
|
159
|
100
|
100
|
|
|
503
|
$params = [$params] if defined $params and not ref $params; |
280
|
159
|
50
|
|
|
|
400
|
$periodic = $opts->{periodic} unless $periodic; |
281
|
159
|
100
|
|
|
|
346
|
my $ws = $periodic ? \%winpersubs : \%winsubs; |
282
|
159
|
50
|
|
|
|
489
|
if ( not exists $ws->{$name}) { |
283
|
0
|
0
|
|
|
|
0
|
my $perstr = $periodic ? 'periodic' : 'symmetric'; |
284
|
0
|
|
|
|
|
0
|
barf "window: Unknown $perstr window '$name'."; |
285
|
|
|
|
|
|
|
} |
286
|
159
|
|
|
|
|
315
|
$self->{name} = $name; |
287
|
159
|
|
|
|
|
248
|
$self->{N} = $N; |
288
|
159
|
|
|
|
|
218
|
$self->{periodic} = $periodic; |
289
|
159
|
|
|
|
|
234
|
$self->{params} = $params; |
290
|
159
|
|
|
|
|
333
|
$self->{code} = $ws->{$name}; |
291
|
159
|
|
|
|
|
280
|
$self->{samples} = undef; |
292
|
159
|
|
|
|
|
665
|
$self->{modfreqs} = undef; |
293
|
159
|
|
|
|
|
973
|
return $self; |
294
|
|
|
|
|
|
|
} |
295
|
|
|
|
|
|
|
|
296
|
|
|
|
|
|
|
=head2 samples |
297
|
|
|
|
|
|
|
|
298
|
|
|
|
|
|
|
=for usage |
299
|
|
|
|
|
|
|
|
300
|
|
|
|
|
|
|
$win->samples(); |
301
|
|
|
|
|
|
|
|
302
|
|
|
|
|
|
|
=for ref |
303
|
|
|
|
|
|
|
|
304
|
|
|
|
|
|
|
Generate and return a reference to the piddle of $N samples for the window C<$win>. |
305
|
|
|
|
|
|
|
This is the real-space representation of the window. |
306
|
|
|
|
|
|
|
|
307
|
|
|
|
|
|
|
The samples are stored in the object C<$win>, but are regenerated |
308
|
|
|
|
|
|
|
every time C is invoked. See the method |
309
|
|
|
|
|
|
|
L below. |
310
|
|
|
|
|
|
|
|
311
|
|
|
|
|
|
|
=for example |
312
|
|
|
|
|
|
|
|
313
|
|
|
|
|
|
|
For example: |
314
|
|
|
|
|
|
|
|
315
|
|
|
|
|
|
|
my $win = new PDL::DSP::Windows(8,'hann'); |
316
|
|
|
|
|
|
|
print $win->samples(), "\n"; |
317
|
|
|
|
|
|
|
|
318
|
|
|
|
|
|
|
=cut |
319
|
|
|
|
|
|
|
|
320
|
|
|
|
|
|
|
sub samples { |
321
|
159
|
|
|
159
|
1
|
220
|
my $self = shift; |
322
|
159
|
100
|
|
|
|
430
|
my @args = defined $self->{params} ? ($self->{N}, @{$self->{params}} ) : ($self->{N}); |
|
50
|
|
|
|
|
129
|
|
323
|
159
|
|
|
|
|
504
|
$self->{samples} = $self->{code}->(@args); |
324
|
|
|
|
|
|
|
} |
325
|
|
|
|
|
|
|
|
326
|
|
|
|
|
|
|
=head2 modfreqs |
327
|
|
|
|
|
|
|
|
328
|
|
|
|
|
|
|
=for usage |
329
|
|
|
|
|
|
|
|
330
|
|
|
|
|
|
|
$win->modfreqs(); |
331
|
|
|
|
|
|
|
|
332
|
|
|
|
|
|
|
=for ref |
333
|
|
|
|
|
|
|
|
334
|
|
|
|
|
|
|
Generate and return a reference to the piddle of the modulus of the |
335
|
|
|
|
|
|
|
fourier transform of the samples for the window C<$win>. |
336
|
|
|
|
|
|
|
|
337
|
|
|
|
|
|
|
These values are stored in the object C<$win>, but are regenerated |
338
|
|
|
|
|
|
|
every time C is invoked. See the method |
339
|
|
|
|
|
|
|
L below. |
340
|
|
|
|
|
|
|
|
341
|
|
|
|
|
|
|
=head3 options |
342
|
|
|
|
|
|
|
|
343
|
|
|
|
|
|
|
=over |
344
|
|
|
|
|
|
|
|
345
|
|
|
|
|
|
|
=item min_bins => MIN |
346
|
|
|
|
|
|
|
|
347
|
|
|
|
|
|
|
This sets the minimum number of frequency bins. |
348
|
|
|
|
|
|
|
Default 1000. If necessary, the piddle of window samples |
349
|
|
|
|
|
|
|
are padded with zeros before the fourier transform is performed. |
350
|
|
|
|
|
|
|
|
351
|
|
|
|
|
|
|
=back |
352
|
|
|
|
|
|
|
|
353
|
|
|
|
|
|
|
=cut |
354
|
|
|
|
|
|
|
|
355
|
|
|
|
|
|
|
sub modfreqs { |
356
|
2
|
|
|
2
|
1
|
15
|
my ($self,$inopts) = @_; |
357
|
2
|
|
|
|
|
11
|
my %opts = iparse( { min_bins => 1000 }, ifhref($inopts)); |
358
|
2
|
|
|
|
|
376
|
my $data = $self->get('samples'); |
359
|
2
|
|
|
|
|
9
|
my $n = $data->nelem; |
360
|
2
|
50
|
|
|
|
10
|
my $fn = $n > $opts{min_bins} ? 2 * $n : $opts{min_bins}; |
361
|
2
|
|
|
|
|
4
|
$n--; |
362
|
2
|
|
|
|
|
6
|
my $freq = zeroes($fn); |
363
|
2
|
|
|
|
|
109
|
$freq->slice("0:$n") .= $data; |
364
|
2
|
|
|
|
|
69
|
PDL::FFT::realfft($freq); |
365
|
2
|
|
|
|
|
351
|
my $real = zeros($freq); |
366
|
2
|
|
|
|
|
107
|
my $img = zeros($freq); |
367
|
2
|
|
|
|
|
103
|
my $mid = ($freq->nelem)/2 - 1; |
368
|
2
|
|
|
|
|
4
|
my $mid1 = $mid + 1; |
369
|
2
|
|
|
|
|
10
|
$real->slice("0:$mid") .= $freq->slice("$mid:0:-1"); |
370
|
2
|
|
|
|
|
88
|
$real->slice("$mid1:-1") .= $freq->slice("0:$mid"); |
371
|
2
|
|
|
|
|
82
|
$img->slice("0:$mid") .= $freq->slice("-1:$mid1:-1"); |
372
|
2
|
|
|
|
|
83
|
$img->slice("$mid1:-1") .= $freq->slice("$mid1:-1"); |
373
|
2
|
|
|
|
|
152
|
my $mod = $real**2 + $img**2; |
374
|
2
|
|
|
|
|
48
|
$self->{modfreqs} = $mod; |
375
|
2
|
|
|
|
|
21
|
return $mod; |
376
|
|
|
|
|
|
|
} |
377
|
|
|
|
|
|
|
|
378
|
|
|
|
|
|
|
=head2 get |
379
|
|
|
|
|
|
|
|
380
|
|
|
|
|
|
|
=for usage |
381
|
|
|
|
|
|
|
|
382
|
|
|
|
|
|
|
my $windata = $win->get('samples'); |
383
|
|
|
|
|
|
|
|
384
|
|
|
|
|
|
|
=for ref |
385
|
|
|
|
|
|
|
|
386
|
|
|
|
|
|
|
Get an attribute (or list of attributes) of the window C<$win>. |
387
|
|
|
|
|
|
|
If attribute C is requested, then the samples are created with the |
388
|
|
|
|
|
|
|
method L if they don't exist. |
389
|
|
|
|
|
|
|
|
390
|
|
|
|
|
|
|
=for example |
391
|
|
|
|
|
|
|
|
392
|
|
|
|
|
|
|
For example: |
393
|
|
|
|
|
|
|
|
394
|
|
|
|
|
|
|
my $win = new PDL::DSP::Windows(8,'hann'); |
395
|
|
|
|
|
|
|
print $win->get('samples'), "\n"; |
396
|
|
|
|
|
|
|
|
397
|
|
|
|
|
|
|
=cut |
398
|
|
|
|
|
|
|
|
399
|
|
|
|
|
|
|
sub get { |
400
|
16
|
|
|
16
|
1
|
23
|
my $self = shift; |
401
|
16
|
|
|
|
|
20
|
my @res; |
402
|
16
|
|
|
|
|
37
|
foreach (@_) { |
403
|
16
|
50
|
33
|
|
|
110
|
$self->samples() if $_ eq 'samples' and not defined $self->{samples}; |
404
|
16
|
50
|
33
|
|
|
36683
|
$self->freqs() if $_ eq 'modfreqs' and not defined $self->{modfreqs}; |
405
|
16
|
|
|
|
|
51
|
push @res, $self->{$_}; |
406
|
|
|
|
|
|
|
}; |
407
|
16
|
50
|
|
|
|
57
|
return wantarray ? @res : $res[0]; |
408
|
|
|
|
|
|
|
} |
409
|
|
|
|
|
|
|
|
410
|
|
|
|
|
|
|
=head2 get_samples |
411
|
|
|
|
|
|
|
|
412
|
|
|
|
|
|
|
=for usage |
413
|
|
|
|
|
|
|
|
414
|
|
|
|
|
|
|
my $windata = $win->get_samples |
415
|
|
|
|
|
|
|
|
416
|
|
|
|
|
|
|
=for ref |
417
|
|
|
|
|
|
|
|
418
|
|
|
|
|
|
|
Return a reference to the pdl of samples for the Window instance C<$win>. |
419
|
|
|
|
|
|
|
The samples will be generated with the method L if and only if |
420
|
|
|
|
|
|
|
they have not yet been generated. |
421
|
|
|
|
|
|
|
|
422
|
|
|
|
|
|
|
=cut |
423
|
|
|
|
|
|
|
|
424
|
|
|
|
|
|
|
sub get_samples { |
425
|
0
|
|
|
0
|
1
|
0
|
my $self = shift; |
426
|
0
|
0
|
|
|
|
0
|
$self->{samples} ? $self->{samples} : $self->samples; |
427
|
|
|
|
|
|
|
} |
428
|
|
|
|
|
|
|
|
429
|
|
|
|
|
|
|
=head2 get_modfreqs |
430
|
|
|
|
|
|
|
|
431
|
|
|
|
|
|
|
=for usage |
432
|
|
|
|
|
|
|
|
433
|
|
|
|
|
|
|
my $winfreqs = $win->get_modfreqs; |
434
|
|
|
|
|
|
|
my $winfreqs = $win->get_modfreqs({OPTS}); |
435
|
|
|
|
|
|
|
|
436
|
|
|
|
|
|
|
=for ref |
437
|
|
|
|
|
|
|
|
438
|
|
|
|
|
|
|
Return a reference to the pdl of the frequency response (modulus of the DFT) |
439
|
|
|
|
|
|
|
for the Window instance C<$win>. |
440
|
|
|
|
|
|
|
|
441
|
|
|
|
|
|
|
Options are passed to the method L. |
442
|
|
|
|
|
|
|
The data are created with L |
443
|
|
|
|
|
|
|
if they don't exist. The data are also created even |
444
|
|
|
|
|
|
|
if they already exist if options are supplied. Otherwise |
445
|
|
|
|
|
|
|
the cached data are returned. |
446
|
|
|
|
|
|
|
|
447
|
|
|
|
|
|
|
=head3 options |
448
|
|
|
|
|
|
|
|
449
|
|
|
|
|
|
|
=over |
450
|
|
|
|
|
|
|
|
451
|
|
|
|
|
|
|
=item min_bins => MIN |
452
|
|
|
|
|
|
|
|
453
|
|
|
|
|
|
|
This sets the minimum number of frequency bins. See |
454
|
|
|
|
|
|
|
L. Default 1000. |
455
|
|
|
|
|
|
|
|
456
|
|
|
|
|
|
|
=back |
457
|
|
|
|
|
|
|
|
458
|
|
|
|
|
|
|
=cut |
459
|
|
|
|
|
|
|
|
460
|
|
|
|
|
|
|
sub get_modfreqs { |
461
|
0
|
|
|
0
|
1
|
0
|
my $self = shift; |
462
|
0
|
|
|
|
|
0
|
my ($in_opts) = @_; |
463
|
0
|
0
|
0
|
|
|
0
|
($self->{modfreqs} and not $in_opts) ? $self->{modfreqs} : $self->modfreqs($in_opts); |
464
|
|
|
|
|
|
|
} |
465
|
|
|
|
|
|
|
|
466
|
|
|
|
|
|
|
=head2 get_params |
467
|
|
|
|
|
|
|
|
468
|
|
|
|
|
|
|
=for usage |
469
|
|
|
|
|
|
|
|
470
|
|
|
|
|
|
|
my $params = $win->get_params |
471
|
|
|
|
|
|
|
|
472
|
|
|
|
|
|
|
=for ref |
473
|
|
|
|
|
|
|
|
474
|
|
|
|
|
|
|
Create a new array containing the parameter values for the instance C<$win> |
475
|
|
|
|
|
|
|
and return a reference to the array. |
476
|
|
|
|
|
|
|
Note that not all window types take parameters. |
477
|
|
|
|
|
|
|
|
478
|
|
|
|
|
|
|
=cut |
479
|
|
|
|
|
|
|
|
480
|
|
|
|
|
|
|
sub get_params { |
481
|
0
|
|
|
0
|
1
|
0
|
my $self = shift; |
482
|
0
|
|
|
|
|
0
|
$self->{params}; |
483
|
|
|
|
|
|
|
} |
484
|
|
|
|
|
|
|
|
485
|
|
|
|
|
|
|
sub get_N { |
486
|
0
|
|
|
0
|
0
|
0
|
my $self = shift; |
487
|
0
|
|
|
|
|
0
|
$self->{N}; |
488
|
|
|
|
|
|
|
} |
489
|
|
|
|
|
|
|
|
490
|
|
|
|
|
|
|
=head2 get_name |
491
|
|
|
|
|
|
|
|
492
|
|
|
|
|
|
|
=for usage |
493
|
|
|
|
|
|
|
|
494
|
|
|
|
|
|
|
print $win->get_name , "\n"; |
495
|
|
|
|
|
|
|
|
496
|
|
|
|
|
|
|
=for ref |
497
|
|
|
|
|
|
|
|
498
|
|
|
|
|
|
|
Return a name suitable for printing associated with the window $win. This is |
499
|
|
|
|
|
|
|
something like the name used in the documentation for the particular |
500
|
|
|
|
|
|
|
window function. This is static data and does not depend on the instance. |
501
|
|
|
|
|
|
|
|
502
|
|
|
|
|
|
|
=cut |
503
|
|
|
|
|
|
|
|
504
|
|
|
|
|
|
|
sub get_name { |
505
|
0
|
|
|
0
|
1
|
0
|
my $self = shift; |
506
|
0
|
|
|
|
|
0
|
my $wd = $window_definitions{$self->{name}}; |
507
|
0
|
0
|
|
|
|
0
|
return $wd->{pfn} . ' window' if $wd->{pfn}; |
508
|
0
|
0
|
0
|
|
|
0
|
return $wd->{fn} . ' window' if $wd->{fn} and not $wd->{fn} =~ /^\*/; |
509
|
0
|
0
|
|
|
|
0
|
return $wd->{fn} if $wd->{fn}; |
510
|
0
|
|
|
|
|
0
|
return ucfirst($self->{name}) . ' window'; |
511
|
|
|
|
|
|
|
} |
512
|
|
|
|
|
|
|
|
513
|
|
|
|
|
|
|
sub get_param_names { |
514
|
0
|
|
|
0
|
0
|
0
|
my $self = shift; |
515
|
0
|
|
|
|
|
0
|
my $wd = $window_definitions{$self->{name}}; |
516
|
0
|
0
|
|
|
|
0
|
$wd->{params} ? ref($wd->{params}) ? $wd->{params} : [$wd->{params}] : undef; |
|
|
0
|
|
|
|
|
|
517
|
|
|
|
|
|
|
} |
518
|
|
|
|
|
|
|
|
519
|
|
|
|
|
|
|
sub format_param_vals { |
520
|
0
|
|
|
0
|
0
|
0
|
my $self = shift; |
521
|
0
|
|
|
|
|
0
|
my $p = $self->get('params'); |
522
|
0
|
0
|
|
|
|
0
|
return '' unless $p; |
523
|
0
|
|
|
|
|
0
|
my $names = $self->get_param_names; |
524
|
0
|
|
|
|
|
0
|
my @p = @$p; |
525
|
0
|
|
|
|
|
0
|
my @names = @$names; |
526
|
0
|
0
|
|
|
|
0
|
return '' unless $names; |
527
|
0
|
|
|
|
|
0
|
my @s; |
528
|
0
|
|
|
|
|
0
|
map { s/^\$// } @names; |
|
0
|
|
|
|
|
0
|
|
529
|
0
|
|
|
|
|
0
|
foreach (@p) { |
530
|
0
|
|
|
|
|
0
|
push @s, (shift @names) . ' = ' . $_; |
531
|
|
|
|
|
|
|
} |
532
|
0
|
|
|
|
|
0
|
join(', ', @s); |
533
|
|
|
|
|
|
|
} |
534
|
|
|
|
|
|
|
|
535
|
|
|
|
|
|
|
sub format_plot_param_vals { |
536
|
0
|
|
|
0
|
0
|
0
|
my $self = shift; |
537
|
0
|
|
|
|
|
0
|
my $ps = $self->format_param_vals; |
538
|
0
|
0
|
|
|
|
0
|
return '' unless $ps; |
539
|
0
|
|
|
|
|
0
|
': ' . $ps; |
540
|
|
|
|
|
|
|
} |
541
|
|
|
|
|
|
|
|
542
|
|
|
|
|
|
|
=head2 plot |
543
|
|
|
|
|
|
|
|
544
|
|
|
|
|
|
|
=for usage |
545
|
|
|
|
|
|
|
|
546
|
|
|
|
|
|
|
$win->plot; |
547
|
|
|
|
|
|
|
|
548
|
|
|
|
|
|
|
=for ref |
549
|
|
|
|
|
|
|
|
550
|
|
|
|
|
|
|
Plot the samples. Currently, only PDL::Graphics::Gnuplot is supported. |
551
|
|
|
|
|
|
|
The default display type is used. |
552
|
|
|
|
|
|
|
|
553
|
|
|
|
|
|
|
=cut |
554
|
|
|
|
|
|
|
|
555
|
|
|
|
|
|
|
sub plot { |
556
|
0
|
|
|
0
|
1
|
0
|
my $self = shift; |
557
|
0
|
0
|
|
|
|
0
|
barf "PDL::DSP::Windows::plot Gnuplot not available!" unless $HAVE_GNUPLOT; |
558
|
0
|
|
|
|
|
0
|
my $w = $self->get('samples'); |
559
|
0
|
|
|
|
|
0
|
my $title = $self->get_name() .$self->format_plot_param_vals; |
560
|
0
|
|
|
|
|
0
|
PDL::Graphics::Gnuplot::plot( title => $title, xlabel => 'Time (samples)', |
561
|
|
|
|
|
|
|
ylabel => 'amplitude', $w ); |
562
|
0
|
|
|
|
|
0
|
return $self; |
563
|
|
|
|
|
|
|
} |
564
|
|
|
|
|
|
|
|
565
|
|
|
|
|
|
|
=head2 plot_freq |
566
|
|
|
|
|
|
|
|
567
|
|
|
|
|
|
|
=for usage |
568
|
|
|
|
|
|
|
|
569
|
|
|
|
|
|
|
Can be called like this |
570
|
|
|
|
|
|
|
|
571
|
|
|
|
|
|
|
$win->plot_freq; |
572
|
|
|
|
|
|
|
|
573
|
|
|
|
|
|
|
|
574
|
|
|
|
|
|
|
Or this |
575
|
|
|
|
|
|
|
|
576
|
|
|
|
|
|
|
$win->plot_freq( {ordinate => ORDINATE }); |
577
|
|
|
|
|
|
|
|
578
|
|
|
|
|
|
|
|
579
|
|
|
|
|
|
|
=for ref |
580
|
|
|
|
|
|
|
|
581
|
|
|
|
|
|
|
Plot the frequency response (magnitude of the DFT of the window samples). |
582
|
|
|
|
|
|
|
The response is plotted in dB, and the frequency |
583
|
|
|
|
|
|
|
(by default) as a fraction of the Nyquist frequency. |
584
|
|
|
|
|
|
|
Currently, only PDL::Graphics::Gnuplot is supported. |
585
|
|
|
|
|
|
|
The default display type is used. |
586
|
|
|
|
|
|
|
|
587
|
|
|
|
|
|
|
=head3 options |
588
|
|
|
|
|
|
|
|
589
|
|
|
|
|
|
|
=over |
590
|
|
|
|
|
|
|
|
591
|
|
|
|
|
|
|
=item coord => COORD |
592
|
|
|
|
|
|
|
|
593
|
|
|
|
|
|
|
This sets the units of frequency of the co-ordinate axis. |
594
|
|
|
|
|
|
|
C must be one of C, for |
595
|
|
|
|
|
|
|
fraction of the nyquist frequency (range C<-1,1>), |
596
|
|
|
|
|
|
|
C, for fraction of the sampling frequncy (range |
597
|
|
|
|
|
|
|
C<-.5,.5>), or C for frequency bin number (range |
598
|
|
|
|
|
|
|
C<0,$N-1>). The default value is C. |
599
|
|
|
|
|
|
|
|
600
|
|
|
|
|
|
|
=item min_bins => MIN |
601
|
|
|
|
|
|
|
|
602
|
|
|
|
|
|
|
This sets the minimum number of frequency bins. See |
603
|
|
|
|
|
|
|
L. Default 1000. |
604
|
|
|
|
|
|
|
|
605
|
|
|
|
|
|
|
=back |
606
|
|
|
|
|
|
|
|
607
|
|
|
|
|
|
|
=cut |
608
|
|
|
|
|
|
|
|
609
|
|
|
|
|
|
|
sub plot_freq { |
610
|
0
|
|
|
0
|
1
|
0
|
my $self = shift; |
611
|
0
|
|
|
|
|
0
|
my $opt = new PDL::Options( |
612
|
|
|
|
|
|
|
{ |
613
|
|
|
|
|
|
|
coord => 'nyquist', |
614
|
|
|
|
|
|
|
min_bins => 1000 |
615
|
|
|
|
|
|
|
}); |
616
|
0
|
0
|
|
|
|
0
|
my $iopts = @_ ? shift : {}; |
617
|
0
|
|
|
|
|
0
|
my $opts = $opt->options($iopts); |
618
|
0
|
0
|
|
|
|
0
|
barf "PDL::DSP::Windows::plot Gnuplot not available!" unless $HAVE_GNUPLOT; |
619
|
0
|
|
|
|
|
0
|
my $mf = $self->get_modfreqs({ min_bins => $opts->{min_bins}}); |
620
|
0
|
|
|
|
|
0
|
$mf /= $mf->max; |
621
|
0
|
|
|
|
|
0
|
my $param_str = $self->format_plot_param_vals; |
622
|
0
|
|
|
|
|
0
|
my $title = $self->get_name() . $param_str |
623
|
|
|
|
|
|
|
. ', frequency response. ENBW=' . sprintf("%2.3f",$self->enbw); |
624
|
0
|
|
|
|
|
0
|
my $coord = $opts->{coord}; |
625
|
0
|
|
|
|
|
0
|
my ($coordinate_range,$xlab); |
626
|
0
|
0
|
|
|
|
0
|
if ($coord eq 'nyquist') { |
|
|
0
|
|
|
|
|
|
|
|
0
|
|
|
|
|
|
627
|
0
|
|
|
|
|
0
|
$coordinate_range = 1; |
628
|
0
|
|
|
|
|
0
|
$xlab = 'Fraction of Nyquist frequency'; |
629
|
|
|
|
|
|
|
} |
630
|
|
|
|
|
|
|
elsif ($coord eq 'sample') { |
631
|
0
|
|
|
|
|
0
|
$coordinate_range = .5; |
632
|
0
|
|
|
|
|
0
|
$xlab = 'Fraction of sampling freqeuncy'; |
633
|
|
|
|
|
|
|
} |
634
|
|
|
|
|
|
|
elsif ($coord eq 'bin') { |
635
|
0
|
|
|
|
|
0
|
$coordinate_range = ($self->get_N)/2; |
636
|
0
|
|
|
|
|
0
|
$xlab = 'bin'; |
637
|
|
|
|
|
|
|
} |
638
|
|
|
|
|
|
|
else { |
639
|
0
|
|
|
|
|
0
|
barf "plot_freq: Unknown ordinate unit specification $coord"; |
640
|
|
|
|
|
|
|
} |
641
|
0
|
|
|
|
|
0
|
my $coordinates = zeroes($mf)->xlinvals(-$coordinate_range,$coordinate_range); |
642
|
0
|
|
|
|
|
0
|
my $ylab = 'freqeuncy response (dB)'; |
643
|
0
|
|
|
|
|
0
|
PDL::Graphics::Gnuplot::plot(title => $title, |
644
|
|
|
|
|
|
|
xmin => -$coordinate_range, xmax => $coordinate_range, |
645
|
|
|
|
|
|
|
xlabel => $xlab, ylabel => $ylab, |
646
|
|
|
|
|
|
|
with => 'line', $coordinates, 20 * log10($mf) ); |
647
|
0
|
|
|
|
|
0
|
return $self; |
648
|
|
|
|
|
|
|
} |
649
|
|
|
|
|
|
|
|
650
|
|
|
|
|
|
|
=head2 enbw |
651
|
|
|
|
|
|
|
|
652
|
|
|
|
|
|
|
=for usage |
653
|
|
|
|
|
|
|
|
654
|
|
|
|
|
|
|
$win->enbw; |
655
|
|
|
|
|
|
|
|
656
|
|
|
|
|
|
|
=for ref |
657
|
|
|
|
|
|
|
|
658
|
|
|
|
|
|
|
Compute and return the equivalent noise bandwidth of the window. |
659
|
|
|
|
|
|
|
|
660
|
|
|
|
|
|
|
=cut |
661
|
|
|
|
|
|
|
|
662
|
|
|
|
|
|
|
sub enbw { |
663
|
14
|
|
|
14
|
1
|
19
|
my $self = shift; |
664
|
14
|
|
|
|
|
39
|
my $w = $self->get('samples'); # hmm have to quote samples here |
665
|
14
|
|
|
|
|
7878
|
($w->nelem) * ($w**2)->sum / ($w->sum)**2; |
666
|
|
|
|
|
|
|
} |
667
|
|
|
|
|
|
|
|
668
|
|
|
|
|
|
|
=head2 coherent_gain |
669
|
|
|
|
|
|
|
|
670
|
|
|
|
|
|
|
=for usage |
671
|
|
|
|
|
|
|
|
672
|
|
|
|
|
|
|
$win->coherent_gain; |
673
|
|
|
|
|
|
|
|
674
|
|
|
|
|
|
|
=for ref |
675
|
|
|
|
|
|
|
|
676
|
|
|
|
|
|
|
Compute and return the coherent gain (the dc gain) of the window. |
677
|
|
|
|
|
|
|
This is just the average of the samples. |
678
|
|
|
|
|
|
|
|
679
|
|
|
|
|
|
|
=cut |
680
|
|
|
|
|
|
|
|
681
|
|
|
|
|
|
|
sub coherent_gain { |
682
|
0
|
|
|
0
|
1
|
0
|
my $self = shift; |
683
|
0
|
|
|
|
|
0
|
my $w = $self->get('samples'); |
684
|
0
|
|
|
|
|
0
|
$w->sum / $w->nelem; |
685
|
|
|
|
|
|
|
} |
686
|
|
|
|
|
|
|
|
687
|
|
|
|
|
|
|
|
688
|
|
|
|
|
|
|
=head2 process_gain |
689
|
|
|
|
|
|
|
|
690
|
|
|
|
|
|
|
=for usage |
691
|
|
|
|
|
|
|
|
692
|
|
|
|
|
|
|
$win->coherent_gain; |
693
|
|
|
|
|
|
|
|
694
|
|
|
|
|
|
|
=for ref |
695
|
|
|
|
|
|
|
|
696
|
|
|
|
|
|
|
Compute and return the processing gain (the dc gain) of the window. |
697
|
|
|
|
|
|
|
This is just the multiplicative inverse of the C. |
698
|
|
|
|
|
|
|
|
699
|
|
|
|
|
|
|
=cut |
700
|
|
|
|
|
|
|
|
701
|
|
|
|
|
|
|
sub process_gain { |
702
|
0
|
|
|
0
|
1
|
0
|
my $self = shift; |
703
|
0
|
|
|
|
|
0
|
1/$self->enbw(); |
704
|
|
|
|
|
|
|
} |
705
|
|
|
|
|
|
|
|
706
|
|
|
|
|
|
|
# not quite correct for some reason. |
707
|
|
|
|
|
|
|
# Seems like 10*log10(this) / 1.154 |
708
|
|
|
|
|
|
|
# gives the correct answer in decibels |
709
|
|
|
|
|
|
|
|
710
|
|
|
|
|
|
|
=head2 scallop_loss |
711
|
|
|
|
|
|
|
|
712
|
|
|
|
|
|
|
=for usage |
713
|
|
|
|
|
|
|
|
714
|
|
|
|
|
|
|
$win->scallop_loss; |
715
|
|
|
|
|
|
|
|
716
|
|
|
|
|
|
|
=for ref |
717
|
|
|
|
|
|
|
|
718
|
|
|
|
|
|
|
**BROKEN**. |
719
|
|
|
|
|
|
|
Compute and return the scalloping loss of the window. |
720
|
|
|
|
|
|
|
|
721
|
|
|
|
|
|
|
=cut |
722
|
|
|
|
|
|
|
|
723
|
|
|
|
|
|
|
sub scallop_loss { |
724
|
0
|
|
|
0
|
1
|
0
|
my ($w) = @_; |
725
|
|
|
|
|
|
|
# my $x = (sequence($w) - ($w->nelem/2)) * (PI/$w->nelem); |
726
|
0
|
|
|
|
|
0
|
my $x = sequence($w) * (PI/$w->nelem); |
727
|
0
|
|
|
|
|
0
|
sqrt( (($w*cos($x))->sum)**2 + (($w*sin($x))->sum)**2 ) / |
728
|
|
|
|
|
|
|
$w->sum; |
729
|
|
|
|
|
|
|
} |
730
|
|
|
|
|
|
|
|
731
|
|
|
|
|
|
|
=head1 WINDOW FUNCTIONS |
732
|
|
|
|
|
|
|
|
733
|
|
|
|
|
|
|
These window-function subroutines return a pdl of $N samples. For most |
734
|
|
|
|
|
|
|
windows, there are a symmetric and a periodic version. The |
735
|
|
|
|
|
|
|
symmetric versions are functions of $N points, uniformly |
736
|
|
|
|
|
|
|
spaced, and taking values from x_lo through x_hi. Here, a |
737
|
|
|
|
|
|
|
periodic function of C< $N > points is equivalent to its |
738
|
|
|
|
|
|
|
symmetric counterpart of C<$N+1> points, with the final |
739
|
|
|
|
|
|
|
point omitted. The name of a periodic window-function subroutine is the |
740
|
|
|
|
|
|
|
same as that for the corresponding symmetric function, except it |
741
|
|
|
|
|
|
|
has the suffix C<_per>. The descriptions below describe the |
742
|
|
|
|
|
|
|
symmetric version of each window. |
743
|
|
|
|
|
|
|
|
744
|
|
|
|
|
|
|
The term 'Blackman-Harris family' is meant to include the Hamming family |
745
|
|
|
|
|
|
|
and the Blackman family. These are functions of sums of cosines. |
746
|
|
|
|
|
|
|
|
747
|
|
|
|
|
|
|
Unless otherwise noted, the arguments in the cosines of all symmetric |
748
|
|
|
|
|
|
|
window functions are multiples of C<$N> numbers uniformly spaced |
749
|
|
|
|
|
|
|
from C<0> through C<2 pi>. |
750
|
|
|
|
|
|
|
|
751
|
|
|
|
|
|
|
=cut |
752
|
|
|
|
|
|
|
|
753
|
|
|
|
|
|
|
sub bartlett { |
754
|
2
|
50
|
|
2
|
1
|
8
|
barf "bartlett: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
755
|
2
|
|
|
|
|
3
|
my ($N) = @_; |
756
|
2
|
|
|
|
|
8
|
1 - abs (zeroes($N)->xlinvals(-1,1)); |
757
|
|
|
|
|
|
|
} |
758
|
|
|
|
|
|
|
$window_definitions{bartlett} = { |
759
|
|
|
|
|
|
|
alias => [ 'fejer'], |
760
|
|
|
|
|
|
|
}; |
761
|
|
|
|
|
|
|
$winsubs{bartlett} = \&bartlett; |
762
|
|
|
|
|
|
|
|
763
|
|
|
|
|
|
|
sub bartlett_per { |
764
|
2
|
50
|
|
2
|
0
|
17
|
barf "bartlett: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
765
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
766
|
2
|
|
|
|
|
9
|
1 - abs (zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N)); |
767
|
|
|
|
|
|
|
} |
768
|
|
|
|
|
|
|
$window_definitions{bartlett} = { |
769
|
|
|
|
|
|
|
alias => [ 'fejer'], |
770
|
|
|
|
|
|
|
}; |
771
|
|
|
|
|
|
|
$winpersubs{bartlett}= \&bartlett_per; |
772
|
|
|
|
|
|
|
|
773
|
|
|
|
|
|
|
sub bartlett_hann { |
774
|
3
|
50
|
|
3
|
1
|
10
|
barf "bartlett_hann: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
775
|
3
|
|
|
|
|
6
|
my ($N) = @_; |
776
|
3
|
|
|
|
|
13
|
0.62 - 0.48 * abs (zeroes($N)->xlinvals(-0.5,0.5)) + 0.38* (cos(zeroes($N)->xlinvals(-(PI),PI))); |
777
|
|
|
|
|
|
|
} |
778
|
|
|
|
|
|
|
$window_definitions{bartlett_hann} = { |
779
|
|
|
|
|
|
|
fn => q!Bartlett-Hann!, |
780
|
|
|
|
|
|
|
alias => [ 'Modified Bartlett-Hann'], |
781
|
|
|
|
|
|
|
}; |
782
|
|
|
|
|
|
|
$winsubs{bartlett_hann} = \&bartlett_hann; |
783
|
|
|
|
|
|
|
|
784
|
|
|
|
|
|
|
sub bartlett_hann_per { |
785
|
2
|
50
|
|
2
|
0
|
10
|
barf "bartlett_hann: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
786
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
787
|
2
|
|
|
|
|
8
|
0.62 - 0.48 * abs (zeroes($N)->xlinvals(-0.5, (-0.5+0.5*($N-1))/$N)) + 0.38* (cos(zeroes($N)->xlinvals(-(PI), (-(PI)+PI*($N-1))/$N))); |
788
|
|
|
|
|
|
|
} |
789
|
|
|
|
|
|
|
$window_definitions{bartlett_hann} = { |
790
|
|
|
|
|
|
|
fn => q!Bartlett-Hann!, |
791
|
|
|
|
|
|
|
alias => [ 'Modified Bartlett-Hann'], |
792
|
|
|
|
|
|
|
}; |
793
|
|
|
|
|
|
|
$winpersubs{bartlett_hann}= \&bartlett_hann_per; |
794
|
|
|
|
|
|
|
|
795
|
|
|
|
|
|
|
sub blackman { |
796
|
3
|
50
|
|
3
|
1
|
11
|
barf "blackman: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
797
|
3
|
|
|
|
|
6
|
my ($N) = @_; |
798
|
3
|
|
|
|
|
12
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
799
|
|
|
|
|
|
|
|
800
|
3
|
|
|
|
|
2036
|
(0.34) + ($cx * ((-0.5) + ($cx * (0.16)))); |
801
|
|
|
|
|
|
|
} |
802
|
|
|
|
|
|
|
$window_definitions{blackman} = { |
803
|
|
|
|
|
|
|
fn => q!'classic' Blackman!, |
804
|
|
|
|
|
|
|
}; |
805
|
|
|
|
|
|
|
$winsubs{blackman} = \&blackman; |
806
|
|
|
|
|
|
|
|
807
|
|
|
|
|
|
|
sub blackman_per { |
808
|
2
|
50
|
|
2
|
0
|
8
|
barf "blackman: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
809
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
810
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
811
|
|
|
|
|
|
|
|
812
|
2
|
|
|
|
|
333
|
(0.34) + ($cx * ((-0.5) + ($cx * (0.16)))); |
813
|
|
|
|
|
|
|
} |
814
|
|
|
|
|
|
|
$window_definitions{blackman} = { |
815
|
|
|
|
|
|
|
fn => q!'classic' Blackman!, |
816
|
|
|
|
|
|
|
}; |
817
|
|
|
|
|
|
|
$winpersubs{blackman}= \&blackman_per; |
818
|
|
|
|
|
|
|
|
819
|
|
|
|
|
|
|
sub blackman_bnh { |
820
|
2
|
50
|
|
2
|
1
|
9
|
barf "blackman_bnh: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
821
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
822
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
823
|
|
|
|
|
|
|
|
824
|
2
|
|
|
|
|
319
|
(0.3461008) + ($cx * ((-0.4973406) + ($cx * (0.1565586)))); |
825
|
|
|
|
|
|
|
} |
826
|
|
|
|
|
|
|
$window_definitions{blackman_bnh} = { |
827
|
|
|
|
|
|
|
pfn => q!Blackman-Harris (bnh)!, |
828
|
|
|
|
|
|
|
fn => q!*An improved version of the 3-term Blackman-Harris window given by Nuttall (Ref 2, p. 89).!, |
829
|
|
|
|
|
|
|
}; |
830
|
|
|
|
|
|
|
$winsubs{blackman_bnh} = \&blackman_bnh; |
831
|
|
|
|
|
|
|
|
832
|
|
|
|
|
|
|
sub blackman_bnh_per { |
833
|
2
|
50
|
|
2
|
0
|
9
|
barf "blackman_bnh: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
834
|
2
|
|
|
|
|
6
|
my ($N) = @_; |
835
|
2
|
|
|
|
|
8
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
836
|
|
|
|
|
|
|
|
837
|
2
|
|
|
|
|
423
|
(0.3461008) + ($cx * ((-0.4973406) + ($cx * (0.1565586)))); |
838
|
|
|
|
|
|
|
} |
839
|
|
|
|
|
|
|
$window_definitions{blackman_bnh} = { |
840
|
|
|
|
|
|
|
pfn => q!Blackman-Harris (bnh)!, |
841
|
|
|
|
|
|
|
fn => q!*An improved version of the 3-term Blackman-Harris window given by Nuttall (Ref 2, p. 89).!, |
842
|
|
|
|
|
|
|
}; |
843
|
|
|
|
|
|
|
$winpersubs{blackman_bnh}= \&blackman_bnh_per; |
844
|
|
|
|
|
|
|
|
845
|
|
|
|
|
|
|
sub blackman_ex { |
846
|
2
|
50
|
|
2
|
1
|
7
|
barf "blackman_ex: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
847
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
848
|
2
|
|
|
|
|
9
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
849
|
|
|
|
|
|
|
|
850
|
2
|
|
|
|
|
318
|
(0.349742046431642) + ($cx * ((-0.496560619088564) + ($cx * (0.153697334479794)))); |
851
|
|
|
|
|
|
|
} |
852
|
|
|
|
|
|
|
$window_definitions{blackman_ex} = { |
853
|
|
|
|
|
|
|
fn => q!'exact' Blackman!, |
854
|
|
|
|
|
|
|
}; |
855
|
|
|
|
|
|
|
$winsubs{blackman_ex} = \&blackman_ex; |
856
|
|
|
|
|
|
|
|
857
|
|
|
|
|
|
|
sub blackman_ex_per { |
858
|
2
|
50
|
|
2
|
0
|
8
|
barf "blackman_ex: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
859
|
2
|
|
|
|
|
6
|
my ($N) = @_; |
860
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
861
|
|
|
|
|
|
|
|
862
|
2
|
|
|
|
|
309
|
(0.349742046431642) + ($cx * ((-0.496560619088564) + ($cx * (0.153697334479794)))); |
863
|
|
|
|
|
|
|
} |
864
|
|
|
|
|
|
|
$window_definitions{blackman_ex} = { |
865
|
|
|
|
|
|
|
fn => q!'exact' Blackman!, |
866
|
|
|
|
|
|
|
}; |
867
|
|
|
|
|
|
|
$winpersubs{blackman_ex}= \&blackman_ex_per; |
868
|
|
|
|
|
|
|
|
869
|
|
|
|
|
|
|
sub blackman_gen { |
870
|
2
|
50
|
|
2
|
1
|
8
|
barf "blackman_gen: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
871
|
2
|
|
|
|
|
5
|
my ($N,$alpha) = @_; |
872
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
873
|
|
|
|
|
|
|
|
874
|
2
|
|
|
|
|
317
|
(.5 - $alpha) + ($cx * ((-.5) + ($cx * ($alpha)))); |
875
|
|
|
|
|
|
|
} |
876
|
|
|
|
|
|
|
$window_definitions{blackman_gen} = { |
877
|
|
|
|
|
|
|
pfn => q!General classic Blackman!, |
878
|
|
|
|
|
|
|
fn => q!*A single parameter family of the 3-term Blackman window. !, |
879
|
|
|
|
|
|
|
params => [ '$alpha'], |
880
|
|
|
|
|
|
|
}; |
881
|
|
|
|
|
|
|
$winsubs{blackman_gen} = \&blackman_gen; |
882
|
|
|
|
|
|
|
|
883
|
|
|
|
|
|
|
sub blackman_gen_per { |
884
|
2
|
50
|
|
2
|
0
|
10
|
barf "blackman_gen: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
885
|
2
|
|
|
|
|
6
|
my ($N,$alpha) = @_; |
886
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
887
|
|
|
|
|
|
|
|
888
|
2
|
|
|
|
|
323
|
(.5 - $alpha) + ($cx * ((-.5) + ($cx * ($alpha)))); |
889
|
|
|
|
|
|
|
} |
890
|
|
|
|
|
|
|
$window_definitions{blackman_gen} = { |
891
|
|
|
|
|
|
|
pfn => q!General classic Blackman!, |
892
|
|
|
|
|
|
|
fn => q!*A single parameter family of the 3-term Blackman window. !, |
893
|
|
|
|
|
|
|
params => [ '$alpha'], |
894
|
|
|
|
|
|
|
}; |
895
|
|
|
|
|
|
|
$winpersubs{blackman_gen}= \&blackman_gen_per; |
896
|
|
|
|
|
|
|
|
897
|
|
|
|
|
|
|
sub blackman_gen3 { |
898
|
2
|
50
|
|
2
|
1
|
8
|
barf "blackman_gen3: 4 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 4; |
899
|
2
|
|
|
|
|
5
|
my ($N,$a0,$a1,$a2) = @_; |
900
|
2
|
|
|
|
|
10
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
901
|
|
|
|
|
|
|
|
902
|
2
|
|
|
|
|
363
|
($a0 - $a2) + ($cx * ((-$a1) + ($cx * (2*$a2)))); |
903
|
|
|
|
|
|
|
} |
904
|
|
|
|
|
|
|
$window_definitions{blackman_gen3} = { |
905
|
|
|
|
|
|
|
fn => q!*The general form of the Blackman family. !, |
906
|
|
|
|
|
|
|
params => [ '$a0','$a1','$a2'], |
907
|
|
|
|
|
|
|
}; |
908
|
|
|
|
|
|
|
$winsubs{blackman_gen3} = \&blackman_gen3; |
909
|
|
|
|
|
|
|
|
910
|
|
|
|
|
|
|
sub blackman_gen3_per { |
911
|
2
|
50
|
|
2
|
0
|
9
|
barf "blackman_gen3: 4 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 4; |
912
|
2
|
|
|
|
|
5
|
my ($N,$a0,$a1,$a2) = @_; |
913
|
2
|
|
|
|
|
10
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
914
|
|
|
|
|
|
|
|
915
|
2
|
|
|
|
|
375
|
($a0 - $a2) + ($cx * ((-$a1) + ($cx * (2*$a2)))); |
916
|
|
|
|
|
|
|
} |
917
|
|
|
|
|
|
|
$window_definitions{blackman_gen3} = { |
918
|
|
|
|
|
|
|
fn => q!*The general form of the Blackman family. !, |
919
|
|
|
|
|
|
|
params => [ '$a0','$a1','$a2'], |
920
|
|
|
|
|
|
|
}; |
921
|
|
|
|
|
|
|
$winpersubs{blackman_gen3}= \&blackman_gen3_per; |
922
|
|
|
|
|
|
|
|
923
|
|
|
|
|
|
|
sub blackman_gen4 { |
924
|
2
|
50
|
|
2
|
1
|
8
|
barf "blackman_gen4: 5 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 5; |
925
|
2
|
|
|
|
|
5
|
my ($N,$a0,$a1,$a2,$a3) = @_; |
926
|
2
|
|
|
|
|
9
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
927
|
|
|
|
|
|
|
|
928
|
2
|
|
|
|
|
323
|
($a0 - $a2) + ($cx * ((-$a1 + 3 * $a3) + ($cx * (2*$a2 + $cx * (-4*$a3) )))); |
929
|
|
|
|
|
|
|
} |
930
|
|
|
|
|
|
|
$window_definitions{blackman_gen4} = { |
931
|
|
|
|
|
|
|
fn => q!*The general 4-term Blackman-Harris window. !, |
932
|
|
|
|
|
|
|
params => [ '$a0','$a1','$a2','$a3'], |
933
|
|
|
|
|
|
|
}; |
934
|
|
|
|
|
|
|
$winsubs{blackman_gen4} = \&blackman_gen4; |
935
|
|
|
|
|
|
|
|
936
|
|
|
|
|
|
|
sub blackman_gen4_per { |
937
|
2
|
50
|
|
2
|
0
|
10
|
barf "blackman_gen4: 5 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 5; |
938
|
2
|
|
|
|
|
6
|
my ($N,$a0,$a1,$a2,$a3) = @_; |
939
|
2
|
|
|
|
|
8
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
940
|
|
|
|
|
|
|
|
941
|
2
|
|
|
|
|
334
|
($a0 - $a2) + ($cx * ((-$a1 + 3 * $a3) + ($cx * (2*$a2 + $cx * (-4*$a3) )))); |
942
|
|
|
|
|
|
|
} |
943
|
|
|
|
|
|
|
$window_definitions{blackman_gen4} = { |
944
|
|
|
|
|
|
|
fn => q!*The general 4-term Blackman-Harris window. !, |
945
|
|
|
|
|
|
|
params => [ '$a0','$a1','$a2','$a3'], |
946
|
|
|
|
|
|
|
}; |
947
|
|
|
|
|
|
|
$winpersubs{blackman_gen4}= \&blackman_gen4_per; |
948
|
|
|
|
|
|
|
|
949
|
|
|
|
|
|
|
sub blackman_gen5 { |
950
|
2
|
50
|
|
2
|
1
|
8
|
barf "blackman_gen5: 6 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 6; |
951
|
2
|
|
|
|
|
4
|
my ($N,$a0,$a1,$a2,$a3,$a4) = @_; |
952
|
2
|
|
|
|
|
9
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
953
|
|
|
|
|
|
|
|
954
|
2
|
|
|
|
|
325
|
($a0 - $a2 + $a4) + ($cx * ((-$a1 + 3 * $a3) + ($cx * (2*$a2 -8*$a4 + $cx * (-4*$a3 +$cx *(8*$a4)) )))); |
955
|
|
|
|
|
|
|
} |
956
|
|
|
|
|
|
|
$window_definitions{blackman_gen5} = { |
957
|
|
|
|
|
|
|
fn => q!*The general 5-term Blackman-Harris window. !, |
958
|
|
|
|
|
|
|
params => [ '$a0','$a1','$a2','$a3','$a4'], |
959
|
|
|
|
|
|
|
}; |
960
|
|
|
|
|
|
|
$winsubs{blackman_gen5} = \&blackman_gen5; |
961
|
|
|
|
|
|
|
|
962
|
|
|
|
|
|
|
sub blackman_gen5_per { |
963
|
2
|
50
|
|
2
|
0
|
9
|
barf "blackman_gen5: 6 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 6; |
964
|
2
|
|
|
|
|
7
|
my ($N,$a0,$a1,$a2,$a3,$a4) = @_; |
965
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
966
|
|
|
|
|
|
|
|
967
|
2
|
|
|
|
|
309
|
($a0 - $a2 + $a4) + ($cx * ((-$a1 + 3 * $a3) + ($cx * (2*$a2 -8*$a4 + $cx * (-4*$a3 +$cx *(8*$a4)) )))); |
968
|
|
|
|
|
|
|
} |
969
|
|
|
|
|
|
|
$window_definitions{blackman_gen5} = { |
970
|
|
|
|
|
|
|
fn => q!*The general 5-term Blackman-Harris window. !, |
971
|
|
|
|
|
|
|
params => [ '$a0','$a1','$a2','$a3','$a4'], |
972
|
|
|
|
|
|
|
}; |
973
|
|
|
|
|
|
|
$winpersubs{blackman_gen5}= \&blackman_gen5_per; |
974
|
|
|
|
|
|
|
|
975
|
|
|
|
|
|
|
sub blackman_harris { |
976
|
2
|
50
|
|
2
|
1
|
8
|
barf "blackman_harris: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
977
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
978
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
979
|
|
|
|
|
|
|
|
980
|
2
|
|
|
|
|
342
|
(0.343103) + ($cx * ((-0.49755) + ($cx * (0.15844)))); |
981
|
|
|
|
|
|
|
} |
982
|
|
|
|
|
|
|
$window_definitions{blackman_harris} = { |
983
|
|
|
|
|
|
|
fn => q!Blackman-Harris!, |
984
|
|
|
|
|
|
|
alias => [ 'Minimum three term (sample) Blackman-Harris'], |
985
|
|
|
|
|
|
|
}; |
986
|
|
|
|
|
|
|
$winsubs{blackman_harris} = \&blackman_harris; |
987
|
|
|
|
|
|
|
|
988
|
|
|
|
|
|
|
sub blackman_harris_per { |
989
|
2
|
50
|
|
2
|
0
|
9
|
barf "blackman_harris: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
990
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
991
|
2
|
|
|
|
|
10
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
992
|
|
|
|
|
|
|
|
993
|
2
|
|
|
|
|
363
|
(0.343103) + ($cx * ((-0.49755) + ($cx * (0.15844)))); |
994
|
|
|
|
|
|
|
} |
995
|
|
|
|
|
|
|
$window_definitions{blackman_harris} = { |
996
|
|
|
|
|
|
|
fn => q!Blackman-Harris!, |
997
|
|
|
|
|
|
|
alias => [ 'Minimum three term (sample) Blackman-Harris'], |
998
|
|
|
|
|
|
|
}; |
999
|
|
|
|
|
|
|
$winpersubs{blackman_harris}= \&blackman_harris_per; |
1000
|
|
|
|
|
|
|
|
1001
|
|
|
|
|
|
|
sub blackman_harris4 { |
1002
|
4
|
50
|
|
4
|
1
|
14
|
barf "blackman_harris4: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1003
|
4
|
|
|
|
|
8
|
my ($N) = @_; |
1004
|
4
|
|
|
|
|
16
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
1005
|
|
|
|
|
|
|
|
1006
|
4
|
|
|
|
|
1969
|
(0.21747) + ($cx * ((-0.45325) + ($cx * (0.28256 + $cx * (-0.04672) )))); |
1007
|
|
|
|
|
|
|
} |
1008
|
|
|
|
|
|
|
$window_definitions{blackman_harris4} = { |
1009
|
|
|
|
|
|
|
fn => q!minimum (sidelobe) four term Blackman-Harris!, |
1010
|
|
|
|
|
|
|
alias => [ 'Blackman-Harris'], |
1011
|
|
|
|
|
|
|
}; |
1012
|
|
|
|
|
|
|
$winsubs{blackman_harris4} = \&blackman_harris4; |
1013
|
|
|
|
|
|
|
|
1014
|
|
|
|
|
|
|
sub blackman_harris4_per { |
1015
|
2
|
50
|
|
2
|
0
|
8
|
barf "blackman_harris4: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1016
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1017
|
2
|
|
|
|
|
9
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1018
|
|
|
|
|
|
|
|
1019
|
2
|
|
|
|
|
311
|
(0.21747) + ($cx * ((-0.45325) + ($cx * (0.28256 + $cx * (-0.04672) )))); |
1020
|
|
|
|
|
|
|
} |
1021
|
|
|
|
|
|
|
$window_definitions{blackman_harris4} = { |
1022
|
|
|
|
|
|
|
fn => q!minimum (sidelobe) four term Blackman-Harris!, |
1023
|
|
|
|
|
|
|
alias => [ 'Blackman-Harris'], |
1024
|
|
|
|
|
|
|
}; |
1025
|
|
|
|
|
|
|
$winpersubs{blackman_harris4}= \&blackman_harris4_per; |
1026
|
|
|
|
|
|
|
|
1027
|
|
|
|
|
|
|
sub blackman_nuttall { |
1028
|
3
|
50
|
|
3
|
1
|
10
|
barf "blackman_nuttall: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1029
|
3
|
|
|
|
|
5
|
my ($N) = @_; |
1030
|
3
|
|
|
|
|
12
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
1031
|
|
|
|
|
|
|
|
1032
|
3
|
|
|
|
|
475
|
(0.2269824) + ($cx * ((-0.4572542) + ($cx * (0.273199 + $cx * (-0.0425644) )))); |
1033
|
|
|
|
|
|
|
} |
1034
|
|
|
|
|
|
|
$window_definitions{blackman_nuttall} = { |
1035
|
|
|
|
|
|
|
fn => q!Blackman-Nuttall!, |
1036
|
|
|
|
|
|
|
}; |
1037
|
|
|
|
|
|
|
$winsubs{blackman_nuttall} = \&blackman_nuttall; |
1038
|
|
|
|
|
|
|
|
1039
|
|
|
|
|
|
|
sub blackman_nuttall_per { |
1040
|
2
|
50
|
|
2
|
0
|
8
|
barf "blackman_nuttall: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1041
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
1042
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1043
|
|
|
|
|
|
|
|
1044
|
2
|
|
|
|
|
309
|
(0.2269824) + ($cx * ((-0.4572542) + ($cx * (0.273199 + $cx * (-0.0425644) )))); |
1045
|
|
|
|
|
|
|
} |
1046
|
|
|
|
|
|
|
$window_definitions{blackman_nuttall} = { |
1047
|
|
|
|
|
|
|
fn => q!Blackman-Nuttall!, |
1048
|
|
|
|
|
|
|
}; |
1049
|
|
|
|
|
|
|
$winpersubs{blackman_nuttall}= \&blackman_nuttall_per; |
1050
|
|
|
|
|
|
|
|
1051
|
|
|
|
|
|
|
sub bohman { |
1052
|
4
|
50
|
|
4
|
1
|
14
|
barf "bohman: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1053
|
4
|
|
|
|
|
8
|
my ($N) = @_; |
1054
|
4
|
|
|
|
|
16
|
my $x = abs((zeroes($N)->xlinvals(-1,1))); |
1055
|
4
|
|
|
|
|
909
|
(1-$x)*cos(PI*$x) +(1/PI)*sin(PI*$x); |
1056
|
|
|
|
|
|
|
} |
1057
|
|
|
|
|
|
|
$window_definitions{bohman} = { |
1058
|
|
|
|
|
|
|
}; |
1059
|
|
|
|
|
|
|
$winsubs{bohman} = \&bohman; |
1060
|
|
|
|
|
|
|
|
1061
|
|
|
|
|
|
|
sub bohman_per { |
1062
|
2
|
50
|
|
2
|
0
|
9
|
barf "bohman: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1063
|
2
|
|
|
|
|
3
|
my ($N) = @_; |
1064
|
2
|
|
|
|
|
7
|
my $x = abs((zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N))); |
1065
|
2
|
|
|
|
|
289
|
(1-$x)*cos(PI*$x) +(1/PI)*sin(PI*$x); |
1066
|
|
|
|
|
|
|
} |
1067
|
|
|
|
|
|
|
$window_definitions{bohman} = { |
1068
|
|
|
|
|
|
|
}; |
1069
|
|
|
|
|
|
|
$winpersubs{bohman}= \&bohman_per; |
1070
|
|
|
|
|
|
|
|
1071
|
|
|
|
|
|
|
sub cauchy { |
1072
|
3
|
50
|
|
3
|
1
|
11
|
barf "cauchy: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1073
|
3
|
|
|
|
|
7
|
my ($N,$alpha) = @_; |
1074
|
3
|
|
|
|
|
13
|
1 / (1 + ((zeroes($N)->xlinvals(-1,1)) * $alpha)**2); |
1075
|
|
|
|
|
|
|
} |
1076
|
|
|
|
|
|
|
$window_definitions{cauchy} = { |
1077
|
|
|
|
|
|
|
params => [ '$alpha'], |
1078
|
|
|
|
|
|
|
alias => [ 'Abel','Poisson'], |
1079
|
|
|
|
|
|
|
}; |
1080
|
|
|
|
|
|
|
$winsubs{cauchy} = \&cauchy; |
1081
|
|
|
|
|
|
|
|
1082
|
|
|
|
|
|
|
sub cauchy_per { |
1083
|
2
|
50
|
|
2
|
0
|
11
|
barf "cauchy: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1084
|
2
|
|
|
|
|
4
|
my ($N,$alpha) = @_; |
1085
|
2
|
|
|
|
|
10
|
1 / (1 + ((zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N)) * $alpha)**2); |
1086
|
|
|
|
|
|
|
} |
1087
|
|
|
|
|
|
|
$window_definitions{cauchy} = { |
1088
|
|
|
|
|
|
|
params => [ '$alpha'], |
1089
|
|
|
|
|
|
|
alias => [ 'Abel','Poisson'], |
1090
|
|
|
|
|
|
|
}; |
1091
|
|
|
|
|
|
|
$winpersubs{cauchy}= \&cauchy_per; |
1092
|
|
|
|
|
|
|
|
1093
|
|
|
|
|
|
|
sub chebyshev { |
1094
|
2
|
50
|
|
2
|
1
|
6
|
barf "chebyshev: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1095
|
2
|
|
|
|
|
4
|
my ($N,$at) = @_; |
1096
|
|
|
|
|
|
|
|
1097
|
2
|
|
|
|
|
3
|
my ($M,$M1,$pos,$pos1); |
1098
|
0
|
|
|
|
|
0
|
my $cw; |
1099
|
2
|
|
|
|
|
63
|
my $beta = cosh(1/($N-1) * acosh(1/(10**(-$at/20)))); |
1100
|
2
|
|
|
|
|
72
|
my $k = sequence($N); |
1101
|
2
|
|
|
|
|
189
|
my $x = $beta * cos(PI*$k/$N); |
1102
|
2
|
|
|
|
|
114
|
$cw = chebpoly($N-1,$x); |
1103
|
2
|
100
|
|
|
|
8
|
if ( $N % 2 ) { # odd |
1104
|
1
|
|
|
|
|
3
|
$M1 = ($N+1)/2; |
1105
|
1
|
|
|
|
|
3
|
$M = $M1 - 1; |
1106
|
1
|
|
|
|
|
2
|
$pos = 0; |
1107
|
1
|
|
|
|
|
1
|
$pos1 = 1; |
1108
|
1
|
|
|
|
|
5
|
PDL::FFT::realfft($cw); |
1109
|
|
|
|
|
|
|
} |
1110
|
|
|
|
|
|
|
else { # half-sample delay (even order) |
1111
|
1
|
|
|
|
|
4
|
my $arg = PI/$N * sequence($N); |
1112
|
1
|
|
|
|
|
102
|
my $cw_im = $cw * sin($arg); |
1113
|
1
|
|
|
|
|
25
|
$cw *= cos($arg); |
1114
|
1
|
|
|
|
|
27
|
PDL::FFT::fftnd($cw,$cw_im); |
1115
|
1
|
|
|
|
|
113
|
$M1 = $N/2; |
1116
|
1
|
|
|
|
|
2
|
$M = $M1-1; |
1117
|
1
|
|
|
|
|
2
|
$pos = 1; |
1118
|
1
|
|
|
|
|
7
|
$pos1 = 0; |
1119
|
|
|
|
|
|
|
} |
1120
|
2
|
|
|
|
|
133
|
$cw /= ($cw->at($pos)); |
1121
|
2
|
|
|
|
|
45
|
my $cwout = zeroes($N); |
1122
|
2
|
|
|
|
|
96
|
$cwout->slice("0:$M") .= $cw->slice("$M:0:-1"); |
1123
|
2
|
|
|
|
|
88
|
$cwout->slice("$M1:-1") .= $cw->slice("$pos1:$M"); |
1124
|
2
|
|
|
|
|
90
|
$cwout /= max($cwout); |
1125
|
2
|
|
|
|
|
148
|
$cwout; |
1126
|
|
|
|
|
|
|
; |
1127
|
|
|
|
|
|
|
} |
1128
|
|
|
|
|
|
|
$window_definitions{chebyshev} = { |
1129
|
|
|
|
|
|
|
params => [ '$at'], |
1130
|
|
|
|
|
|
|
alias => [ 'Dolph-Chebyshev'], |
1131
|
|
|
|
|
|
|
}; |
1132
|
|
|
|
|
|
|
$winsubs{chebyshev} = \&chebyshev; |
1133
|
|
|
|
|
|
|
|
1134
|
|
|
|
|
|
|
sub cos_alpha { |
1135
|
5
|
50
|
|
5
|
1
|
14
|
barf "cos_alpha: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1136
|
5
|
|
|
|
|
11
|
my ($N,$alpha) = @_; |
1137
|
5
|
|
|
|
|
17
|
(sin(zeroes($N)->xlinvals(0,PI)))**$alpha ; |
1138
|
|
|
|
|
|
|
} |
1139
|
|
|
|
|
|
|
$window_definitions{cos_alpha} = { |
1140
|
|
|
|
|
|
|
params => [ '$alpha'], |
1141
|
|
|
|
|
|
|
alias => [ 'Power-of-cosine'], |
1142
|
|
|
|
|
|
|
}; |
1143
|
|
|
|
|
|
|
$winsubs{cos_alpha} = \&cos_alpha; |
1144
|
|
|
|
|
|
|
|
1145
|
|
|
|
|
|
|
sub cos_alpha_per { |
1146
|
2
|
50
|
|
2
|
0
|
8
|
barf "cos_alpha: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1147
|
2
|
|
|
|
|
5
|
my ($N,$alpha) = @_; |
1148
|
2
|
|
|
|
|
9
|
(sin(zeroes($N)->xlinvals(0, PI*($N-1)/$N)))**$alpha ; |
1149
|
|
|
|
|
|
|
} |
1150
|
|
|
|
|
|
|
$window_definitions{cos_alpha} = { |
1151
|
|
|
|
|
|
|
params => [ '$alpha'], |
1152
|
|
|
|
|
|
|
alias => [ 'Power-of-cosine'], |
1153
|
|
|
|
|
|
|
}; |
1154
|
|
|
|
|
|
|
$winpersubs{cos_alpha}= \&cos_alpha_per; |
1155
|
|
|
|
|
|
|
|
1156
|
|
|
|
|
|
|
sub cosine { |
1157
|
3
|
50
|
|
3
|
1
|
10
|
barf "cosine: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1158
|
3
|
|
|
|
|
5
|
my ($N) = @_; |
1159
|
3
|
|
|
|
|
11
|
(sin(zeroes($N)->xlinvals(0,PI))); |
1160
|
|
|
|
|
|
|
} |
1161
|
|
|
|
|
|
|
$window_definitions{cosine} = { |
1162
|
|
|
|
|
|
|
alias => [ 'sine'], |
1163
|
|
|
|
|
|
|
}; |
1164
|
|
|
|
|
|
|
$winsubs{cosine} = \&cosine; |
1165
|
|
|
|
|
|
|
|
1166
|
|
|
|
|
|
|
sub cosine_per { |
1167
|
2
|
50
|
|
2
|
0
|
7
|
barf "cosine: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1168
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1169
|
2
|
|
|
|
|
8
|
(sin(zeroes($N)->xlinvals(0, PI*($N-1)/$N))); |
1170
|
|
|
|
|
|
|
} |
1171
|
|
|
|
|
|
|
$window_definitions{cosine} = { |
1172
|
|
|
|
|
|
|
alias => [ 'sine'], |
1173
|
|
|
|
|
|
|
}; |
1174
|
|
|
|
|
|
|
$winpersubs{cosine}= \&cosine_per; |
1175
|
|
|
|
|
|
|
|
1176
|
|
|
|
|
|
|
sub dpss { |
1177
|
0
|
0
|
|
0
|
1
|
0
|
barf "dpss: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1178
|
0
|
|
|
|
|
0
|
my ($N,$beta) = @_; |
1179
|
|
|
|
|
|
|
|
1180
|
0
|
0
|
|
|
|
0
|
barf 'dpss: PDL::LinearAlgebra not installed.' unless $HAVE_LinearAlgebra; |
1181
|
0
|
0
|
0
|
|
|
0
|
barf "dpss: $beta not between 0 and $N." unless |
1182
|
|
|
|
|
|
|
$beta >= 0 and $beta <= $N; |
1183
|
0
|
|
|
|
|
0
|
$beta /= ($N/2); |
1184
|
0
|
|
|
|
|
0
|
my $k = sequence($N); |
1185
|
0
|
|
|
|
|
0
|
my $s = sin(PI*$beta*$k)/$k; |
1186
|
0
|
|
|
|
|
0
|
$s->slice('0') .= $beta; |
1187
|
0
|
|
|
|
|
0
|
my ($ev,$e) = eigens_sym(PDL::LinearAlgebra::Special::mtoeplitz($s)); |
1188
|
0
|
|
|
|
|
0
|
my $i = $e->maximum_ind; |
1189
|
0
|
|
|
|
|
0
|
$ev->slice("($i)")->copy; |
1190
|
|
|
|
|
|
|
; |
1191
|
|
|
|
|
|
|
} |
1192
|
|
|
|
|
|
|
$window_definitions{dpss} = { |
1193
|
|
|
|
|
|
|
fn => q!Digital Prolate Spheroidal Sequence (DPSS)!, |
1194
|
|
|
|
|
|
|
params => [ '$beta'], |
1195
|
|
|
|
|
|
|
alias => [ 'sleppian'], |
1196
|
|
|
|
|
|
|
}; |
1197
|
|
|
|
|
|
|
$winsubs{dpss} = \&dpss; |
1198
|
|
|
|
|
|
|
|
1199
|
|
|
|
|
|
|
sub dpss_per { |
1200
|
0
|
0
|
|
0
|
0
|
0
|
barf "dpss: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1201
|
0
|
|
|
|
|
0
|
my ($N,$beta) = @_; |
1202
|
0
|
|
|
|
|
0
|
$N++; |
1203
|
|
|
|
|
|
|
|
1204
|
0
|
0
|
|
|
|
0
|
barf 'dpss: PDL::LinearAlgebra not installed.' unless $HAVE_LinearAlgebra; |
1205
|
0
|
0
|
0
|
|
|
0
|
barf "dpss: $beta not between 0 and $N." unless |
1206
|
|
|
|
|
|
|
$beta >= 0 and $beta <= $N; |
1207
|
0
|
|
|
|
|
0
|
$beta /= ($N/2); |
1208
|
0
|
|
|
|
|
0
|
my $k = sequence($N); |
1209
|
0
|
|
|
|
|
0
|
my $s = sin(PI*$beta*$k)/$k; |
1210
|
0
|
|
|
|
|
0
|
$s->slice('0') .= $beta; |
1211
|
0
|
|
|
|
|
0
|
my ($ev,$e) = eigens_sym(PDL::LinearAlgebra::Special::mtoeplitz($s)); |
1212
|
0
|
|
|
|
|
0
|
my $i = $e->maximum_ind; |
1213
|
0
|
|
|
|
|
0
|
$ev->slice("($i),0:-2")->copy; |
1214
|
|
|
|
|
|
|
; |
1215
|
|
|
|
|
|
|
} |
1216
|
|
|
|
|
|
|
$window_definitions{dpss} = { |
1217
|
|
|
|
|
|
|
fn => q!Digital Prolate Spheroidal Sequence (DPSS)!, |
1218
|
|
|
|
|
|
|
params => [ '$beta'], |
1219
|
|
|
|
|
|
|
alias => [ 'sleppian'], |
1220
|
|
|
|
|
|
|
}; |
1221
|
|
|
|
|
|
|
$winpersubs{dpss}= \&dpss_per; |
1222
|
|
|
|
|
|
|
|
1223
|
|
|
|
|
|
|
sub exponential { |
1224
|
2
|
50
|
|
2
|
1
|
8
|
barf "exponential: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1225
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
1226
|
2
|
|
|
|
|
7
|
2 ** (1 - abs (zeroes($N)->xlinvals(-1,1))) - 1; |
1227
|
|
|
|
|
|
|
} |
1228
|
|
|
|
|
|
|
$window_definitions{exponential} = { |
1229
|
|
|
|
|
|
|
}; |
1230
|
|
|
|
|
|
|
$winsubs{exponential} = \&exponential; |
1231
|
|
|
|
|
|
|
|
1232
|
|
|
|
|
|
|
sub exponential_per { |
1233
|
2
|
50
|
|
2
|
0
|
9
|
barf "exponential: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1234
|
2
|
|
|
|
|
3
|
my ($N) = @_; |
1235
|
2
|
|
|
|
|
9
|
2 ** (1 - abs (zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N))) - 1; |
1236
|
|
|
|
|
|
|
} |
1237
|
|
|
|
|
|
|
$window_definitions{exponential} = { |
1238
|
|
|
|
|
|
|
}; |
1239
|
|
|
|
|
|
|
$winpersubs{exponential}= \&exponential_per; |
1240
|
|
|
|
|
|
|
|
1241
|
|
|
|
|
|
|
sub flattop { |
1242
|
4
|
50
|
|
4
|
1
|
13
|
barf "flattop: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1243
|
4
|
|
|
|
|
7
|
my ($N) = @_; |
1244
|
4
|
|
|
|
|
16
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
1245
|
|
|
|
|
|
|
|
1246
|
4
|
|
|
|
|
1951
|
(-0.05473684) + ($cx * ((-0.165894739) + ($cx * (0.498947372 + $cx * (-0.334315788 +$cx *(0.055578944)) )))); |
1247
|
|
|
|
|
|
|
} |
1248
|
|
|
|
|
|
|
$window_definitions{flattop} = { |
1249
|
|
|
|
|
|
|
fn => q!flat top!, |
1250
|
|
|
|
|
|
|
}; |
1251
|
|
|
|
|
|
|
$winsubs{flattop} = \&flattop; |
1252
|
|
|
|
|
|
|
|
1253
|
|
|
|
|
|
|
sub flattop_per { |
1254
|
2
|
50
|
|
2
|
0
|
8
|
barf "flattop: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1255
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1256
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1257
|
|
|
|
|
|
|
|
1258
|
2
|
|
|
|
|
317
|
(-0.05473684) + ($cx * ((-0.165894739) + ($cx * (0.498947372 + $cx * (-0.334315788 +$cx *(0.055578944)) )))); |
1259
|
|
|
|
|
|
|
} |
1260
|
|
|
|
|
|
|
$window_definitions{flattop} = { |
1261
|
|
|
|
|
|
|
fn => q!flat top!, |
1262
|
|
|
|
|
|
|
}; |
1263
|
|
|
|
|
|
|
$winpersubs{flattop}= \&flattop_per; |
1264
|
|
|
|
|
|
|
|
1265
|
|
|
|
|
|
|
sub gaussian { |
1266
|
2
|
50
|
|
2
|
1
|
10
|
barf "gaussian: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1267
|
2
|
|
|
|
|
3
|
my ($N,$beta) = @_; |
1268
|
2
|
|
|
|
|
9
|
exp (-0.5 * ($beta * (zeroes($N)->xlinvals(-1,1)) )**2); |
1269
|
|
|
|
|
|
|
} |
1270
|
|
|
|
|
|
|
$window_definitions{gaussian} = { |
1271
|
|
|
|
|
|
|
params => [ '$beta'], |
1272
|
|
|
|
|
|
|
alias => [ 'Weierstrass'], |
1273
|
|
|
|
|
|
|
}; |
1274
|
|
|
|
|
|
|
$winsubs{gaussian} = \&gaussian; |
1275
|
|
|
|
|
|
|
|
1276
|
|
|
|
|
|
|
sub gaussian_per { |
1277
|
2
|
50
|
|
2
|
0
|
9
|
barf "gaussian: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1278
|
2
|
|
|
|
|
3
|
my ($N,$beta) = @_; |
1279
|
2
|
|
|
|
|
9
|
exp (-0.5 * ($beta * (zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N)) )**2); |
1280
|
|
|
|
|
|
|
} |
1281
|
|
|
|
|
|
|
$window_definitions{gaussian} = { |
1282
|
|
|
|
|
|
|
params => [ '$beta'], |
1283
|
|
|
|
|
|
|
alias => [ 'Weierstrass'], |
1284
|
|
|
|
|
|
|
}; |
1285
|
|
|
|
|
|
|
$winpersubs{gaussian}= \&gaussian_per; |
1286
|
|
|
|
|
|
|
|
1287
|
|
|
|
|
|
|
sub hamming { |
1288
|
6
|
50
|
|
6
|
1
|
17
|
barf "hamming: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1289
|
6
|
|
|
|
|
10
|
my ($N) = @_; |
1290
|
6
|
|
|
|
|
21
|
0.54 + -0.46 * (cos(zeroes($N)->xlinvals(0,TPI))); |
1291
|
|
|
|
|
|
|
} |
1292
|
|
|
|
|
|
|
$window_definitions{hamming} = { |
1293
|
|
|
|
|
|
|
}; |
1294
|
|
|
|
|
|
|
$winsubs{hamming} = \&hamming; |
1295
|
|
|
|
|
|
|
|
1296
|
|
|
|
|
|
|
sub hamming_per { |
1297
|
2
|
50
|
|
2
|
0
|
7
|
barf "hamming: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1298
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1299
|
2
|
|
|
|
|
9
|
0.54 + -0.46 * (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1300
|
|
|
|
|
|
|
} |
1301
|
|
|
|
|
|
|
$window_definitions{hamming} = { |
1302
|
|
|
|
|
|
|
}; |
1303
|
|
|
|
|
|
|
$winpersubs{hamming}= \&hamming_per; |
1304
|
|
|
|
|
|
|
|
1305
|
|
|
|
|
|
|
sub hamming_ex { |
1306
|
2
|
50
|
|
2
|
1
|
10
|
barf "hamming_ex: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1307
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
1308
|
2
|
|
|
|
|
36
|
0.53836 + -0.46164 * (cos(zeroes($N)->xlinvals(0,TPI))); |
1309
|
|
|
|
|
|
|
} |
1310
|
|
|
|
|
|
|
$window_definitions{hamming_ex} = { |
1311
|
|
|
|
|
|
|
fn => q!'exact' Hamming!, |
1312
|
|
|
|
|
|
|
}; |
1313
|
|
|
|
|
|
|
$winsubs{hamming_ex} = \&hamming_ex; |
1314
|
|
|
|
|
|
|
|
1315
|
|
|
|
|
|
|
sub hamming_ex_per { |
1316
|
2
|
50
|
|
2
|
0
|
6
|
barf "hamming_ex: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1317
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
1318
|
2
|
|
|
|
|
9
|
0.53836 + -0.46164 * (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1319
|
|
|
|
|
|
|
} |
1320
|
|
|
|
|
|
|
$window_definitions{hamming_ex} = { |
1321
|
|
|
|
|
|
|
fn => q!'exact' Hamming!, |
1322
|
|
|
|
|
|
|
}; |
1323
|
|
|
|
|
|
|
$winpersubs{hamming_ex}= \&hamming_ex_per; |
1324
|
|
|
|
|
|
|
|
1325
|
|
|
|
|
|
|
sub hamming_gen { |
1326
|
2
|
50
|
|
2
|
1
|
8
|
barf "hamming_gen: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1327
|
2
|
|
|
|
|
5
|
my ($N,$a) = @_; |
1328
|
2
|
|
|
|
|
9
|
$a + -(1-$a) * (cos(zeroes($N)->xlinvals(0,TPI))); |
1329
|
|
|
|
|
|
|
} |
1330
|
|
|
|
|
|
|
$window_definitions{hamming_gen} = { |
1331
|
|
|
|
|
|
|
fn => q!general Hamming!, |
1332
|
|
|
|
|
|
|
params => [ '$a'], |
1333
|
|
|
|
|
|
|
}; |
1334
|
|
|
|
|
|
|
$winsubs{hamming_gen} = \&hamming_gen; |
1335
|
|
|
|
|
|
|
|
1336
|
|
|
|
|
|
|
sub hamming_gen_per { |
1337
|
2
|
50
|
|
2
|
0
|
8
|
barf "hamming_gen: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1338
|
2
|
|
|
|
|
4
|
my ($N,$a) = @_; |
1339
|
2
|
|
|
|
|
10
|
$a + -(1-$a) * (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1340
|
|
|
|
|
|
|
} |
1341
|
|
|
|
|
|
|
$window_definitions{hamming_gen} = { |
1342
|
|
|
|
|
|
|
fn => q!general Hamming!, |
1343
|
|
|
|
|
|
|
params => [ '$a'], |
1344
|
|
|
|
|
|
|
}; |
1345
|
|
|
|
|
|
|
$winpersubs{hamming_gen}= \&hamming_gen_per; |
1346
|
|
|
|
|
|
|
|
1347
|
|
|
|
|
|
|
sub hann { |
1348
|
5
|
50
|
|
5
|
1
|
18
|
barf "hann: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1349
|
5
|
|
|
|
|
10
|
my ($N) = @_; |
1350
|
5
|
|
|
|
|
21
|
0.5 + -0.5 * (cos(zeroes($N)->xlinvals(0,TPI))); |
1351
|
|
|
|
|
|
|
} |
1352
|
|
|
|
|
|
|
$window_definitions{hann} = { |
1353
|
|
|
|
|
|
|
alias => [ 'hanning'], |
1354
|
|
|
|
|
|
|
}; |
1355
|
|
|
|
|
|
|
$winsubs{hann} = \&hann; |
1356
|
|
|
|
|
|
|
|
1357
|
|
|
|
|
|
|
sub hann_per { |
1358
|
2
|
50
|
|
2
|
0
|
10
|
barf "hann: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1359
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1360
|
2
|
|
|
|
|
9
|
0.5 + -0.5 * (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1361
|
|
|
|
|
|
|
} |
1362
|
|
|
|
|
|
|
$window_definitions{hann} = { |
1363
|
|
|
|
|
|
|
alias => [ 'hanning'], |
1364
|
|
|
|
|
|
|
}; |
1365
|
|
|
|
|
|
|
$winpersubs{hann}= \&hann_per; |
1366
|
|
|
|
|
|
|
|
1367
|
|
|
|
|
|
|
sub hann_matlab { |
1368
|
1
|
50
|
|
1
|
1
|
6
|
barf "hann_matlab: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1369
|
1
|
|
|
|
|
2
|
my ($N) = @_; |
1370
|
1
|
|
|
|
|
5
|
0.5 - 0.5 * (cos(zeroes($N)->xlinvals(TPI/($N+1),TPI *$N /($N+1)))); |
1371
|
|
|
|
|
|
|
} |
1372
|
|
|
|
|
|
|
$window_definitions{hann_matlab} = { |
1373
|
|
|
|
|
|
|
pfn => q!Hann (matlab)!, |
1374
|
|
|
|
|
|
|
fn => q!*Equivalent to the Hann window of N+2 points, with the endpoints (which are both zero) removed.!, |
1375
|
|
|
|
|
|
|
}; |
1376
|
|
|
|
|
|
|
$winsubs{hann_matlab} = \&hann_matlab; |
1377
|
|
|
|
|
|
|
|
1378
|
|
|
|
|
|
|
sub hann_poisson { |
1379
|
1
|
50
|
|
1
|
1
|
6
|
barf "hann_poisson: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1380
|
1
|
|
|
|
|
2
|
my ($N,$alpha) = @_; |
1381
|
1
|
|
|
|
|
5
|
0.5 * (1 + (cos(zeroes($N)->xlinvals(-(PI),PI)))) * exp (-$alpha * abs (zeroes($N)->xlinvals(-1,1))); |
1382
|
|
|
|
|
|
|
} |
1383
|
|
|
|
|
|
|
$window_definitions{hann_poisson} = { |
1384
|
|
|
|
|
|
|
fn => q!Hann-Poisson!, |
1385
|
|
|
|
|
|
|
params => [ '$alpha'], |
1386
|
|
|
|
|
|
|
}; |
1387
|
|
|
|
|
|
|
$winsubs{hann_poisson} = \&hann_poisson; |
1388
|
|
|
|
|
|
|
|
1389
|
|
|
|
|
|
|
sub hann_poisson_per { |
1390
|
0
|
0
|
|
0
|
0
|
0
|
barf "hann_poisson: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1391
|
0
|
|
|
|
|
0
|
my ($N,$alpha) = @_; |
1392
|
0
|
|
|
|
|
0
|
0.5 * (1 + (cos(zeroes($N)->xlinvals(-(PI), (-(PI)+PI*($N-1))/$N)))) * exp (-$alpha * abs (zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N))); |
1393
|
|
|
|
|
|
|
} |
1394
|
|
|
|
|
|
|
$window_definitions{hann_poisson} = { |
1395
|
|
|
|
|
|
|
fn => q!Hann-Poisson!, |
1396
|
|
|
|
|
|
|
params => [ '$alpha'], |
1397
|
|
|
|
|
|
|
}; |
1398
|
|
|
|
|
|
|
$winpersubs{hann_poisson}= \&hann_poisson_per; |
1399
|
|
|
|
|
|
|
|
1400
|
|
|
|
|
|
|
sub kaiser { |
1401
|
0
|
0
|
|
0
|
1
|
0
|
barf "kaiser: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1402
|
0
|
|
|
|
|
0
|
my ($N,$beta) = @_; |
1403
|
|
|
|
|
|
|
|
1404
|
0
|
0
|
|
|
|
0
|
barf "kaiser: PDL::GSLSF not installed" unless $HAVE_BESSEL; |
1405
|
0
|
|
|
|
|
0
|
$beta *= PI; |
1406
|
0
|
|
|
|
|
0
|
my @n = PDL::GSLSF::BESSEL::gsl_sf_bessel_In ($beta * sqrt(1 - (zeroes($N)->xlinvals(-1,1)) **2),0); |
1407
|
0
|
|
|
|
|
0
|
my @d = PDL::GSLSF::BESSEL::gsl_sf_bessel_In($beta,0); |
1408
|
0
|
|
|
|
|
0
|
(shift @n)/(shift @d); |
1409
|
|
|
|
|
|
|
} |
1410
|
|
|
|
|
|
|
$window_definitions{kaiser} = { |
1411
|
|
|
|
|
|
|
params => [ '$beta'], |
1412
|
|
|
|
|
|
|
alias => [ 'Kaiser-Bessel'], |
1413
|
|
|
|
|
|
|
}; |
1414
|
|
|
|
|
|
|
$winsubs{kaiser} = \&kaiser; |
1415
|
|
|
|
|
|
|
|
1416
|
|
|
|
|
|
|
sub kaiser_per { |
1417
|
0
|
0
|
|
0
|
0
|
0
|
barf "kaiser: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1418
|
0
|
|
|
|
|
0
|
my ($N,$beta) = @_; |
1419
|
|
|
|
|
|
|
|
1420
|
0
|
0
|
|
|
|
0
|
barf "kaiser: PDL::GSLSF not installed" unless $HAVE_BESSEL; |
1421
|
0
|
|
|
|
|
0
|
$beta *= PI; |
1422
|
0
|
|
|
|
|
0
|
my @n = PDL::GSLSF::BESSEL::gsl_sf_bessel_In ($beta * sqrt(1 - (zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N)) **2),0); |
1423
|
0
|
|
|
|
|
0
|
my @d = PDL::GSLSF::BESSEL::gsl_sf_bessel_In($beta,0); |
1424
|
0
|
|
|
|
|
0
|
(shift @n)/(shift @d); |
1425
|
|
|
|
|
|
|
} |
1426
|
|
|
|
|
|
|
$window_definitions{kaiser} = { |
1427
|
|
|
|
|
|
|
params => [ '$beta'], |
1428
|
|
|
|
|
|
|
alias => [ 'Kaiser-Bessel'], |
1429
|
|
|
|
|
|
|
}; |
1430
|
|
|
|
|
|
|
$winpersubs{kaiser}= \&kaiser_per; |
1431
|
|
|
|
|
|
|
|
1432
|
|
|
|
|
|
|
sub lanczos { |
1433
|
3
|
50
|
|
3
|
1
|
11
|
barf "lanczos: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1434
|
3
|
|
|
|
|
9
|
my ($N) = @_; |
1435
|
|
|
|
|
|
|
|
1436
|
3
|
|
|
|
|
11
|
my $x = PI * (zeroes($N)->xlinvals(-1,1)); |
1437
|
3
|
|
|
|
|
807
|
my $res = sin($x)/$x; |
1438
|
3
|
|
|
|
|
1020
|
my $mid; |
1439
|
3
|
100
|
|
|
|
20
|
$mid = int($N/2), $res->slice($mid) .= 1 if $N % 2; |
1440
|
3
|
|
|
|
|
61
|
$res;; |
1441
|
|
|
|
|
|
|
} |
1442
|
|
|
|
|
|
|
$window_definitions{lanczos} = { |
1443
|
|
|
|
|
|
|
alias => [ 'sinc'], |
1444
|
|
|
|
|
|
|
}; |
1445
|
|
|
|
|
|
|
$winsubs{lanczos} = \&lanczos; |
1446
|
|
|
|
|
|
|
|
1447
|
|
|
|
|
|
|
sub lanczos_per { |
1448
|
2
|
50
|
|
2
|
0
|
9
|
barf "lanczos: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1449
|
2
|
|
|
|
|
3
|
my ($N) = @_; |
1450
|
|
|
|
|
|
|
|
1451
|
2
|
|
|
|
|
9
|
my $x = PI * (zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N)); |
1452
|
2
|
|
|
|
|
363
|
my $res = sin($x)/$x; |
1453
|
2
|
|
|
|
|
64
|
my $mid; |
1454
|
2
|
100
|
|
|
|
11
|
$mid = int($N/2), $res->slice($mid) .= 1 unless $N % 2; |
1455
|
2
|
|
|
|
|
49
|
$res;; |
1456
|
|
|
|
|
|
|
} |
1457
|
|
|
|
|
|
|
$window_definitions{lanczos} = { |
1458
|
|
|
|
|
|
|
alias => [ 'sinc'], |
1459
|
|
|
|
|
|
|
}; |
1460
|
|
|
|
|
|
|
$winpersubs{lanczos}= \&lanczos_per; |
1461
|
|
|
|
|
|
|
|
1462
|
|
|
|
|
|
|
sub nuttall { |
1463
|
2
|
50
|
|
2
|
1
|
8
|
barf "nuttall: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1464
|
2
|
|
|
|
|
5
|
my ($N) = @_; |
1465
|
2
|
|
|
|
|
7
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
1466
|
|
|
|
|
|
|
|
1467
|
2
|
|
|
|
|
332
|
(0.2269824) + ($cx * ((-0.4572542) + ($cx * (0.273199 + $cx * (-0.0425644) )))); |
1468
|
|
|
|
|
|
|
} |
1469
|
|
|
|
|
|
|
$window_definitions{nuttall} = { |
1470
|
|
|
|
|
|
|
}; |
1471
|
|
|
|
|
|
|
$winsubs{nuttall} = \&nuttall; |
1472
|
|
|
|
|
|
|
|
1473
|
|
|
|
|
|
|
sub nuttall_per { |
1474
|
2
|
50
|
|
2
|
0
|
7
|
barf "nuttall: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1475
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1476
|
2
|
|
|
|
|
8
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1477
|
|
|
|
|
|
|
|
1478
|
2
|
|
|
|
|
301
|
(0.2269824) + ($cx * ((-0.4572542) + ($cx * (0.273199 + $cx * (-0.0425644) )))); |
1479
|
|
|
|
|
|
|
} |
1480
|
|
|
|
|
|
|
$window_definitions{nuttall} = { |
1481
|
|
|
|
|
|
|
}; |
1482
|
|
|
|
|
|
|
$winpersubs{nuttall}= \&nuttall_per; |
1483
|
|
|
|
|
|
|
|
1484
|
|
|
|
|
|
|
sub nuttall1 { |
1485
|
0
|
0
|
|
0
|
1
|
0
|
barf "nuttall1: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1486
|
0
|
|
|
|
|
0
|
my ($N) = @_; |
1487
|
0
|
|
|
|
|
0
|
my $cx = (cos(zeroes($N)->xlinvals(0,TPI))); |
1488
|
|
|
|
|
|
|
|
1489
|
0
|
|
|
|
|
0
|
(0.211536) + ($cx * ((-0.449584) + ($cx * (0.288464 + $cx * (-0.050416) )))); |
1490
|
|
|
|
|
|
|
} |
1491
|
|
|
|
|
|
|
$window_definitions{nuttall1} = { |
1492
|
|
|
|
|
|
|
pfn => q!Nuttall (v1)!, |
1493
|
|
|
|
|
|
|
fn => q!*A window referred to as the Nuttall window.!, |
1494
|
|
|
|
|
|
|
}; |
1495
|
|
|
|
|
|
|
$winsubs{nuttall1} = \&nuttall1; |
1496
|
|
|
|
|
|
|
|
1497
|
|
|
|
|
|
|
sub nuttall1_per { |
1498
|
0
|
0
|
|
0
|
0
|
0
|
barf "nuttall1: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1499
|
0
|
|
|
|
|
0
|
my ($N) = @_; |
1500
|
0
|
|
|
|
|
0
|
my $cx = (cos(zeroes($N)->xlinvals(0, TPI*($N-1)/$N))); |
1501
|
|
|
|
|
|
|
|
1502
|
0
|
|
|
|
|
0
|
(0.211536) + ($cx * ((-0.449584) + ($cx * (0.288464 + $cx * (-0.050416) )))); |
1503
|
|
|
|
|
|
|
} |
1504
|
|
|
|
|
|
|
$window_definitions{nuttall1} = { |
1505
|
|
|
|
|
|
|
pfn => q!Nuttall (v1)!, |
1506
|
|
|
|
|
|
|
fn => q!*A window referred to as the Nuttall window.!, |
1507
|
|
|
|
|
|
|
}; |
1508
|
|
|
|
|
|
|
$winpersubs{nuttall1}= \&nuttall1_per; |
1509
|
|
|
|
|
|
|
|
1510
|
|
|
|
|
|
|
sub parzen { |
1511
|
3
|
50
|
|
3
|
1
|
11
|
barf "parzen: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1512
|
3
|
|
|
|
|
5
|
my ($N) = @_; |
1513
|
|
|
|
|
|
|
|
1514
|
3
|
|
|
|
|
43
|
my $x = zeroes($N)->xlinvals(-1,1); |
1515
|
3
|
|
|
|
|
842
|
my $x1 = $x->where($x <= -.5); |
1516
|
3
|
|
|
|
|
572
|
my $x2 = $x->where( ($x < .5) & ($x > -.5) ); |
1517
|
3
|
|
|
|
|
305
|
my $x3 = $x->where($x >= .5); |
1518
|
3
|
|
|
|
|
230
|
$x1 .= 2 * (1-abs($x1))**3; |
1519
|
|
|
|
|
|
|
# $x3 .= 2 * (1-abs($x3))**3; |
1520
|
3
|
|
|
|
|
1063
|
$x3 .= $x1->slice('-1:0:-1'); |
1521
|
3
|
|
|
|
|
410
|
$x2 .= 1 - 6 * $x2**2 *(1-abs($x2)); |
1522
|
3
|
|
|
|
|
480
|
return $x; |
1523
|
|
|
|
|
|
|
} |
1524
|
|
|
|
|
|
|
$window_definitions{parzen} = { |
1525
|
|
|
|
|
|
|
alias => [ 'Jackson','Valle-Poussin'], |
1526
|
|
|
|
|
|
|
}; |
1527
|
|
|
|
|
|
|
$winsubs{parzen} = \&parzen; |
1528
|
|
|
|
|
|
|
|
1529
|
|
|
|
|
|
|
sub parzen_per { |
1530
|
2
|
50
|
|
2
|
0
|
8
|
barf "parzen: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1531
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1532
|
|
|
|
|
|
|
|
1533
|
2
|
|
|
|
|
8
|
my $x = zeroes($N)->xlinvals(-1,(-1 + ($N-1))/($N)); |
1534
|
2
|
|
|
|
|
346
|
my $x1 = $x->where($x <= -.5); |
1535
|
2
|
|
|
|
|
146
|
my $x2 = $x->where( ($x < .5) & ($x > -.5) ); |
1536
|
2
|
|
|
|
|
120
|
my $x3 = $x->where($x >= .5); |
1537
|
2
|
|
|
|
|
88
|
$x1 .= 2 * (1-abs($x1))**3; |
1538
|
|
|
|
|
|
|
# $x3 .= 2 * (1-abs($x3))**3; |
1539
|
2
|
|
|
|
|
150
|
$x3 .= $x1->slice('-1:1:-1'); |
1540
|
2
|
|
|
|
|
91
|
$x2 .= 1 - 6 * $x2**2 *(1-abs($x2)); |
1541
|
2
|
|
|
|
|
183
|
return $x; |
1542
|
|
|
|
|
|
|
} |
1543
|
|
|
|
|
|
|
$window_definitions{parzen} = { |
1544
|
|
|
|
|
|
|
alias => [ 'Jackson','Valle-Poussin'], |
1545
|
|
|
|
|
|
|
}; |
1546
|
|
|
|
|
|
|
$winpersubs{parzen}= \&parzen_per; |
1547
|
|
|
|
|
|
|
|
1548
|
|
|
|
|
|
|
sub parzen_octave { |
1549
|
0
|
0
|
|
0
|
1
|
0
|
barf "parzen_octave: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1550
|
0
|
|
|
|
|
0
|
my ($N) = @_; |
1551
|
|
|
|
|
|
|
|
1552
|
0
|
|
|
|
|
0
|
my $L = $N-1; |
1553
|
0
|
|
|
|
|
0
|
my $r = ($L/2); |
1554
|
0
|
|
|
|
|
0
|
my $r4 = ($r/2); |
1555
|
0
|
|
|
|
|
0
|
my $n = sequence(2*$r+1)-$r; |
1556
|
0
|
|
|
|
|
0
|
my $n1 = $n->where(abs($n) <= $r4); |
1557
|
0
|
|
|
|
|
0
|
my $n2 = $n->where($n > $r4); |
1558
|
0
|
|
|
|
|
0
|
my $n3 = $n->where($n < -$r4); |
1559
|
0
|
|
|
|
|
0
|
$n1 .= 1 -6.*(abs($n1)/($N/2))**2 + 6*(abs($n1)/($N/2))**3; |
1560
|
0
|
|
|
|
|
0
|
$n2 .= 2.*(1-abs($n2)/($N/2))**3; |
1561
|
0
|
|
|
|
|
0
|
$n3 .= 2.*(1-abs($n3)/($N/2))**3; |
1562
|
0
|
|
|
|
|
0
|
$n; |
1563
|
|
|
|
|
|
|
; |
1564
|
|
|
|
|
|
|
} |
1565
|
|
|
|
|
|
|
$window_definitions{parzen_octave} = { |
1566
|
|
|
|
|
|
|
fn => q!Parzen!, |
1567
|
|
|
|
|
|
|
}; |
1568
|
|
|
|
|
|
|
$winsubs{parzen_octave} = \&parzen_octave; |
1569
|
|
|
|
|
|
|
|
1570
|
|
|
|
|
|
|
sub poisson { |
1571
|
3
|
50
|
|
3
|
1
|
10
|
barf "poisson: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1572
|
3
|
|
|
|
|
6
|
my ($N,$alpha) = @_; |
1573
|
3
|
|
|
|
|
12
|
exp (-$alpha * abs (zeroes($N)->xlinvals(-1,1))); |
1574
|
|
|
|
|
|
|
} |
1575
|
|
|
|
|
|
|
$window_definitions{poisson} = { |
1576
|
|
|
|
|
|
|
params => [ '$alpha'], |
1577
|
|
|
|
|
|
|
}; |
1578
|
|
|
|
|
|
|
$winsubs{poisson} = \&poisson; |
1579
|
|
|
|
|
|
|
|
1580
|
|
|
|
|
|
|
sub poisson_per { |
1581
|
2
|
50
|
|
2
|
0
|
7
|
barf "poisson: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1582
|
2
|
|
|
|
|
5
|
my ($N,$alpha) = @_; |
1583
|
2
|
|
|
|
|
8
|
exp (-$alpha * abs (zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N))); |
1584
|
|
|
|
|
|
|
} |
1585
|
|
|
|
|
|
|
$window_definitions{poisson} = { |
1586
|
|
|
|
|
|
|
params => [ '$alpha'], |
1587
|
|
|
|
|
|
|
}; |
1588
|
|
|
|
|
|
|
$winpersubs{poisson}= \&poisson_per; |
1589
|
|
|
|
|
|
|
|
1590
|
|
|
|
|
|
|
sub rectangular { |
1591
|
4
|
50
|
|
4
|
1
|
16
|
barf "rectangular: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1592
|
4
|
|
|
|
|
9
|
my ($N) = @_; |
1593
|
4
|
|
|
|
|
48
|
ones($N); |
1594
|
|
|
|
|
|
|
} |
1595
|
|
|
|
|
|
|
$window_definitions{rectangular} = { |
1596
|
|
|
|
|
|
|
alias => [ 'dirichlet','boxcar'], |
1597
|
|
|
|
|
|
|
}; |
1598
|
|
|
|
|
|
|
$winsubs{rectangular} = \&rectangular; |
1599
|
|
|
|
|
|
|
|
1600
|
|
|
|
|
|
|
sub rectangular_per { |
1601
|
2
|
50
|
|
2
|
0
|
8
|
barf "rectangular: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1602
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1603
|
2
|
|
|
|
|
8
|
ones($N); |
1604
|
|
|
|
|
|
|
} |
1605
|
|
|
|
|
|
|
$window_definitions{rectangular} = { |
1606
|
|
|
|
|
|
|
alias => [ 'dirichlet','boxcar'], |
1607
|
|
|
|
|
|
|
}; |
1608
|
|
|
|
|
|
|
$winpersubs{rectangular}= \&rectangular_per; |
1609
|
|
|
|
|
|
|
|
1610
|
|
|
|
|
|
|
sub triangular { |
1611
|
4
|
50
|
|
4
|
1
|
16
|
barf "triangular: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1612
|
4
|
|
|
|
|
7
|
my ($N) = @_; |
1613
|
4
|
|
|
|
|
18
|
1 - abs (zeroes($N)->xlinvals(-($N-1)/$N,($N-1)/$N)); |
1614
|
|
|
|
|
|
|
} |
1615
|
|
|
|
|
|
|
$window_definitions{triangular} = { |
1616
|
|
|
|
|
|
|
}; |
1617
|
|
|
|
|
|
|
$winsubs{triangular} = \&triangular; |
1618
|
|
|
|
|
|
|
|
1619
|
|
|
|
|
|
|
sub triangular_per { |
1620
|
2
|
50
|
|
2
|
0
|
9
|
barf "triangular: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1621
|
2
|
|
|
|
|
4
|
my ($N) = @_; |
1622
|
2
|
|
|
|
|
9
|
1 - abs (zeroes($N)->xlinvals(-$N/($N+1),-1/($N+1)+($N-1)/($N+1))); |
1623
|
|
|
|
|
|
|
} |
1624
|
|
|
|
|
|
|
$window_definitions{triangular} = { |
1625
|
|
|
|
|
|
|
}; |
1626
|
|
|
|
|
|
|
$winpersubs{triangular}= \&triangular_per; |
1627
|
|
|
|
|
|
|
|
1628
|
|
|
|
|
|
|
sub tukey { |
1629
|
4
|
50
|
|
4
|
1
|
13
|
barf "tukey: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1630
|
4
|
|
|
|
|
8
|
my ($N,$alpha) = @_; |
1631
|
|
|
|
|
|
|
|
1632
|
4
|
50
|
33
|
|
|
26
|
barf("tukey: alpha must be between 0 and 1") unless |
1633
|
|
|
|
|
|
|
$alpha >=0 and $alpha <= 1; |
1634
|
4
|
50
|
|
|
|
12
|
return ones($N) if $alpha == 0; |
1635
|
4
|
|
|
|
|
17
|
my $x = zeroes($N)->xlinvals(0,1); |
1636
|
4
|
|
|
|
|
941
|
my $x1 = $x->where($x < $alpha/2); |
1637
|
4
|
|
|
|
|
668
|
my $x2 = $x->where( ($x <= 1-$alpha/2) & ($x >= $alpha/2) ); |
1638
|
4
|
|
|
|
|
377
|
my $x3 = $x->where($x > 1 - $alpha/2); |
1639
|
4
|
|
|
|
|
298
|
$x1 .= 0.5 * ( 1 + cos( PI * (2*$x1/$alpha -1))); |
1640
|
4
|
|
|
|
|
692
|
$x2 .= 1; |
1641
|
4
|
|
|
|
|
238
|
$x3 .= $x1->slice('-1:0:-1'); |
1642
|
4
|
|
|
|
|
250
|
return $x; |
1643
|
|
|
|
|
|
|
} |
1644
|
|
|
|
|
|
|
$window_definitions{tukey} = { |
1645
|
|
|
|
|
|
|
params => [ '$alpha'], |
1646
|
|
|
|
|
|
|
alias => [ 'tapered cosine'], |
1647
|
|
|
|
|
|
|
}; |
1648
|
|
|
|
|
|
|
$winsubs{tukey} = \&tukey; |
1649
|
|
|
|
|
|
|
|
1650
|
|
|
|
|
|
|
sub tukey_per { |
1651
|
2
|
50
|
|
2
|
0
|
8
|
barf "tukey: 2 arguments expected. Got " . scalar(@_) . ' arguments.' unless @_ == 2; |
1652
|
2
|
|
|
|
|
3
|
my ($N,$alpha) = @_; |
1653
|
|
|
|
|
|
|
|
1654
|
2
|
50
|
33
|
|
|
14
|
barf("tukey: alpha must be between 0 and 1") unless |
1655
|
|
|
|
|
|
|
$alpha >=0 and $alpha <= 1; |
1656
|
2
|
50
|
|
|
|
8
|
return ones($N) if $alpha == 0; |
1657
|
2
|
|
|
|
|
8
|
my $x = zeroes($N)->xlinvals(0,($N-1)/$N); |
1658
|
2
|
|
|
|
|
323
|
my $x1 = $x->where($x < $alpha/2); |
1659
|
2
|
|
|
|
|
167
|
my $x2 = $x->where( ($x <= 1-$alpha/2) & ($x >= $alpha/2) ); |
1660
|
2
|
|
|
|
|
131
|
my $x3 = $x->where($x > 1 - $alpha/2); |
1661
|
2
|
|
|
|
|
98
|
$x1 .= 0.5 * ( 1 + cos( PI * (2*$x1/$alpha -1))); |
1662
|
2
|
|
|
|
|
230
|
$x2 .= 1; |
1663
|
2
|
|
|
|
|
47
|
$x3 .= $x1->slice('-1:1:-1'); |
1664
|
2
|
|
|
|
|
92
|
return $x; |
1665
|
|
|
|
|
|
|
} |
1666
|
|
|
|
|
|
|
$window_definitions{tukey} = { |
1667
|
|
|
|
|
|
|
params => [ '$alpha'], |
1668
|
|
|
|
|
|
|
alias => [ 'tapered cosine'], |
1669
|
|
|
|
|
|
|
}; |
1670
|
|
|
|
|
|
|
$winpersubs{tukey}= \&tukey_per; |
1671
|
|
|
|
|
|
|
|
1672
|
|
|
|
|
|
|
sub welch { |
1673
|
3
|
50
|
|
3
|
1
|
12
|
barf "welch: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1674
|
3
|
|
|
|
|
5
|
my ($N) = @_; |
1675
|
3
|
|
|
|
|
13
|
1 - (zeroes($N)->xlinvals(-1,1))**2; |
1676
|
|
|
|
|
|
|
} |
1677
|
|
|
|
|
|
|
$window_definitions{welch} = { |
1678
|
|
|
|
|
|
|
alias => [ 'Riez','Bochner','Parzen','parabolic'], |
1679
|
|
|
|
|
|
|
}; |
1680
|
|
|
|
|
|
|
$winsubs{welch} = \&welch; |
1681
|
|
|
|
|
|
|
|
1682
|
|
|
|
|
|
|
sub welch_per { |
1683
|
2
|
50
|
|
2
|
0
|
10
|
barf "welch: 1 argument expected. Got " . scalar(@_) . ' arguments.' unless @_ == 1; |
1684
|
2
|
|
|
|
|
6
|
my ($N) = @_; |
1685
|
2
|
|
|
|
|
9
|
1 - (zeroes($N)->xlinvals(-1, (-1+1*($N-1))/$N))**2; |
1686
|
|
|
|
|
|
|
} |
1687
|
|
|
|
|
|
|
$window_definitions{welch} = { |
1688
|
|
|
|
|
|
|
alias => [ 'Riez','Bochner','Parzen','parabolic'], |
1689
|
|
|
|
|
|
|
}; |
1690
|
|
|
|
|
|
|
$winpersubs{welch}= \&welch_per; |
1691
|
|
|
|
|
|
|
|
1692
|
|
|
|
|
|
|
|
1693
|
|
|
|
|
|
|
=head1 Symmetric window functions |
1694
|
|
|
|
|
|
|
|
1695
|
|
|
|
|
|
|
|
1696
|
|
|
|
|
|
|
|
1697
|
|
|
|
|
|
|
=head2 bartlett($N) |
1698
|
|
|
|
|
|
|
|
1699
|
|
|
|
|
|
|
The Bartlett window. (Ref 1). Another name for this window is the fejer window. This window is defined by |
1700
|
|
|
|
|
|
|
|
1701
|
|
|
|
|
|
|
1 - abs arr, |
1702
|
|
|
|
|
|
|
|
1703
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
1704
|
|
|
|
|
|
|
See also L. |
1705
|
|
|
|
|
|
|
|
1706
|
|
|
|
|
|
|
|
1707
|
|
|
|
|
|
|
=cut |
1708
|
|
|
|
|
|
|
|
1709
|
|
|
|
|
|
|
|
1710
|
|
|
|
|
|
|
=head2 bartlett_hann($N) |
1711
|
|
|
|
|
|
|
|
1712
|
|
|
|
|
|
|
The Bartlett-Hann window. Another name for this window is the Modified Bartlett-Hann window. This window is defined by |
1713
|
|
|
|
|
|
|
|
1714
|
|
|
|
|
|
|
0.62 - 0.48 * abs arr + 0.38* arr1, |
1715
|
|
|
|
|
|
|
|
1716
|
|
|
|
|
|
|
where the points in arr range from -1/2 through 1/2, and arr1 are the cos of points ranging from -PI through PI. |
1717
|
|
|
|
|
|
|
|
1718
|
|
|
|
|
|
|
|
1719
|
|
|
|
|
|
|
=cut |
1720
|
|
|
|
|
|
|
|
1721
|
|
|
|
|
|
|
|
1722
|
|
|
|
|
|
|
=head2 blackman($N) |
1723
|
|
|
|
|
|
|
|
1724
|
|
|
|
|
|
|
The 'classic' Blackman window. (Ref 1). One of the Blackman-Harris family, with coefficients |
1725
|
|
|
|
|
|
|
|
1726
|
|
|
|
|
|
|
a0 = 0.42, a1 = 0.5, a2 = 0.08 |
1727
|
|
|
|
|
|
|
|
1728
|
|
|
|
|
|
|
|
1729
|
|
|
|
|
|
|
|
1730
|
|
|
|
|
|
|
=cut |
1731
|
|
|
|
|
|
|
|
1732
|
|
|
|
|
|
|
|
1733
|
|
|
|
|
|
|
=head2 blackman_bnh($N) |
1734
|
|
|
|
|
|
|
|
1735
|
|
|
|
|
|
|
The Blackman-Harris (bnh) window. An improved version of the 3-term Blackman-Harris window given by Nuttall (Ref 2, p. 89). One of the Blackman-Harris family, with coefficients |
1736
|
|
|
|
|
|
|
|
1737
|
|
|
|
|
|
|
a0 = 0.4243801, a1 = 0.4973406, a2 = 0.0782793 |
1738
|
|
|
|
|
|
|
|
1739
|
|
|
|
|
|
|
|
1740
|
|
|
|
|
|
|
|
1741
|
|
|
|
|
|
|
=cut |
1742
|
|
|
|
|
|
|
|
1743
|
|
|
|
|
|
|
|
1744
|
|
|
|
|
|
|
=head2 blackman_ex($N) |
1745
|
|
|
|
|
|
|
|
1746
|
|
|
|
|
|
|
The 'exact' Blackman window. (Ref 1). One of the Blackman-Harris family, with coefficients |
1747
|
|
|
|
|
|
|
|
1748
|
|
|
|
|
|
|
a0 = 0.426590713671539, a1 = 0.496560619088564, a2 = 0.0768486672398968 |
1749
|
|
|
|
|
|
|
|
1750
|
|
|
|
|
|
|
|
1751
|
|
|
|
|
|
|
|
1752
|
|
|
|
|
|
|
=cut |
1753
|
|
|
|
|
|
|
|
1754
|
|
|
|
|
|
|
|
1755
|
|
|
|
|
|
|
=head2 blackman_gen($N,$alpha) |
1756
|
|
|
|
|
|
|
|
1757
|
|
|
|
|
|
|
The General classic Blackman window. A single parameter family of the 3-term Blackman window. This window is defined by |
1758
|
|
|
|
|
|
|
|
1759
|
|
|
|
|
|
|
my $cx = arr; |
1760
|
|
|
|
|
|
|
|
1761
|
|
|
|
|
|
|
(.5 - $alpha) + ($cx * ((-.5) + ($cx * ($alpha)))), |
1762
|
|
|
|
|
|
|
|
1763
|
|
|
|
|
|
|
where the points in arr are the cos of points ranging from 0 through 2PI. |
1764
|
|
|
|
|
|
|
|
1765
|
|
|
|
|
|
|
|
1766
|
|
|
|
|
|
|
=cut |
1767
|
|
|
|
|
|
|
|
1768
|
|
|
|
|
|
|
|
1769
|
|
|
|
|
|
|
=head2 blackman_gen3($N,$a0,$a1,$a2) |
1770
|
|
|
|
|
|
|
|
1771
|
|
|
|
|
|
|
The general form of the Blackman family. One of the Blackman-Harris family, with coefficients |
1772
|
|
|
|
|
|
|
|
1773
|
|
|
|
|
|
|
a0 = $a0, a1 = $a1, a2 = $a2 |
1774
|
|
|
|
|
|
|
|
1775
|
|
|
|
|
|
|
|
1776
|
|
|
|
|
|
|
|
1777
|
|
|
|
|
|
|
=cut |
1778
|
|
|
|
|
|
|
|
1779
|
|
|
|
|
|
|
|
1780
|
|
|
|
|
|
|
=head2 blackman_gen4($N,$a0,$a1,$a2,$a3) |
1781
|
|
|
|
|
|
|
|
1782
|
|
|
|
|
|
|
The general 4-term Blackman-Harris window. One of the Blackman-Harris family, with coefficients |
1783
|
|
|
|
|
|
|
|
1784
|
|
|
|
|
|
|
a0 = $a0, a1 = $a1, a2 = $a2, a3 = $a3 |
1785
|
|
|
|
|
|
|
|
1786
|
|
|
|
|
|
|
|
1787
|
|
|
|
|
|
|
|
1788
|
|
|
|
|
|
|
=cut |
1789
|
|
|
|
|
|
|
|
1790
|
|
|
|
|
|
|
|
1791
|
|
|
|
|
|
|
=head2 blackman_gen5($N,$a0,$a1,$a2,$a3,$a4) |
1792
|
|
|
|
|
|
|
|
1793
|
|
|
|
|
|
|
The general 5-term Blackman-Harris window. One of the Blackman-Harris family, with coefficients |
1794
|
|
|
|
|
|
|
|
1795
|
|
|
|
|
|
|
a0 = $a0, a1 = $a1, a2 = $a2, a3 = $a3, a4 = $a4 |
1796
|
|
|
|
|
|
|
|
1797
|
|
|
|
|
|
|
|
1798
|
|
|
|
|
|
|
|
1799
|
|
|
|
|
|
|
=cut |
1800
|
|
|
|
|
|
|
|
1801
|
|
|
|
|
|
|
|
1802
|
|
|
|
|
|
|
=head2 blackman_harris($N) |
1803
|
|
|
|
|
|
|
|
1804
|
|
|
|
|
|
|
The Blackman-Harris window. (Ref 1). One of the Blackman-Harris family, with coefficients |
1805
|
|
|
|
|
|
|
|
1806
|
|
|
|
|
|
|
a0 = 0.422323, a1 = 0.49755, a2 = 0.07922 |
1807
|
|
|
|
|
|
|
|
1808
|
|
|
|
|
|
|
Another name for this window is the Minimum three term (sample) Blackman-Harris window. |
1809
|
|
|
|
|
|
|
|
1810
|
|
|
|
|
|
|
=cut |
1811
|
|
|
|
|
|
|
|
1812
|
|
|
|
|
|
|
|
1813
|
|
|
|
|
|
|
=head2 blackman_harris4($N) |
1814
|
|
|
|
|
|
|
|
1815
|
|
|
|
|
|
|
The minimum (sidelobe) four term Blackman-Harris window. (Ref 1). One of the Blackman-Harris family, with coefficients |
1816
|
|
|
|
|
|
|
|
1817
|
|
|
|
|
|
|
a0 = 0.35875, a1 = 0.48829, a2 = 0.14128, a3 = 0.01168 |
1818
|
|
|
|
|
|
|
|
1819
|
|
|
|
|
|
|
Another name for this window is the Blackman-Harris window. |
1820
|
|
|
|
|
|
|
|
1821
|
|
|
|
|
|
|
=cut |
1822
|
|
|
|
|
|
|
|
1823
|
|
|
|
|
|
|
|
1824
|
|
|
|
|
|
|
=head2 blackman_nuttall($N) |
1825
|
|
|
|
|
|
|
|
1826
|
|
|
|
|
|
|
The Blackman-Nuttall window. One of the Blackman-Harris family, with coefficients |
1827
|
|
|
|
|
|
|
|
1828
|
|
|
|
|
|
|
a0 = 0.3635819, a1 = 0.4891775, a2 = 0.1365995, a3 = 0.0106411 |
1829
|
|
|
|
|
|
|
|
1830
|
|
|
|
|
|
|
|
1831
|
|
|
|
|
|
|
|
1832
|
|
|
|
|
|
|
=cut |
1833
|
|
|
|
|
|
|
|
1834
|
|
|
|
|
|
|
|
1835
|
|
|
|
|
|
|
=head2 bohman($N) |
1836
|
|
|
|
|
|
|
|
1837
|
|
|
|
|
|
|
The Bohman window. (Ref 1). This window is defined by |
1838
|
|
|
|
|
|
|
|
1839
|
|
|
|
|
|
|
my $x = abs(arr); |
1840
|
|
|
|
|
|
|
(1-$x)*cos(PI*$x) +(1/PI)*sin(PI*$x), |
1841
|
|
|
|
|
|
|
|
1842
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
1843
|
|
|
|
|
|
|
|
1844
|
|
|
|
|
|
|
|
1845
|
|
|
|
|
|
|
=cut |
1846
|
|
|
|
|
|
|
|
1847
|
|
|
|
|
|
|
|
1848
|
|
|
|
|
|
|
=head2 cauchy($N,$alpha) |
1849
|
|
|
|
|
|
|
|
1850
|
|
|
|
|
|
|
The Cauchy window. (Ref 1). Other names for this window are: Abel, Poisson. This window is defined by |
1851
|
|
|
|
|
|
|
|
1852
|
|
|
|
|
|
|
1 / (1 + (arr * $alpha)**2), |
1853
|
|
|
|
|
|
|
|
1854
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
1855
|
|
|
|
|
|
|
|
1856
|
|
|
|
|
|
|
|
1857
|
|
|
|
|
|
|
=cut |
1858
|
|
|
|
|
|
|
|
1859
|
|
|
|
|
|
|
|
1860
|
|
|
|
|
|
|
=head2 chebyshev($N,$at) |
1861
|
|
|
|
|
|
|
|
1862
|
|
|
|
|
|
|
The Chebyshev window. The frequency response of this window has C<$at> dB of attenuation in the stop-band. |
1863
|
|
|
|
|
|
|
Another name for this window is the Dolph-Chebyshev window. No periodic version of this window is defined. |
1864
|
|
|
|
|
|
|
This routine gives the same result as the routine B in Octave 3.6.2. |
1865
|
|
|
|
|
|
|
|
1866
|
|
|
|
|
|
|
|
1867
|
|
|
|
|
|
|
=cut |
1868
|
|
|
|
|
|
|
|
1869
|
|
|
|
|
|
|
|
1870
|
|
|
|
|
|
|
=head2 cos_alpha($N,$alpha) |
1871
|
|
|
|
|
|
|
|
1872
|
|
|
|
|
|
|
The Cos_alpha window. (Ref 1). Another name for this window is the Power-of-cosine window. This window is defined by |
1873
|
|
|
|
|
|
|
|
1874
|
|
|
|
|
|
|
arr**$alpha , |
1875
|
|
|
|
|
|
|
|
1876
|
|
|
|
|
|
|
where the points in arr are the sin of points ranging from 0 through PI. |
1877
|
|
|
|
|
|
|
|
1878
|
|
|
|
|
|
|
|
1879
|
|
|
|
|
|
|
=cut |
1880
|
|
|
|
|
|
|
|
1881
|
|
|
|
|
|
|
|
1882
|
|
|
|
|
|
|
=head2 cosine($N) |
1883
|
|
|
|
|
|
|
|
1884
|
|
|
|
|
|
|
The Cosine window. Another name for this window is the sine window. This window is defined by |
1885
|
|
|
|
|
|
|
|
1886
|
|
|
|
|
|
|
arr, |
1887
|
|
|
|
|
|
|
|
1888
|
|
|
|
|
|
|
where the points in arr are the sin of points ranging from 0 through PI. |
1889
|
|
|
|
|
|
|
|
1890
|
|
|
|
|
|
|
|
1891
|
|
|
|
|
|
|
=cut |
1892
|
|
|
|
|
|
|
|
1893
|
|
|
|
|
|
|
|
1894
|
|
|
|
|
|
|
=head2 dpss($N,$beta) |
1895
|
|
|
|
|
|
|
|
1896
|
|
|
|
|
|
|
The Digital Prolate Spheroidal Sequence (DPSS) window. The parameter C<$beta> is the half-width of the mainlobe, measured in frequency bins. This window maximizes the power in the mainlobe for given C<$N> and C<$beta>. |
1897
|
|
|
|
|
|
|
Another name for this window is the sleppian window. |
1898
|
|
|
|
|
|
|
|
1899
|
|
|
|
|
|
|
=cut |
1900
|
|
|
|
|
|
|
|
1901
|
|
|
|
|
|
|
|
1902
|
|
|
|
|
|
|
=head2 exponential($N) |
1903
|
|
|
|
|
|
|
|
1904
|
|
|
|
|
|
|
The Exponential window. This window is defined by |
1905
|
|
|
|
|
|
|
|
1906
|
|
|
|
|
|
|
2 ** (1 - abs arr) - 1, |
1907
|
|
|
|
|
|
|
|
1908
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
1909
|
|
|
|
|
|
|
|
1910
|
|
|
|
|
|
|
|
1911
|
|
|
|
|
|
|
=cut |
1912
|
|
|
|
|
|
|
|
1913
|
|
|
|
|
|
|
|
1914
|
|
|
|
|
|
|
=head2 flattop($N) |
1915
|
|
|
|
|
|
|
|
1916
|
|
|
|
|
|
|
The flat top window. One of the Blackman-Harris family, with coefficients |
1917
|
|
|
|
|
|
|
|
1918
|
|
|
|
|
|
|
a0 = 0.21557895, a1 = 0.41663158, a2 = 0.277263158, a3 = 0.083578947, a4 = 0.006947368 |
1919
|
|
|
|
|
|
|
|
1920
|
|
|
|
|
|
|
|
1921
|
|
|
|
|
|
|
|
1922
|
|
|
|
|
|
|
=cut |
1923
|
|
|
|
|
|
|
|
1924
|
|
|
|
|
|
|
|
1925
|
|
|
|
|
|
|
=head2 gaussian($N,$beta) |
1926
|
|
|
|
|
|
|
|
1927
|
|
|
|
|
|
|
The Gaussian window. (Ref 1). Another name for this window is the Weierstrass window. This window is defined by |
1928
|
|
|
|
|
|
|
|
1929
|
|
|
|
|
|
|
exp (-0.5 * ($beta * arr )**2), |
1930
|
|
|
|
|
|
|
|
1931
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
1932
|
|
|
|
|
|
|
|
1933
|
|
|
|
|
|
|
|
1934
|
|
|
|
|
|
|
=cut |
1935
|
|
|
|
|
|
|
|
1936
|
|
|
|
|
|
|
|
1937
|
|
|
|
|
|
|
=head2 hamming($N) |
1938
|
|
|
|
|
|
|
|
1939
|
|
|
|
|
|
|
The Hamming window. (Ref 1). One of the Blackman-Harris family, with coefficients |
1940
|
|
|
|
|
|
|
|
1941
|
|
|
|
|
|
|
a0 = 0.54, a1 = 0.46 |
1942
|
|
|
|
|
|
|
|
1943
|
|
|
|
|
|
|
|
1944
|
|
|
|
|
|
|
|
1945
|
|
|
|
|
|
|
=cut |
1946
|
|
|
|
|
|
|
|
1947
|
|
|
|
|
|
|
|
1948
|
|
|
|
|
|
|
=head2 hamming_ex($N) |
1949
|
|
|
|
|
|
|
|
1950
|
|
|
|
|
|
|
The 'exact' Hamming window. (Ref 1). One of the Blackman-Harris family, with coefficients |
1951
|
|
|
|
|
|
|
|
1952
|
|
|
|
|
|
|
a0 = 0.53836, a1 = 0.46164 |
1953
|
|
|
|
|
|
|
|
1954
|
|
|
|
|
|
|
|
1955
|
|
|
|
|
|
|
|
1956
|
|
|
|
|
|
|
=cut |
1957
|
|
|
|
|
|
|
|
1958
|
|
|
|
|
|
|
|
1959
|
|
|
|
|
|
|
=head2 hamming_gen($N,$a) |
1960
|
|
|
|
|
|
|
|
1961
|
|
|
|
|
|
|
The general Hamming window. (Ref 1). One of the Blackman-Harris family, with coefficients |
1962
|
|
|
|
|
|
|
|
1963
|
|
|
|
|
|
|
a0 = $a, a1 = (1-$a) |
1964
|
|
|
|
|
|
|
|
1965
|
|
|
|
|
|
|
|
1966
|
|
|
|
|
|
|
|
1967
|
|
|
|
|
|
|
=cut |
1968
|
|
|
|
|
|
|
|
1969
|
|
|
|
|
|
|
|
1970
|
|
|
|
|
|
|
=head2 hann($N) |
1971
|
|
|
|
|
|
|
|
1972
|
|
|
|
|
|
|
The Hann window. (Ref 1). One of the Blackman-Harris family, with coefficients |
1973
|
|
|
|
|
|
|
|
1974
|
|
|
|
|
|
|
a0 = 0.5, a1 = 0.5 |
1975
|
|
|
|
|
|
|
|
1976
|
|
|
|
|
|
|
Another name for this window is the hanning window. See also L. |
1977
|
|
|
|
|
|
|
|
1978
|
|
|
|
|
|
|
|
1979
|
|
|
|
|
|
|
=cut |
1980
|
|
|
|
|
|
|
|
1981
|
|
|
|
|
|
|
|
1982
|
|
|
|
|
|
|
=head2 hann_matlab($N) |
1983
|
|
|
|
|
|
|
|
1984
|
|
|
|
|
|
|
The Hann (matlab) window. Equivalent to the Hann window of N+2 points, with the endpoints (which are both zero) removed. No periodic version of this window is defined. |
1985
|
|
|
|
|
|
|
This window is defined by |
1986
|
|
|
|
|
|
|
|
1987
|
|
|
|
|
|
|
0.5 - 0.5 * arr, |
1988
|
|
|
|
|
|
|
|
1989
|
|
|
|
|
|
|
where the points in arr are the cosine of points ranging from 2PI/($N+1) through 2PI*$N/($N+1). |
1990
|
|
|
|
|
|
|
This routine gives the same result as the routine B in Matlab. |
1991
|
|
|
|
|
|
|
See also L. |
1992
|
|
|
|
|
|
|
|
1993
|
|
|
|
|
|
|
|
1994
|
|
|
|
|
|
|
=cut |
1995
|
|
|
|
|
|
|
|
1996
|
|
|
|
|
|
|
|
1997
|
|
|
|
|
|
|
=head2 hann_poisson($N,$alpha) |
1998
|
|
|
|
|
|
|
|
1999
|
|
|
|
|
|
|
The Hann-Poisson window. (Ref 1). This window is defined by |
2000
|
|
|
|
|
|
|
|
2001
|
|
|
|
|
|
|
0.5 * (1 + arr1) * exp (-$alpha * abs arr), |
2002
|
|
|
|
|
|
|
|
2003
|
|
|
|
|
|
|
where the points in arr range from -1 through 1, and arr1 are the cos of points ranging from -PI through PI. |
2004
|
|
|
|
|
|
|
|
2005
|
|
|
|
|
|
|
|
2006
|
|
|
|
|
|
|
=cut |
2007
|
|
|
|
|
|
|
|
2008
|
|
|
|
|
|
|
|
2009
|
|
|
|
|
|
|
=head2 kaiser($N,$beta) |
2010
|
|
|
|
|
|
|
|
2011
|
|
|
|
|
|
|
The Kaiser window. (Ref 1). The parameter C<$beta> is the approximate half-width of the mainlobe, measured in frequency bins. |
2012
|
|
|
|
|
|
|
Another name for this window is the Kaiser-Bessel window. This window is defined by |
2013
|
|
|
|
|
|
|
|
2014
|
|
|
|
|
|
|
|
2015
|
|
|
|
|
|
|
barf "kaiser: PDL::GSLSF not installed" unless $HAVE_BESSEL; |
2016
|
|
|
|
|
|
|
$beta *= PI; |
2017
|
|
|
|
|
|
|
my @n = PDL::GSLSF::BESSEL::gsl_sf_bessel_In ($beta * sqrt(1 - arr **2),0); |
2018
|
|
|
|
|
|
|
my @d = PDL::GSLSF::BESSEL::gsl_sf_bessel_In($beta,0); |
2019
|
|
|
|
|
|
|
(shift @n)/(shift @d), |
2020
|
|
|
|
|
|
|
|
2021
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
2022
|
|
|
|
|
|
|
|
2023
|
|
|
|
|
|
|
|
2024
|
|
|
|
|
|
|
=cut |
2025
|
|
|
|
|
|
|
|
2026
|
|
|
|
|
|
|
|
2027
|
|
|
|
|
|
|
=head2 lanczos($N) |
2028
|
|
|
|
|
|
|
|
2029
|
|
|
|
|
|
|
The Lanczos window. Another name for this window is the sinc window. This window is defined by |
2030
|
|
|
|
|
|
|
|
2031
|
|
|
|
|
|
|
|
2032
|
|
|
|
|
|
|
my $x = PI * arr; |
2033
|
|
|
|
|
|
|
my $res = sin($x)/$x; |
2034
|
|
|
|
|
|
|
my $mid; |
2035
|
|
|
|
|
|
|
$mid = int($N/2), $res->slice($mid) .= 1 if $N % 2; |
2036
|
|
|
|
|
|
|
$res;, |
2037
|
|
|
|
|
|
|
|
2038
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
2039
|
|
|
|
|
|
|
|
2040
|
|
|
|
|
|
|
|
2041
|
|
|
|
|
|
|
=cut |
2042
|
|
|
|
|
|
|
|
2043
|
|
|
|
|
|
|
|
2044
|
|
|
|
|
|
|
=head2 nuttall($N) |
2045
|
|
|
|
|
|
|
|
2046
|
|
|
|
|
|
|
The Nuttall window. One of the Blackman-Harris family, with coefficients |
2047
|
|
|
|
|
|
|
|
2048
|
|
|
|
|
|
|
a0 = 0.3635819, a1 = 0.4891775, a2 = 0.1365995, a3 = 0.0106411 |
2049
|
|
|
|
|
|
|
|
2050
|
|
|
|
|
|
|
See also L. |
2051
|
|
|
|
|
|
|
|
2052
|
|
|
|
|
|
|
|
2053
|
|
|
|
|
|
|
=cut |
2054
|
|
|
|
|
|
|
|
2055
|
|
|
|
|
|
|
|
2056
|
|
|
|
|
|
|
=head2 nuttall1($N) |
2057
|
|
|
|
|
|
|
|
2058
|
|
|
|
|
|
|
The Nuttall (v1) window. A window referred to as the Nuttall window. One of the Blackman-Harris family, with coefficients |
2059
|
|
|
|
|
|
|
|
2060
|
|
|
|
|
|
|
a0 = 0.355768, a1 = 0.487396, a2 = 0.144232, a3 = 0.012604 |
2061
|
|
|
|
|
|
|
|
2062
|
|
|
|
|
|
|
This routine gives the same result as the routine B in Octave 3.6.2. |
2063
|
|
|
|
|
|
|
See also L. |
2064
|
|
|
|
|
|
|
|
2065
|
|
|
|
|
|
|
|
2066
|
|
|
|
|
|
|
=cut |
2067
|
|
|
|
|
|
|
|
2068
|
|
|
|
|
|
|
|
2069
|
|
|
|
|
|
|
=head2 parzen($N) |
2070
|
|
|
|
|
|
|
|
2071
|
|
|
|
|
|
|
The Parzen window. (Ref 1). Other names for this window are: Jackson, Valle-Poussin. This function disagrees with the Octave subroutine B, but agrees with Ref. 1. |
2072
|
|
|
|
|
|
|
See also L. |
2073
|
|
|
|
|
|
|
|
2074
|
|
|
|
|
|
|
|
2075
|
|
|
|
|
|
|
=cut |
2076
|
|
|
|
|
|
|
|
2077
|
|
|
|
|
|
|
|
2078
|
|
|
|
|
|
|
=head2 parzen_octave($N) |
2079
|
|
|
|
|
|
|
|
2080
|
|
|
|
|
|
|
The Parzen window. No periodic version of this window is defined. |
2081
|
|
|
|
|
|
|
This routine gives the same result as the routine B in Octave 3.6.2. |
2082
|
|
|
|
|
|
|
See also L. |
2083
|
|
|
|
|
|
|
|
2084
|
|
|
|
|
|
|
|
2085
|
|
|
|
|
|
|
=cut |
2086
|
|
|
|
|
|
|
|
2087
|
|
|
|
|
|
|
|
2088
|
|
|
|
|
|
|
=head2 poisson($N,$alpha) |
2089
|
|
|
|
|
|
|
|
2090
|
|
|
|
|
|
|
The Poisson window. (Ref 1). This window is defined by |
2091
|
|
|
|
|
|
|
|
2092
|
|
|
|
|
|
|
exp (-$alpha * abs arr), |
2093
|
|
|
|
|
|
|
|
2094
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
2095
|
|
|
|
|
|
|
|
2096
|
|
|
|
|
|
|
|
2097
|
|
|
|
|
|
|
=cut |
2098
|
|
|
|
|
|
|
|
2099
|
|
|
|
|
|
|
|
2100
|
|
|
|
|
|
|
=head2 rectangular($N) |
2101
|
|
|
|
|
|
|
|
2102
|
|
|
|
|
|
|
The Rectangular window. (Ref 1). Other names for this window are: dirichlet, boxcar. |
2103
|
|
|
|
|
|
|
|
2104
|
|
|
|
|
|
|
=cut |
2105
|
|
|
|
|
|
|
|
2106
|
|
|
|
|
|
|
|
2107
|
|
|
|
|
|
|
=head2 triangular($N) |
2108
|
|
|
|
|
|
|
|
2109
|
|
|
|
|
|
|
The Triangular window. This window is defined by |
2110
|
|
|
|
|
|
|
|
2111
|
|
|
|
|
|
|
1 - abs arr, |
2112
|
|
|
|
|
|
|
|
2113
|
|
|
|
|
|
|
where the points in arr range from -$N/($N-1) through $N/($N-1). |
2114
|
|
|
|
|
|
|
See also L. |
2115
|
|
|
|
|
|
|
|
2116
|
|
|
|
|
|
|
|
2117
|
|
|
|
|
|
|
=cut |
2118
|
|
|
|
|
|
|
|
2119
|
|
|
|
|
|
|
|
2120
|
|
|
|
|
|
|
=head2 tukey($N,$alpha) |
2121
|
|
|
|
|
|
|
|
2122
|
|
|
|
|
|
|
The Tukey window. (Ref 1). Another name for this window is the tapered cosine window. |
2123
|
|
|
|
|
|
|
|
2124
|
|
|
|
|
|
|
=cut |
2125
|
|
|
|
|
|
|
|
2126
|
|
|
|
|
|
|
|
2127
|
|
|
|
|
|
|
=head2 welch($N) |
2128
|
|
|
|
|
|
|
|
2129
|
|
|
|
|
|
|
The Welch window. (Ref 1). Other names for this window are: Riez, Bochner, Parzen, parabolic. This window is defined by |
2130
|
|
|
|
|
|
|
|
2131
|
|
|
|
|
|
|
1 - arr**2, |
2132
|
|
|
|
|
|
|
|
2133
|
|
|
|
|
|
|
where the points in arr range from -1 through 1. |
2134
|
|
|
|
|
|
|
|
2135
|
|
|
|
|
|
|
|
2136
|
|
|
|
|
|
|
=cut |
2137
|
|
|
|
|
|
|
|
2138
|
|
|
|
|
|
|
|
2139
|
|
|
|
|
|
|
# Maxima code to convert between powers of cos and multiple angles in cos |
2140
|
|
|
|
|
|
|
#grind(trigsimp(trigexpand(a0 - a1*cos(x) +a2*cos(2*x) -a3*cos(3*x) + a4*cos(4*x) -a5*cos(5*x) +a6*cos(6*x)))); |
2141
|
|
|
|
|
|
|
# |
2142
|
|
|
|
|
|
|
#32*a6*cos(x)^6-16*a5*cos(x)^5+(8*a4-48*a6)*cos(x)^4+(20*a5-4*a3)*cos(x)^3 |
2143
|
|
|
|
|
|
|
# +(18*a6-8*a4+2*a2)*cos(x)^2+(-5*a5+3*a3-a1)*cos(x)-a6+a4-a2+a0$ |
2144
|
|
|
|
|
|
|
|
2145
|
|
|
|
|
|
|
#(%i37) grind(trigsimp(trigreduce(c0 + c1*cos(x) |
2146
|
|
|
|
|
|
|
# + c2*cos(x)^2 + c3*cos(x)^3 + c4*cos(x)^4 + c5*cos(x)^5 + c6*cos(x)^6))); |
2147
|
|
|
|
|
|
|
# |
2148
|
|
|
|
|
|
|
#(c6*cos(6*x)+2*c5*cos(5*x)+(6*c6+4*c4)*cos(4*x)+(10*c5+8*c3)*cos(3*x) |
2149
|
|
|
|
|
|
|
# +(15*c6+16*c4+16*c2)*cos(2*x)+(20*c5+24*c3+32*c1)*cos(x)+10*c6+12*c4+16*c2+32*c0) /32$ |
2150
|
|
|
|
|
|
|
|
2151
|
|
|
|
|
|
|
=head1 AUXILIARY SUBROUTINES |
2152
|
|
|
|
|
|
|
|
2153
|
|
|
|
|
|
|
These subroutines are used internally, but are also available for export. |
2154
|
|
|
|
|
|
|
|
2155
|
|
|
|
|
|
|
=head2 cos_mult_to_pow |
2156
|
|
|
|
|
|
|
|
2157
|
|
|
|
|
|
|
Convert Blackman-Harris coefficients. The BH windows are usually defined via coefficients |
2158
|
|
|
|
|
|
|
for cosines of integer multiples of an argument. The same windows may be written instead |
2159
|
|
|
|
|
|
|
as terms of powers of cosines of the same argument. These may be computed faster as they |
2160
|
|
|
|
|
|
|
replace evaluation of cosines with multiplications. |
2161
|
|
|
|
|
|
|
This subroutine is used internally to implement the Blackman-Harris |
2162
|
|
|
|
|
|
|
family of windows more efficiently. |
2163
|
|
|
|
|
|
|
|
2164
|
|
|
|
|
|
|
This subroutine takes between 1 and 7 numeric arguments a0, a1, ... |
2165
|
|
|
|
|
|
|
|
2166
|
|
|
|
|
|
|
It converts the coefficients of this |
2167
|
|
|
|
|
|
|
|
2168
|
|
|
|
|
|
|
a0 - a1 cos(arg) + a2 cos( 2 * arg) - a3 cos( 3 * arg) + ... |
2169
|
|
|
|
|
|
|
|
2170
|
|
|
|
|
|
|
To the cofficients of this |
2171
|
|
|
|
|
|
|
|
2172
|
|
|
|
|
|
|
c0 + c1 cos(arg) + c2 cos(arg)**2 + c3 cos(arg)**3 + ... |
2173
|
|
|
|
|
|
|
|
2174
|
|
|
|
|
|
|
=head2 cos_pow_to_mult |
2175
|
|
|
|
|
|
|
|
2176
|
|
|
|
|
|
|
This function is the inverse of L. |
2177
|
|
|
|
|
|
|
|
2178
|
|
|
|
|
|
|
This subroutine takes between 1 and 7 numeric arguments c0, c1, ... |
2179
|
|
|
|
|
|
|
|
2180
|
|
|
|
|
|
|
It converts the coefficients of this |
2181
|
|
|
|
|
|
|
|
2182
|
|
|
|
|
|
|
c0 + c1 cos(arg) + c2 cos(arg)**2 + c3 cos(arg)**3 + ... |
2183
|
|
|
|
|
|
|
|
2184
|
|
|
|
|
|
|
To the cofficients of this |
2185
|
|
|
|
|
|
|
|
2186
|
|
|
|
|
|
|
a0 - a1 cos(arg) + a2 cos( 2 * arg) - a3 cos( 3 * arg) + ... |
2187
|
|
|
|
|
|
|
|
2188
|
|
|
|
|
|
|
=cut |
2189
|
|
|
|
|
|
|
|
2190
|
|
|
|
|
|
|
sub cos_pow_to_mult { |
2191
|
0
|
|
|
0
|
1
|
0
|
my( @cin ) = @_; |
2192
|
0
|
0
|
|
|
|
0
|
barf "cos_pow_to_mult: number of args not less than 8." if @cin > 7; |
2193
|
0
|
|
|
|
|
0
|
my $ex = 7 - @cin; |
2194
|
0
|
|
|
|
|
0
|
my @c = (@cin, (0) x $ex); |
2195
|
0
|
|
|
|
|
0
|
my (@as) = ( |
2196
|
|
|
|
|
|
|
10*$c[6]+12*$c[4]+16*$c[2]+32*$c[0], 20*$c[5]+24*$c[3]+32*$c[1], |
2197
|
|
|
|
|
|
|
15*$c[6]+16*$c[4]+16*$c[2], 10*$c[5]+8*$c[3], 6*$c[6]+4*$c[4], 2*$c[5], $c[6]); |
2198
|
0
|
|
|
|
|
0
|
foreach (1..$ex) {pop (@as)} |
|
0
|
|
|
|
|
0
|
|
2199
|
0
|
|
|
|
|
0
|
my $sign = -1; |
2200
|
0
|
|
|
|
|
0
|
foreach (@as) { $_ /= (-$sign*32); $sign *= -1 } |
|
0
|
|
|
|
|
0
|
|
|
0
|
|
|
|
|
0
|
|
2201
|
0
|
|
|
|
|
0
|
@as; |
2202
|
|
|
|
|
|
|
} |
2203
|
|
|
|
|
|
|
|
2204
|
|
|
|
|
|
|
=head2 chebpoly |
2205
|
|
|
|
|
|
|
|
2206
|
|
|
|
|
|
|
=for usage |
2207
|
|
|
|
|
|
|
|
2208
|
|
|
|
|
|
|
chebpoly($n,$x) |
2209
|
|
|
|
|
|
|
|
2210
|
|
|
|
|
|
|
=for ref |
2211
|
|
|
|
|
|
|
|
2212
|
|
|
|
|
|
|
Returns the value of the C<$n>-th order Chebyshev polynomial at point C<$x>. |
2213
|
|
|
|
|
|
|
$n and $x may be scalar numbers, pdl's, or array refs. However, |
2214
|
|
|
|
|
|
|
at least one of $n and $x must be a scalar number. |
2215
|
|
|
|
|
|
|
|
2216
|
|
|
|
|
|
|
All mixtures of pdls and scalars could be handled much more |
2217
|
|
|
|
|
|
|
easily as a PP routine. But, at this point PDL::DSP::Windows |
2218
|
|
|
|
|
|
|
is pure perl/pdl, requiring no C/Fortran compiler. |
2219
|
|
|
|
|
|
|
|
2220
|
|
|
|
|
|
|
=cut |
2221
|
|
|
|
|
|
|
|
2222
|
|
|
|
|
|
|
sub chebpoly { |
2223
|
5
|
50
|
|
5
|
1
|
2410
|
barf 'chebpoly: Two arguments expected. Got ' .scalar(@_) ."\n" unless @_==2; |
2224
|
5
|
|
|
|
|
11
|
my ($n,$x) = @_; |
2225
|
5
|
100
|
|
|
|
13
|
if (ref($x)) { |
2226
|
4
|
|
|
|
|
14
|
$x = topdl($x); |
2227
|
4
|
50
|
|
|
|
69
|
barf "chebpoly: neither $n nor $x is a scalar number" if ref($n); |
2228
|
4
|
|
|
|
|
13
|
my $tn = zeroes($x); |
2229
|
4
|
|
|
|
|
207
|
my ($ind1,$ind2); |
2230
|
4
|
|
|
|
|
12
|
($ind1,$ind2) = which_both(abs($x) <= 1); |
2231
|
4
|
|
|
|
|
384
|
$tn->index($ind1) .= cos($n*(acos($x->index($ind1)))); |
2232
|
4
|
|
|
|
|
285
|
$tn->index($ind2) .= cosh($n*(acosh($x->index($ind2)))); |
2233
|
4
|
|
|
|
|
197
|
return $tn; |
2234
|
|
|
|
|
|
|
} |
2235
|
|
|
|
|
|
|
else { |
2236
|
1
|
50
|
|
|
|
5
|
$n = topdl($n) if ref($n); |
2237
|
1
|
50
|
|
|
|
5
|
return cos($n*(acos($x))) if abs($x) <= 1; |
2238
|
1
|
|
|
|
|
13
|
return cosh($n*(acosh($x))); |
2239
|
|
|
|
|
|
|
} |
2240
|
|
|
|
|
|
|
} |
2241
|
|
|
|
|
|
|
|
2242
|
|
|
|
|
|
|
|
2243
|
|
|
|
|
|
|
sub cos_mult_to_pow { |
2244
|
0
|
|
|
0
|
1
|
|
my( @ain ) = @_; |
2245
|
0
|
0
|
|
|
|
|
barf("cos_mult_to_pow: number of args not less than 8.") if @ain > 7; |
2246
|
0
|
|
|
|
|
|
my $ex = 7 - @ain; |
2247
|
0
|
|
|
|
|
|
my @a = (@ain, (0) x $ex); |
2248
|
0
|
|
|
|
|
|
my (@cs) = ( |
2249
|
|
|
|
|
|
|
-$a[6]+$a[4]-$a[2]+$a[0], -5*$a[5]+3*$a[3]-$a[1], 18*$a[6]-8*$a[4]+2*$a[2], 20*$a[5]-4*$a[3], |
2250
|
|
|
|
|
|
|
8*$a[4]-48*$a[6], -16*$a[5], 32*$a[6]); |
2251
|
0
|
|
|
|
|
|
foreach (1..$ex) {pop (@cs)} |
|
0
|
|
|
|
|
|
|
2252
|
0
|
|
|
|
|
|
@cs; |
2253
|
|
|
|
|
|
|
} |
2254
|
|
|
|
|
|
|
|
2255
|
|
|
|
|
|
|
=head1 REFERENCES |
2256
|
|
|
|
|
|
|
|
2257
|
|
|
|
|
|
|
=over |
2258
|
|
|
|
|
|
|
|
2259
|
|
|
|
|
|
|
=item 1 |
2260
|
|
|
|
|
|
|
|
2261
|
|
|
|
|
|
|
Harris, F.J. C, |
2262
|
|
|
|
|
|
|
I, 1978, vol 66, pp 51-83. |
2263
|
|
|
|
|
|
|
|
2264
|
|
|
|
|
|
|
=item 2 |
2265
|
|
|
|
|
|
|
|
2266
|
|
|
|
|
|
|
Nuttall, A.H. C, I, |
2267
|
|
|
|
|
|
|
1981, vol. ASSP-29, pp. 84-91. |
2268
|
|
|
|
|
|
|
|
2269
|
|
|
|
|
|
|
=back |
2270
|
|
|
|
|
|
|
|
2271
|
|
|
|
|
|
|
=head1 AUTHOR |
2272
|
|
|
|
|
|
|
|
2273
|
|
|
|
|
|
|
John Lapeyre, C<< >> |
2274
|
|
|
|
|
|
|
|
2275
|
|
|
|
|
|
|
=head1 ACKNOWLEDGMENTS |
2276
|
|
|
|
|
|
|
|
2277
|
|
|
|
|
|
|
For study and comparison, the author used documents or output from: |
2278
|
|
|
|
|
|
|
Thomas Cokelaer's spectral analysis software; Julius O Smith III's |
2279
|
|
|
|
|
|
|
Spectral Audio Signal Processing web pages; André Carezia's |
2280
|
|
|
|
|
|
|
chebwin.m Octave code; Other code in the Octave signal package. |
2281
|
|
|
|
|
|
|
|
2282
|
|
|
|
|
|
|
=head1 LICENSE AND COPYRIGHT |
2283
|
|
|
|
|
|
|
|
2284
|
|
|
|
|
|
|
Copyright 2012 John Lapeyre. |
2285
|
|
|
|
|
|
|
|
2286
|
|
|
|
|
|
|
This program is free software; you can redistribute it and/or modify it |
2287
|
|
|
|
|
|
|
under the terms of either: the GNU General Public License as published |
2288
|
|
|
|
|
|
|
by the Free Software Foundation; or the Artistic License. |
2289
|
|
|
|
|
|
|
|
2290
|
|
|
|
|
|
|
See http://dev.perl.org/licenses/ for more information. |
2291
|
|
|
|
|
|
|
|
2292
|
|
|
|
|
|
|
This software is neither licensed nor distributed by The MathWorks, Inc., |
2293
|
|
|
|
|
|
|
maker and liscensor of MATLAB. |
2294
|
|
|
|
|
|
|
|
2295
|
|
|
|
|
|
|
=cut |
2296
|
|
|
|
|
|
|
|
2297
|
|
|
|
|
|
|
1; # End of PDL::DSP::Windows.pm\n"; |
2298
|
|
|
|
|
|
|
|