Chris@12
|
1
|
Chris@12
|
2 program silvet;
|
Chris@12
|
3
|
Chris@12
|
4 { prepareTimeFrequency } = load timefreq;
|
Chris@13
|
5 { loadTemplates, extractRanges } = load templates;
|
Chris@14
|
6
|
Chris@14
|
7 em = load em;
|
Chris@19
|
8 em1 = load em_onecolumn;
|
Chris@12
|
9
|
Chris@12
|
10 mat = load may.matrix;
|
Chris@12
|
11 vec = load may.vector;
|
Chris@14
|
12 plot = load may.plot;
|
Chris@12
|
13
|
Chris@12
|
14 templates = loadTemplates ();
|
Chris@12
|
15
|
Chris@13
|
16 ranges = extractRanges templates;
|
Chris@13
|
17
|
Chris@13
|
18 eprintln "\nWe have \(length (keys templates)) instruments:";
|
Chris@13
|
19 for (sort (keys templates)) do k:
|
Chris@14
|
20 eprintln " * \(k) \(mat.size templates[k]) range \(ranges[k].lowest) -> \(ranges[k].highest)";
|
Chris@13
|
21 done;
|
Chris@12
|
22 eprintln "";
|
Chris@12
|
23
|
Chris@12
|
24 columns = prepareTimeFrequency "test.wav";
|
Chris@12
|
25
|
Chris@19
|
26 /*
|
Chris@12
|
27 height = if empty? columns then 0 else vec.length (head columns) fi;
|
Chris@12
|
28
|
Chris@13
|
29 chunkSize = { rows = height, columns = 100 };
|
Chris@13
|
30
|
Chris@14
|
31 emdata = em.initialise ranges templates 88 chunkSize;
|
Chris@13
|
32
|
Chris@14
|
33 eprintln "initialised EM data: overall pitch range \(emdata.lowest) -> \(emdata.highest)";
|
Chris@13
|
34
|
Chris@12
|
35 chunkify cols =
|
Chris@12
|
36 if empty? cols then []
|
Chris@12
|
37 else
|
Chris@13
|
38 (mat.resizedTo chunkSize
|
Chris@13
|
39 (mat.fromColumns (take chunkSize.columns cols)))
|
Chris@13
|
40 :. \(chunkify (drop chunkSize.columns cols));
|
Chris@12
|
41 fi;
|
Chris@12
|
42
|
Chris@12
|
43 chunks = chunkify columns;
|
Chris@12
|
44
|
Chris@12
|
45 eprintln "we have \(length chunks) chunks of size \(mat.size (head chunks))";
|
Chris@12
|
46
|
Chris@14
|
47 eprintln "attempting one expectation phase...";
|
Chris@13
|
48
|
Chris@14
|
49 error = em.performExpectation emdata (head chunks);
|
Chris@14
|
50
|
Chris@14
|
51 eprintln "done, result has dimension \(mat.size error)";
|
Chris@14
|
52
|
Chris@14
|
53 eprintln "attempting one maximisation phase...";
|
Chris@14
|
54
|
Chris@14
|
55 newP = em.performMaximisation emdata (head chunks) error;
|
Chris@14
|
56
|
Chris@14
|
57 eprintln "done";
|
Chris@14
|
58
|
Chris@14
|
59 \() (plot.plot [ Grid (head chunks) ]);
|
Chris@14
|
60 \() (plot.plot [ Grid error ]);
|
Chris@14
|
61
|
Chris@14
|
62 \() (plot.plot [ Grid newP ]);
|
Chris@19
|
63 */
|
Chris@19
|
64
|
Chris@19
|
65 em1data = em1.initialise ranges templates 88;
|
Chris@19
|
66
|
Chris@21
|
67 col = head (drop ((length columns) / 2) columns);
|
Chris@19
|
68
|
Chris@20
|
69 \() (plot.plot [ Caption "Source frequency distribution", Vector col ]);
|
Chris@19
|
70
|
Chris@20
|
71 sourceGrid d =
|
Chris@20
|
72 mat.fromColumns (map do k: d.sources[k] done (sort (keys d.sources)));
|
Chris@20
|
73
|
Chris@20
|
74 \() (plot.plot [ Caption "Source distribution beforehand", Grid (sourceGrid em1data)]);
|
Chris@20
|
75
|
Chris@20
|
76 oneIteration em1data col n =
|
Chris@19
|
77 (q = em1.performExpectation em1data col;
|
Chris@20
|
78 // \() (plot.plot [ Caption "Frequency distribution and output of E-step for iteration \(n)", Vector col, Vector q ]);
|
Chris@19
|
79 newdata = em1.performMaximisation em1data col q;
|
Chris@21
|
80 if (n % 6 == 0) then
|
Chris@20
|
81 \() (plot.plot [ Caption "Pitch distribution before and after M-step update for iteration \(n)", Vector (em1data.pitches), Vector (newdata.pitches) ]);
|
Chris@20
|
82 \() (plot.plot [ Caption "Source distribution after M-step update for iteration \(n)", Grid (sourceGrid newdata) ]);
|
Chris@20
|
83 fi;
|
Chris@19
|
84 newdata);
|
Chris@19
|
85
|
Chris@21
|
86 iterations = 12;
|
Chris@19
|
87
|
Chris@20
|
88 var d = em1data;
|
Chris@20
|
89
|
Chris@20
|
90 for [1..iterations] do i:
|
Chris@20
|
91 d := oneIteration d col i;
|
Chris@20
|
92 done;
|
Chris@14
|
93
|
Chris@21
|
94 var sounding = [];
|
Chris@21
|
95
|
Chris@21
|
96 println "pitch distribution: \(d.pitches)";
|
Chris@21
|
97
|
Chris@21
|
98 for [d.lowest .. d.highest] do p:
|
Chris@21
|
99 if (vec.at d.pitches p) > 0.05 then
|
Chris@21
|
100 sounding := sounding ++ [p];
|
Chris@21
|
101 fi;
|
Chris@21
|
102 done;
|
Chris@21
|
103
|
Chris@21
|
104 println "Sounding: \(sounding)";
|
Chris@21
|
105
|
Chris@21
|
106 toNote p = (array ["A","A#","B","C","C#","D","D#","E","F","F#","G","G#"])[p % 12];
|
Chris@21
|
107
|
Chris@21
|
108 println "Notes: \(map toNote sounding)";
|
Chris@21
|
109
|
Chris@21
|
110 var instruments = [];
|
Chris@21
|
111 for sounding do p:
|
Chris@21
|
112 var best = "";
|
Chris@21
|
113 var bestp = 0;
|
Chris@21
|
114 for d.instruments do i:
|
Chris@21
|
115 if vec.at d.sources[i] p > bestp then
|
Chris@21
|
116 bestp := vec.at d.sources[i] p;
|
Chris@21
|
117 best := i;
|
Chris@21
|
118 fi;
|
Chris@21
|
119 done;
|
Chris@21
|
120 if bestp > 0 then
|
Chris@21
|
121 instruments := instruments ++ [best];
|
Chris@21
|
122 else
|
Chris@21
|
123 instruments := instruments ++ ["unknown"];
|
Chris@21
|
124 fi;
|
Chris@21
|
125 done;
|
Chris@21
|
126
|
Chris@21
|
127 println "Instruments: \(instruments)";
|
Chris@21
|
128
|
Chris@21
|
129
|
Chris@14
|
130 ();
|
Chris@14
|
131
|
Chris@14
|
132
|