annotate yeti/silvet.yeti @ 21:8e61ec97b34e

I think this may be the right arrangement. Print out some results.
author Chris Cannam
date Fri, 28 Mar 2014 12:52:11 +0000
parents 982aa1197a7e
children 782b910394f3
rev   line source
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