1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
#1 "fftw3SD.ml"
open Bigarray
open Fftw3_utils
type float_elt = Bigarray.float64_elt
type complex_elt = Bigarray.complex64_elt
let float = Bigarray.float64
let complex = Bigarray.complex64
type 'a fftw_plan
type c2c
type r2c
type c2r
type r2r
type dir = Forward | Backward
type measure = Estimate | Measure | Patient | Exhaustive
type r2r_kind =
| R2HC | HC2R | DHT
| REDFT00 | REDFT01 | REDFT10 | REDFT11
| RODFT00 | RODFT01 | RODFT10 | RODFT11
exception Failure of string
let is_c_layout m =
(Genarray.layout m = (Obj.magic c_layout : 'a layout))
(** External declarations
***********************************************************************)
type 'l complex_array = (Complex.t, Bigarray.complex64_elt, 'l) Genarray.t
type 'l float_array = (float, Bigarray.float64_elt, 'l) Genarray.t
external fftw_exec : 'a fftw_plan -> unit = "fftw_ocaml_execute" [@@noalloc]
external exec_dft : c2c fftw_plan -> 'l complex_array -> 'l complex_array
-> unit = "fftw_ocaml_execute_dft" [@@noalloc]
external exec_split_dft : c2c fftw_plan -> 'l float_array -> 'l float_array ->
'l float_array -> 'l float_array -> unit
= "fftw_ocaml_execute_split_dft" [@@noalloc]
external exec_dft_r2c : r2c fftw_plan -> 'l float_array -> 'l complex_array
-> unit
= "fftw_ocaml_execute_dft_r2c" [@@noalloc]
external exec_split_dft_r2c : r2c fftw_plan -> 'l float_array ->
'l float_array -> 'l float_array -> unit
= "fftw_ocaml_execute_split_dft_r2c" [@@noalloc]
external exec_dft_c2r : c2r fftw_plan -> 'l complex_array -> 'l float_array
-> unit
= "fftw_ocaml_execute_dft_c2r" [@@noalloc]
external exec_split_dft_c2r : c2r fftw_plan -> 'l float_array -> 'l float_array
-> 'l float_array -> unit
= "fftw_ocaml_execute_split_dft_c2r" [@@noalloc]
external exec_r2r : r2r fftw_plan -> 'l float_array -> 'l float_array
-> unit
= "fftw_ocaml_execute_r2r" [@@noalloc]
external guru_dft :
'l complex_array ->
'l complex_array ->
int ->
int ->
int ->
int ->
int array ->
int array ->
int array ->
int array ->
int array ->
int array
-> c2c fftw_plan
= "fftw_ocaml_guru_dft_bc" "fftw_ocaml_guru_dft"
external guru_r2c :
'l float_array ->
'l complex_array ->
int ->
int ->
int ->
int array ->
int array ->
int array ->
int array ->
int array ->
int array
-> r2c fftw_plan
= "fftw_ocaml_guru_r2c_bc" "fftw_ocaml_guru_r2c"
external guru_c2r :
'l complex_array ->
'l float_array ->
int ->
int ->
int ->
int array ->
int array ->
int array ->
int array ->
int array ->
int array
-> c2r fftw_plan
= "fftw_ocaml_guru_c2r_bc" "fftw_ocaml_guru_c2r"
external guru_r2r :
'l float_array ->
'l float_array ->
r2r_kind array ->
int ->
int ->
int ->
int array ->
int array ->
int array ->
int array ->
int array ->
int array
-> r2r fftw_plan
= "fftw_ocaml_guru_r2r_bc" "fftw_ocaml_guru_r2r"
(** Plans on the OCaml side
***********************************************************************)
type genarray
external genarray : (_,_,_) Genarray.t -> genarray = "%identity"
type 'a plan = {
plan: 'a fftw_plan;
i : genarray;
offseto : int;
strideo : int array;
no : int array;
o : genarray;
}
let sign_of_dir = function
| Forward -> -1
| Backward -> 1
let flags meas unaligned ~destroy_input : int =
let f = match meas with
| Measure -> 0
| Exhaustive -> 8
| Patient -> 32
| Estimate -> 64 in
let f = if unaligned then f lor 2 else f in
if destroy_input then f lor 1 else f lor 16
(** {2 Execution of plans}
***********************************************************************)
let exec p =
fftw_exec p.plan
module Guru = struct
let dft plan i o =
exec_dft plan i o
let split_dft plan ri ii ro io =
exec_split_dft plan ri ii ro io
end
(** {2 Creating plans}
***********************************************************************)
module Genarray = struct
external create: ('a, 'b) Bigarray.kind -> 'c Bigarray.layout ->
int array -> ('a, 'b, 'c) Bigarray.Genarray.t
= "fftw3_ocaml_ba_create"
type 'l complex_array = (Complex.t, Bigarray.complex64_elt, 'l) Genarray.t
type 'l float_array = (float, Bigarray.float64_elt, 'l) Genarray.t
type coord = int array
let apply name mk_plan hm_n hmi ?ni ofsi inci i hmo ?no ofso inco o
~logical_dims =
let make offseti offseto n stridei strideo hm_ni hm_stridei hm_strideo =
let p = (mk_plan offseti offseto n stridei strideo
hm_ni hm_stridei hm_strideo) in
{ plan = p;
i = genarray i;
offseto = offseto;
strideo = strideo;
no = n;
o = genarray o;
} in
(if is_c_layout i then Fftw3_geomC.apply else Fftw3_geomF.apply)
name make hm_n hmi ?ni ofsi inci i hmo ?no ofso inco o ~logical_dims
let dft_name = "Fftw3.D.Genarray.dft"
let dft dir ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?inci (i: 'l complex_array)
?(howmanyo=[]) ?no ?ofso ?inco (o: 'l complex_array) =
apply dft_name ~logical_dims:Geom.logical_c2c
(guru_dft i o (sign_of_dir dir) (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci i howmanyo ?no ofso inco o
let r2c_name = "Fftw3.D.Genarray.r2c"
let r2c ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?inci (i: 'l float_array)
?(howmanyo=[]) ?no ?ofso ?inco (o: 'l complex_array) =
apply r2c_name ~logical_dims:Geom.logical_r2c
(guru_r2c i o (flags meas unaligned ~destroy_input))
howmany_n howmanyi ofsi ?ni inci i howmanyo ?no ofso inco o
let c2r_name = "Fftw3.D.Genarray.c2r"
let c2r ?(meas=Measure)
?(destroy_input=true) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?inci (i: 'l complex_array)
?(howmanyo=[]) ?no ?ofso ?inco (o: 'l float_array) =
apply c2r_name ~logical_dims:Geom.logical_c2r
(guru_c2r i o (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci i howmanyo ?no ofso inco o
let r2r_name = "Fftw3.D.Genarray.r2r"
let r2r kind ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?inci (i: 'l float_array)
?(howmanyo=[]) ?no ?ofso ?inco (o: 'l float_array) =
apply r2r_name ~logical_dims:Geom.logical_r2r
(guru_r2r i o kind (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci i howmanyo ?no ofso inco o
end
module Array1 = struct
external array1_of_ba : ('a,'b,'c) Bigarray.Genarray.t -> ('a,'b,'c) Array1.t
= "%identity"
let create kind layout dim =
array1_of_ba(Genarray.create kind layout [|dim|])
let of_array kind layout data =
let ba = create kind layout (Array.length data) in
let ofs = if layout = (Obj.magic c_layout : 'a layout) then 0 else 1 in
for i = 0 to Array.length data - 1 do ba.{i + ofs} <- data.(i) done;
ba
type 'l complex_array = (Complex.t, Bigarray.complex64_elt, 'l) Array1.t
type 'l float_array = (float, Bigarray.float64_elt, 'l) Array1.t
let apply name make_plan hm_n hmi ?ni ofsi inci i hmo ?no ofso inco o
~logical_dims =
let hmi = List.map (fun v -> [| v |]) hmi in
let ni = option_map (fun n -> [| n |]) ni in
let ofsi = option_map (fun n -> [| n |]) ofsi in
let inci = Some [| inci |] in
let hmo = List.map (fun v -> [| v |]) hmo in
let no = option_map (fun n -> [| n |]) no in
let ofso = option_map (fun n -> [| n |]) ofso in
let inco = Some [| inco |] in
Genarray.apply name make_plan
hm_n hmi ?ni ofsi inci i hmo ?no ofso inco o ~logical_dims
let dft_name = "Fftw3.D.Array1.dft"
let dft dir ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=1) (i: 'l complex_array)
?(howmanyo=[]) ?no ?ofso ?(inco=1) (o: 'l complex_array) =
let gi = genarray_of_array1 i
and go = genarray_of_array1 o in
apply dft_name ~logical_dims:Geom.logical_c2c
(guru_dft gi go (sign_of_dir dir) (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let r2c_name = "Fftw3.D.Array1.r2c"
let r2c ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=1) (i: 'l float_array)
?(howmanyo=[]) ?no ?ofso ?(inco=1) (o: 'l complex_array) =
let gi = genarray_of_array1 i
and go = genarray_of_array1 o in
apply r2c_name ~logical_dims:Geom.logical_r2c
(guru_r2c gi go (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let c2r_name = "Fftw3.D.Array1.c2r"
let c2r ?(meas=Measure)
?(destroy_input=true) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=1) (i: 'l complex_array)
?(howmanyo=[]) ?no ?ofso ?(inco=1) (o: 'l float_array) =
let gi = genarray_of_array1 i
and go = genarray_of_array1 o in
apply c2r_name ~logical_dims:Geom.logical_c2r
(guru_c2r gi go (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let r2r_name = "Fftw3.D.Array1.r2r"
let r2r kind ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=1) (i: 'l float_array)
?(howmanyo=[]) ?no ?ofso ?(inco=1) (o: 'l float_array) =
let gi = genarray_of_array1 i
and go = genarray_of_array1 o in
let kind = [| kind |] in
apply r2r_name ~logical_dims:Geom.logical_r2r
(guru_r2r gi go kind (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
end
module Array2 = struct
external array2_of_ba : ('a,'b,'c) Bigarray.Genarray.t -> ('a,'b,'c) Array2.t
= "%identity"
let create kind layout dim1 dim2 =
array2_of_ba(Genarray.create kind layout [|dim1; dim2|])
type 'l complex_array = (Complex.t, Bigarray.complex64_elt, 'l) Array2.t
type 'l float_array = (float, Bigarray.float64_elt, 'l) Array2.t
type coord = int * int
let apply name make_plan hm_n hmi ?ni ofsi (inci1,inci2) i
hmo ?no ofso (inco1,inco2) o ~logical_dims =
let hmi = List.map (fun (d1,d2) -> [| d1; d2 |]) hmi in
let ni = option_map (fun (n1,n2) -> [| n1; n2 |]) ni in
let ofsi = option_map (fun (n1,n2) -> [| n1; n2 |]) ofsi in
let inci = Some [| inci1; inci2 |] in
let hmo = List.map (fun (d1,d2) -> [| d1; d2 |]) hmo in
let no = option_map (fun (n1,n2) -> [| n1; n2 |]) no in
let ofso = option_map (fun (n1,n2) -> [| n1; n2 |]) ofso in
let inco = Some [| inco1; inco2 |] in
Genarray.apply name make_plan
hm_n hmi ?ni ofsi inci i hmo ?no ofso inco o ~logical_dims
let dft_name = "Fftw3.D.Array2.dft"
let dft dir ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=(1,1)) (i: 'l complex_array)
?(howmanyo=[]) ?no ?ofso ?(inco=(1,1)) (o: 'l complex_array) =
let gi = genarray_of_array2 i
and go = genarray_of_array2 o in
apply dft_name ~logical_dims:Geom.logical_c2c
(guru_dft gi go (sign_of_dir dir) (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let r2c_name = "Fftw3.D.Array2.r2c"
let r2c ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=(1,1)) (i: 'l float_array)
?(howmanyo=[]) ?no ?ofso ?(inco=(1,1)) (o: 'l complex_array) =
let gi = genarray_of_array2 i
and go = genarray_of_array2 o in
apply r2c_name ~logical_dims:Geom.logical_r2c
(guru_r2c gi go (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let c2r_name = "Fftw3.D.Array2.c2r"
let c2r ?(meas=Measure)
?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=(1,1)) (i: 'l complex_array)
?(howmanyo=[]) ?no ?ofso ?(inco=(1,1)) (o: 'l float_array) =
let gi = genarray_of_array2 i
and go = genarray_of_array2 o in
apply c2r_name ~logical_dims:Geom.logical_c2r
(guru_c2r gi go (flags meas unaligned ~destroy_input:true))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let r2r_name = "Fftw3.D.Array2.r2r"
let r2r (kind1,kind2) ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=(1,1)) (i: 'l float_array)
?(howmanyo=[]) ?no ?ofso ?(inco=(1,1)) (o: 'l float_array) =
let gi = genarray_of_array2 i
and go = genarray_of_array2 o in
let kind = [| kind1; kind2 |] in
apply r2r_name ~logical_dims:Geom.logical_r2r
(guru_r2r gi go kind (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
end
module Array3 = struct
external array3_of_ba : ('a,'b,'c) Bigarray.Genarray.t -> ('a,'b,'c) Array3.t
= "%identity"
let create kind layout dim1 dim2 dim3 =
array3_of_ba(Genarray.create kind layout [|dim1; dim2; dim3|])
type 'l complex_array = (Complex.t, Bigarray.complex64_elt, 'l) Array3.t
type 'l float_array = (float, Bigarray.float64_elt, 'l) Array3.t
type coord = int * int * int
let apply name make_plan hm_n hmi ?ni ofsi (inci1,inci2,inci3) i
hmo ?no ofso (inco1,inco2,inco3) o ~logical_dims =
let hmi = List.map (fun (d1,d2,d3) -> [| d1; d2; d3 |]) hmi in
let ni = option_map (fun (n1,n2,n3) -> [| n1; n2; n3 |]) ni in
let ofsi = option_map (fun (n1,n2,n3) -> [| n1; n2; n3 |]) ofsi in
let inci = Some [| inci1; inci2; inci3 |] in
let hmo = List.map (fun (d1,d2,d3) -> [| d1; d2; d3 |]) hmo in
let no = option_map (fun (n1,n2,n3) -> [| n1; n2; n3 |]) no in
let ofso = option_map (fun (n1,n2,n3) -> [| n1; n2; n3 |]) ofso in
let inco = Some [| inco1; inco2; inco3 |] in
Genarray.apply name make_plan
hm_n hmi ?ni ofsi inci i hmo ?no ofso inco o ~logical_dims
let dft_name = "Fftw3.D.Array3.dft"
let dft dir ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=(1,1,1)) (i: 'l complex_array)
?(howmanyo=[]) ?no ?ofso ?(inco=(1,1,1)) (o: 'l complex_array) =
let gi = genarray_of_array3 i
and go = genarray_of_array3 o in
apply dft_name ~logical_dims:Geom.logical_c2c
(guru_dft gi go (sign_of_dir dir) (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let r2c_name = "Fftw3.D.Array3.r2c"
let r2c ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=(1,1,1)) (i: 'l float_array)
?(howmanyo=[]) ?no ?ofso ?(inco=(1,1,1)) (o: 'l complex_array) =
let gi = genarray_of_array3 i
and go = genarray_of_array3 o in
apply r2c_name ~logical_dims:Geom.logical_r2c
(guru_r2c gi go (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let c2r_name = "Fftw3.D.Array3.c2r"
let c2r ?(meas=Measure)
?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=(1,1,1)) (i: 'l complex_array)
?(howmanyo=[]) ?no ?ofso ?(inco=(1,1,1)) (o: 'l float_array) =
let gi = genarray_of_array3 i
and go = genarray_of_array3 o in
apply c2r_name ~logical_dims:Geom.logical_c2r
(guru_c2r gi go (flags meas unaligned ~destroy_input:true))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
let r2r_name = "Fftw3.D.Array3.r2r"
let r2r (kind1,kind2,kind3) ?(meas=Measure)
?(destroy_input=false) ?(unaligned=false) ?(howmany_n=[| |])
?(howmanyi=[]) ?ni ?ofsi ?(inci=(1,1,1)) (i: 'l float_array)
?(howmanyo=[]) ?no ?ofso ?(inco=(1,1,1)) (o: 'l float_array) =
let gi = genarray_of_array3 i
and go = genarray_of_array3 o in
let kind = [| kind1; kind2; kind3 |] in
apply r2r_name ~logical_dims:Geom.logical_r2r
(guru_r2r gi go kind (flags meas unaligned ~destroy_input))
howmany_n howmanyi ?ni ofsi inci gi howmanyo ?no ofso inco go
end