77 [scicloj.tableplot.v1.plotly :as plotly]
88 [tech.v3.tensor :as tensor]
99 [tablecloth.column.api :as tcc]
10- [fastmath.stats :as stats]))
11-
10+ [fastmath.stats :as stats]
11+ [tech.v3.datatype :as dtype]))
1212
1313; ; ## Data source
1414; ;
@@ -86,7 +86,6 @@ wav-format
8686 (aset-short audio-shorts i (.getShort bb 0 )))
8787 audio-shorts))
8888
89-
9089(def wav-shorts
9190 (with-open [wav-stream (io/input-stream violin-file-path)]
9291 (audio-data wav-stream)))
@@ -95,29 +94,56 @@ wav-format
9594 ; ; one of the two stereo channels
9695 (dfn// (->> wav-shorts
9796 (partition 2 )
98- (map first))
97+ (map first)
98+ dtype/->double-array)
9999 32768.0 ))
100100
101- ; ; ## Listening to the data back as audio
101+ (def wav-ds
102+ (tc/dataset
103+ {:time (dfn// (range )
104+ (double sample-rate))
105+ :sample wav-samples}))
102106
107+ wav-ds
103108
109+ (-> wav-ds
110+ (tc/select-rows #(<= 0 (:time %) 1 ))
111+ (plotly/layer-line {:=x :time
112+ :=y :sample }))
113+
114+ (-> wav-ds
115+ (tc/select-rows #(<= 1 (:time %) 1.06 ))
116+ (plotly/layer-line {:=x :time
117+ :=y :sample }))
118+
119+ ; ; ## Listening to the data back as audio
104120
105121^kind/audio
106122{:samples wav-samples
107123 :sample-rate sample-rate}
108124
125+ ^kind/audio
126+ {:samples (->> (-> wav-ds
127+ (tc/select-rows #(<= 1 (:time %) 1.06 ))
128+ :sample )
129+ (repeat 30 )
130+ (apply concat))
131+ :sample-rate sample-rate}
132+
109133; ; ## Computing the power spectrum of the data
110134
111135(import 'com.github.psambit9791.jdsp.transform.DiscreteFourier)
112136
113137(count wav-samples)
114138
115- (def first-part
116- (take sample-rate ; just a second
117- wav-samples))
139+ (def some-part
140+ (->> wav-samples
141+ (drop sample-rate) ; drop a second
142+ (take (* sample-rate 0.06 )) ; take 6% of a second
143+ ))
118144
119145(def power-spectrum
120- (let [fft (-> first -part
146+ (let [fft (-> some -part
121147 double-array
122148 DiscreteFourier.)]
123149 (.transform fft)
@@ -126,9 +152,12 @@ wav-format
126152
127153(count power-spectrum)
128154
155+ (def Nyquist-freq (* 0.5 sample-rate))
156+
129157(def power-spectrum-ds
130158 (tc/dataset
131- {:freq (range )
159+ {:freq (dfn// (range )
160+ (/ (count power-spectrum) Nyquist-freq))
132161 :power power-spectrum}))
133162
134163(-> power-spectrum-ds
@@ -145,7 +174,8 @@ wav-format
145174 (let [peaks (-> power-spectrum
146175 FindPeak.
147176 .detectPeaks)]
148- (-> (tc/dataset {:freq (.getPeaks peaks)
177+ (-> (tc/dataset {:freq (dfn// (.getPeaks peaks)
178+ (/ (count power-spectrum) Nyquist-freq ))
149179 :power (.getHeights peaks)
150180 :prominence (.getProminence peaks)})
151181 (tc/order-by [:prominence ] :desc )
@@ -161,6 +191,7 @@ wav-format
161191; ; ## Synthesizing from the peaks
162192
163193(require '[clojure.math :as math])
194+ (require '[tech.v3.datatype :as dtype])
164195
165196(def components-dataset
166197 (let [duration 1
@@ -191,130 +222,118 @@ wav-format
191222 :=y :value
192223 :=color :$column}))
193224
225+ (defn normalize [values]
226+ (tcc// values
227+ (tcc/reduce-max (tcc/abs values))))
194228
195- (def violin -dataset
229+ (def combined -dataset
196230 (-> components-dataset
197- (tc/add-column :violin
231+ (tc/add-column :combined
198232 (fn [ds]
199233 (-> ds
200234 (tc/drop-columns [:time ])
201235 vals
202- (->> (apply tcc/+)))))))
236+ (->> (apply tcc/+))
237+ normalize)))))
203238
204- (-> violin -dataset
239+ (-> combined -dataset
205240 (tc/head 200 )
206241 (plotly/layer-line {:=x :time
207- :=y :violin }))
242+ :=y :combined }))
208243
209244(defn audio [samples]
210245 (with-meta
211246 {:samples samples
212247 :sample-rate sample-rate}
213248 {:kind/audio true }))
214249
215- (defn normalize [values]
216- (tcc// values
217- (tcc/reduce-max (tcc/abs values))))
250+ (-> some-part
251+ audio)
218252
219- (-> violin-dataset
220- :violin
221- normalize
253+ (-> combined-dataset
254+ :combined
222255 audio)
223256
224- ; ; ## Power spectrum by part
225-
226- (require '[tech.v3.tensor :as tensor])
227-
228- (def parts
229- (let [n-parts 25 ]
230- (-> wav-samples
231- tensor/->tensor
232- (tensor/reshape [n-parts
233- (/ (count wav-samples) n-parts)])
234- (tensor/slice 1 ))))
235-
236-
237- (def spectra
238- (->> parts
239- (pmap (fn [part]
240- (let [fft (-> part
241- double-array
242- DiscreteFourier.)]
243- (.transform fft)
244- ; ; Get magnitude spectrum (power)
245- (.getMagnitude fft true ))))))
246-
247- (def spectra-ds
248- (->> spectra
249- (map (fn [s]
250- (tc/dataset
251- {:freq (range )
252- :power s})))))
253-
254- (->> spectra-ds
255- (map (fn [ds]
256- (-> ds
257- (plotly/base {:=height 300
258- :=width 900 })
259- (plotly/layer-line {:=x :freq
260- :=y :power })
261- plotly/plot
262- (assoc-in [:layout :xaxis :range ] [0 1000 ])))))
263-
264- ; ; ## Peaks by part
265-
266- (def all-peaks-ds
267- (->> spectra
268- (map (fn [spectrum]
269- (let [peaks (-> spectrum
270- FindPeak.
271- .detectPeaks)]
272- (-> (tc/dataset {:freq (.getPeaks peaks)
273- :power (.getHeights peaks)
274- :prominence (.getProminence peaks)})
275- (tc/order-by [:prominence ] :desc )
276- (tc/head n-peaks)))))))
277-
278- (map (fn [spectrum-ds peaks-ds]
279- (-> spectrum-ds
280- (plotly/base {:=x :freq
281- :=y :power })
282- plotly/layer-line
283- (plotly/layer-point {:=dataset peaks-ds})))
284- spectra-ds
285- all-peaks-ds)
257+ ; ; ## Spectogram (WIP)
258+
259+ (import 'com.github.psambit9791.jdsp.windows.Hanning)
260+ (require '[tech.v3.dataset.tensor :as ds-tensor])
261+
262+ (let [window-size 0.03 ; seconds
263+ window-samples (* sample-rate window-size)
264+
265+ ; ; Create Hanning window to reduce spectral leakage
266+ hanning (-> window-samples
267+ Hanning.
268+ .getWindow
269+ dtype/as-reader)
270+
271+ ; ; Overlap parameters
272+ overlap-fraction 0.05
273+ hop (int (* overlap-fraction window-samples))
274+
275+ n-samples (count wav-samples)
276+
277+ ; ; Calculate how many windows we can extract
278+ n-windows (-> n-samples
279+ (- window-samples)
280+ (quot hop))
281+
282+ ; ; Create sliding windows: [time × windows × channels]
283+ windows-shape [window-samples n-windows]
284+ windows (tensor/compute-tensor
285+ windows-shape
286+ (fn [i j]
287+ ; ; Extract sample at time i, window j, channel k
288+ (wav-samples (+ (* j hop) i)))
289+ :float32 )
290+
291+ ; ; Apply Hanning window to each segment
292+ hanninged-windows (tensor/compute-tensor
293+ windows-shape
294+ (fn [i j]
295+ (* (windows i j)
296+ (hanning i)))
297+ :float32 )
298+
299+ ; ; ;; Compute power spectrum for each window
300+ spectogram (-> hanninged-windows
301+ ; ; Rearrange to [windows × time]
302+ (tensor/transpose [1 0 ])
303+ (tensor/slice 1 )
304+ ; ; Apply FFT to each window
305+ (->> (pmap (fn [window]
306+ (let [fft (-> window
307+ dtype/as-reader
308+ double-array
309+ DiscreteFourier.)]
310+ (.transform fft)
311+ ; ; Get magnitude spectrum (power)
312+ (.getMagnitude fft true )))))
313+ tensor/->tensor
314+ ; ; Reshape to [windows × frequencies]
315+ (#(tensor/reshape %
316+ [n-windows
317+ (-> % first count)]))
318+ ; ; Transpose to [frequencies × windows] for plotting
319+ (tensor/transpose [1 0 ]))]
320+ ; ; (-> windows
321+ ; ; (dfn/* 100)
322+ ; ; plotly/imshow)
323+
324+ ; ; (-> hanninged-windows
325+ ; ; (dfn/* 100)
326+ ; ; plotly/imshow)
327+
328+ (-> spectogram
329+ (#(tensor/broadcast % (cons 3 (dtype/shape %))))
330+ (dfn/* 50 )
331+ (tensor/transpose [1 2 0 ])
332+ plotly/imshow))
333+
334+
335+
336+
286337
287- ; ; ## Synthesizing from the peaks
288338
289- (def combined-violin-ds
290- (->> all-peaks-ds
291- (map (fn [peaks-ds]
292- (let [duration (/ 3.0 25 )
293- num-samples (* duration sample-rate)
294- time (dtype/make-reader :float32
295- num-samples
296- (/ idx sample-rate))]
297- (-> (->> (tc/rows peaks-ds :as-maps )
298- ; ; For each peak, generate a sine wave
299- (map-indexed (fn [i {:keys [freq power]}]
300- [(str " wave" i)
301- (-> time
302- (dfn/* (* 2 math/PI freq))
303- dfn/sin
304- (dfn/* power))]))
305- (into {:time time})
306- tc/dataset)
307- (tc/add-column :violin
308- (fn [ds]
309- (-> ds
310- (tc/drop-columns [:time ])
311- vals
312- (->> (apply tcc/+)))))
313- (tc/select-columns [:violin ])))))
314- (apply tc/concat)))
315-
316- (-> combined-violin-ds
317- :violin
318- normalize
319- audio)
320339
0 commit comments