Flake: use new window functions (partial_tukey, punchout_tukey).

This commit is contained in:
Grigory Chudov
2014-08-26 23:48:16 -04:00
parent 5e784881f1
commit 52879ed70c
15 changed files with 1314 additions and 214 deletions

View File

@@ -6,12 +6,14 @@
public int* samples;
public uint done_fixed;
public LpcContext[] lpc_ctx;
public LpcSubframeInfo sf;
public ALACSubframeInfo()
{
best = new ALACSubframe();
lpc_ctx = new LpcContext[Alac.MAX_LPC_WINDOWS];
for (int i = 0; i < Alac.MAX_LPC_WINDOWS; i++)
sf = new LpcSubframeInfo();
lpc_ctx = new LpcContext[lpc.MAX_LPC_WINDOWS];
for (int i = 0; i < lpc.MAX_LPC_WINDOWS; i++)
lpc_ctx[i] = new LpcContext();
}
@@ -21,7 +23,8 @@
best.residual = r;
best.size = AudioSamples.UINT32_MAX;
best.order = 0;
for (int iWindow = 0; iWindow < Alac.MAX_LPC_WINDOWS; iWindow++)
sf.Reset();
for (int iWindow = 0; iWindow < lpc.MAX_LPC_WINDOWS; iWindow++)
lpc_ctx[iWindow].Reset();
done_fixed = 0;
}

View File

@@ -36,7 +36,7 @@ namespace CUETools.Codecs.ALAC
public class ALACWriterSettings: AudioEncoderSettings
{
public ALACWriterSettings()
: base("0 1 2 3 4 5 6 7 8 9 10", "3")
: base("0 1 2 3 4 5 6 7 8 9 10", "5")
{
}
@@ -91,7 +91,8 @@ namespace CUETools.Codecs.ALAC
int[] verifyBuffer;
int[] residualBuffer;
float[] windowBuffer;
int samplesInBuffer = 0;
LpcWindowSection[,] windowSections;
int samplesInBuffer = 0;
int m_blockSize = 0;
int _totalSize = 0;
@@ -122,7 +123,8 @@ namespace CUETools.Codecs.ALAC
samplesBuffer = new int[Alac.MAX_BLOCKSIZE * (Settings.PCM.ChannelCount == 2 ? 5 : Settings.PCM.ChannelCount)];
residualBuffer = new int[Alac.MAX_BLOCKSIZE * (Settings.PCM.ChannelCount == 2 ? 6 : Settings.PCM.ChannelCount + 1)];
windowBuffer = new float[Alac.MAX_BLOCKSIZE * 2 * Alac.MAX_LPC_WINDOWS];
windowBuffer = new float[Alac.MAX_BLOCKSIZE * 2 * lpc.MAX_LPC_WINDOWS];
windowSections = new LpcWindowSection[lpc.MAX_LPC_WINDOWS, lpc.MAX_LPC_SECTIONS];
eparams.set_defaults(m_settings.EncoderModeIndex);
@@ -875,7 +877,11 @@ namespace CUETools.Codecs.ALAC
LpcContext lpc_ctx = frame.subframes[ch].lpc_ctx[iWindow];
lpc_ctx.GetReflection(max_order, smp, n, frame.window_buffer + iWindow * Alac.MAX_BLOCKSIZE * 2);
fixed (LpcWindowSection* sections = &windowSections[iWindow, 0])
lpc_ctx.GetReflection(
frame.subframes[ch].sf, max_order, smp,
frame.window_buffer + iWindow * Alac.MAX_BLOCKSIZE * 2, sections,
Settings.PCM.BitsPerSample * 2 + BitReader.log2i(frame.blocksize) >= 61);
lpc_ctx.ComputeLPC(lpcs);
lpc_ctx.SortOrdersAkaike(frame.blocksize, eparams.estimation_depth, min_order, max_order, 5.0, 1.0/18);
for (i = 0; i < eparams.estimation_depth && i < max_order; i++)
@@ -945,20 +951,36 @@ namespace CUETools.Codecs.ALAC
{
case WindowMethod.Estimate:
{
int best_window = -1;
double best_error = 0;
int order = 2;
for (int i = 0; i < _windowcount; i++)
int order = 4;
float* err = stackalloc float[lpc.MAX_LPC_ORDER];
for (int i = 0; i < _windowcount; i++)
{
frame.subframes[ch].lpc_ctx[i].GetReflection(order, frame.subframes[ch].samples, frame.blocksize, frame.window_buffer + i * Alac.MAX_BLOCKSIZE * 2);
double err = frame.subframes[ch].lpc_ctx[i].prediction_error[order - 1] / frame.subframes[ch].lpc_ctx[i].autocorr_values[0];
if (best_window == -1 || best_error > err)
{
best_window = i;
best_error = err;
}
LpcContext lpc_ctx = frame.subframes[ch].lpc_ctx[i];
fixed (LpcWindowSection* sections = &windowSections[i, 0])
lpc_ctx.GetReflection(
frame.subframes[ch].sf, order,
frame.subframes[ch].samples,
frame.window_buffer + i * Alac.MAX_BLOCKSIZE * 2, sections,
Settings.PCM.BitsPerSample * 2 + BitReader.log2i(frame.blocksize) >= 61);
lpc_ctx.SortOrdersAkaike(frame.blocksize, 1, 1, order, 4.5, 0.0);
err[i] = (float)(lpc_ctx.Akaike(frame.blocksize, lpc_ctx.best_orders[0], 4.5, 0.0) - frame.blocksize * Math.Log(lpc_ctx.autocorr_values[0]) / 2);
}
return best_window;
int* best_windows = stackalloc int[lpc.MAX_LPC_ORDER];
for (int i = 0; i < _windowcount; i++)
best_windows[i] = i;
for (int i = 0; i < _windowcount; i++)
{
for (int j = i + 1; j < _windowcount; j++)
{
if (err[best_windows[i]] > err[best_windows[j]])
{
int tmp = best_windows[j];
best_windows[j] = best_windows[i];
best_windows[i] = tmp;
}
}
}
return best_windows[0];
}
case WindowMethod.Evaluate:
encode_residual_pass1(frame, ch, -1);
@@ -978,10 +1000,16 @@ namespace CUETools.Codecs.ALAC
case StereoMethod.Estimate:
for (int ch = 0; ch < subframes; ch++)
{
LpcContext lpc_ctx = frame.subframes[ch].lpc_ctx[0];
int stereo_order = Math.Min(8, eparams.max_prediction_order);
int iWindow = 0;
LpcContext lpc_ctx = frame.subframes[ch].lpc_ctx[iWindow];
int stereo_order = Math.Min(8, eparams.max_prediction_order);
double alpha = 1.5; // 4.5 + eparams.max_prediction_order / 10.0;
lpc_ctx.GetReflection(stereo_order, frame.subframes[ch].samples, frame.blocksize, frame.window_buffer);
fixed (LpcWindowSection* sections = &windowSections[iWindow, 0])
lpc_ctx.GetReflection(
frame.subframes[ch].sf, stereo_order,
frame.subframes[ch].samples,
frame.window_buffer + iWindow * Alac.MAX_BLOCKSIZE * 2, sections,
Settings.PCM.BitsPerSample * 2 + BitReader.log2i(frame.blocksize) >= 61);
lpc_ctx.SortOrdersAkaike(frame.blocksize, 1, 1, stereo_order, alpha, 0);
frame.subframes[ch].best.size = (uint)Math.Max(0, lpc_ctx.Akaike(frame.blocksize, lpc_ctx.best_orders[0], alpha, 0));
}
@@ -1060,13 +1088,17 @@ namespace CUETools.Codecs.ALAC
unsafe void calculate_window(float * window, window_function func, WindowFunction flag)
{
if ((eparams.window_function & flag) == 0 || _windowcount == Alac.MAX_LPC_WINDOWS)
if ((eparams.window_function & flag) == 0 || _windowcount == lpc.MAX_LPC_WINDOWS)
return;
int sz = _windowsize;
float* pos = window + _windowcount * Alac.MAX_BLOCKSIZE * 2;
do
{
func(pos, sz);
windowSections[_windowcount, 0].setData(0, sz);
for (int j = 1; j < lpc.MAX_LPC_SECTIONS; j++)
windowSections[_windowcount, j].setZero(sz, sz);
func(pos, sz);
break;
if ((sz & 1) != 0)
break;
pos += sz;
@@ -1091,9 +1123,33 @@ namespace CUETools.Codecs.ALAC
calculate_window(window, lpc.window_tukey, WindowFunction.Tukey);
calculate_window(window, lpc.window_hann, WindowFunction.Hann);
calculate_window(window, lpc.window_flattop, WindowFunction.Flattop);
if (_windowcount == 0)
int tukey_parts = 2;
double overlap = -0.3;
double overlap_units = overlap / (1.0 - overlap);
for (int m = 0; m < tukey_parts; m++)
calculate_window(window, (w, wsz) =>
{
lpc.window_punchout_tukey(w, wsz, 0.1,
m / (tukey_parts + overlap_units),
(m + 1 + overlap_units) / (tukey_parts + overlap_units));
}, WindowFunction.PartialTukey);
tukey_parts = 3;
overlap = -0.1;
//overlap = 0.1;
overlap_units = overlap / (1.0 - overlap);
for (int m = 0; m < tukey_parts; m++)
calculate_window(window, (w, wsz) =>
{
lpc.window_punchout_tukey(w, wsz, 0.1,
m / (tukey_parts + overlap_units),
(m + 1 + overlap_units) / (tukey_parts + overlap_units));
}, WindowFunction.PunchoutTukey);
if (_windowcount == 0)
throw new Exception("invalid windowfunction");
}
fixed (LpcWindowSection* sections = &windowSections[0, 0])
LpcWindowSection.Detect(_windowcount, window, Alac.MAX_BLOCKSIZE * 2, _windowsize, sections);
}
frame.window_buffer = window;
int bps = Settings.PCM.BitsPerSample + Settings.PCM.ChannelCount - 1;
@@ -1793,11 +1849,13 @@ namespace CUETools.Codecs.ALAC
case 4:
stereo_method = StereoMethod.Estimate;
window_method = WindowMethod.Estimate;
window_function = WindowFunction.PartialTukey | WindowFunction.PunchoutTukey;
max_prediction_order = 8;
break;
case 5:
stereo_method = StereoMethod.Estimate;
window_method = WindowMethod.Estimate;
window_function = WindowFunction.PartialTukey | WindowFunction.PunchoutTukey;
break;
case 6:
stereo_method = StereoMethod.Estimate;

View File

@@ -28,7 +28,6 @@ namespace CUETools.Codecs.ALAC
public const int MAX_RICE_PARAM = 14;
public const int MAX_PARTITION_ORDER = 8;
public const int MAX_PARTITIONS = 1 << MAX_PARTITION_ORDER;
public const int MAX_LPC_WINDOWS = 4;
public const uint UINT32_MAX = 0xffffffff;

View File

@@ -7,6 +7,8 @@
Hann = 4,
Flattop = 8,
Bartlett = 16,
TukFlat = 10
TukFlat = 10,
PartialTukey = 32,
PunchoutTukey = 64,
}
}

View File

@@ -270,6 +270,8 @@ namespace CUETools.Codecs.FLACCL
string _path;
long _position;
public const int MAX_LPC_WINDOWS = 2;
// number of audio channels
// valid values are 1 to 8
int channels, ch_code;
@@ -1090,7 +1092,7 @@ namespace CUETools.Codecs.FLACCL
unsafe void calculate_window(FLACCLTask task, window_function func, WindowFunction flag)
{
if ((eparams.window_function & flag) == 0 || task.nWindowFunctions == lpc.MAX_LPC_WINDOWS)
if ((eparams.window_function & flag) == 0 || task.nWindowFunctions == MAX_LPC_WINDOWS)
return;
func(((float*)task.clWindowFunctionsPtr) + task.nWindowFunctions * task.frameSize, task.frameSize);
@@ -2508,10 +2510,10 @@ namespace CUETools.Codecs.FLACCL
int residualBufferLen = sizeof(int) * MAX_CHANNELSIZE * channels; // need to adjust residualOffset?
int partitionsLen = sizeof(int) * ((writer.Settings.PCM.BitsPerSample > 16 ? 31 : 15) * 2 << 8) * channels * MAX_FRAMES;
int riceParamsLen = sizeof(int) * (4 << 8) * channels * MAX_FRAMES;
int autocorLen = sizeof(float) * (MAX_ORDER + 1) * lpc.MAX_LPC_WINDOWS * channelsCount * MAX_FRAMES;
int autocorLen = sizeof(float) * (MAX_ORDER + 1) * FLACCLWriter.MAX_LPC_WINDOWS * channelsCount * MAX_FRAMES;
int lpcDataLen = autocorLen * 32;
int resOutLen = sizeof(int) * channelsCount * (lpc.MAX_LPC_WINDOWS * lpc.MAX_LPC_ORDER + 8) * MAX_FRAMES;
int wndLen = sizeof(float) * MAX_CHANNELSIZE /** 2*/ * lpc.MAX_LPC_WINDOWS;
int resOutLen = sizeof(int) * channelsCount * (FLACCLWriter.MAX_LPC_WINDOWS * lpc.MAX_LPC_ORDER + 8) * MAX_FRAMES;
int wndLen = sizeof(float) * MAX_CHANNELSIZE /** 2*/ * FLACCLWriter.MAX_LPC_WINDOWS;
int selectedLen = sizeof(int) * 32 * channelsCount * MAX_FRAMES;
int riceLen = sizeof(int) * channels * MAX_CHANNELSIZE;

View File

@@ -11,6 +11,7 @@
public int frame_number;
public FlacSubframe current;
public float* window_buffer;
public int nSeg = 0;
public BitWriter writer = null;
public int writer_offset = 0;

View File

@@ -7,6 +7,7 @@ namespace CUETools.Codecs.FLAKE
public FlacSubframeInfo()
{
best = new FlacSubframe();
sf = new LpcSubframeInfo();
lpc_ctx = new LpcContext[lpc.MAX_LPC_WINDOWS];
for (int i = 0; i < lpc.MAX_LPC_WINDOWS; i++)
lpc_ctx[i] = new LpcContext();
@@ -22,8 +23,10 @@ namespace CUETools.Codecs.FLAKE
best.residual = r;
best.type = SubframeType.Verbatim;
best.size = AudioSamples.UINT32_MAX;
sf.Reset();
for (int iWindow = 0; iWindow < lpc.MAX_LPC_WINDOWS; iWindow++)
lpc_ctx[iWindow].Reset();
//sf.obits = obits;
done_fixed = 0;
}
@@ -33,5 +36,6 @@ namespace CUETools.Codecs.FLAKE
public int* samples;
public uint done_fixed;
public LpcContext[] lpc_ctx;
public LpcSubframeInfo sf;
};
}

File diff suppressed because it is too large Load Diff

View File

@@ -6,6 +6,6 @@ namespace CUETools.Codecs.FLAKE
Independent = 0,
Estimate = 1,
Evaluate = 2,
Search = 3
Search = 3,
}
}

View File

@@ -2,11 +2,17 @@
{
public enum WindowFunction
{
None = 0,
Welch = 1,
Tukey = 2,
Hann = 4,
Flattop = 8,
Bartlett = 16,
TukeyFlattop = 10
TukeyFlattop = 10,
PartialTukey = 32,
TukeyPartialTukey = 34,
TukeyFlattopPartialTukey = 42,
PunchoutTukey = 64,
TukeyFlattopPartialTukeyPunchoutTukey = 106,
}
}

View File

@@ -2,8 +2,16 @@
{
public enum WindowMethod
{
Estimate = 0,
Evaluate = 1,
Search = 2
Evaluate,
Search,
Estimate,
Estimate2,
Estimate3,
EstimateN,
Evaluate2,
Evaluate2N,
Evaluate3,
Evaluate3N,
EvaluateN,
}
}

View File

@@ -73,6 +73,19 @@ namespace CUETools.Codecs
}
}
public bool HasDefaultValuesForMode(int index)
{
bool res = true;
foreach (PropertyDescriptor property in TypeDescriptor.GetProperties(this))
foreach (var attribute in property.Attributes)
if (attribute is DefaultValueForModeAttribute)
{
var defaultValueForMode = attribute as DefaultValueForModeAttribute;
res &= (int)property.GetValue(this) == defaultValueForMode.m_values[index];
}
return res;
}
[Browsable(false)]
[XmlIgnore]
public AudioPCMConfig PCM

View File

@@ -5,8 +5,9 @@ namespace CUETools.Codecs
public class lpc
{
public const int MAX_LPC_ORDER = 32;
public const int MAX_LPC_WINDOWS = 2;
public const int MAX_LPC_WINDOWS = 8;
public const int MAX_LPC_PRECISIONS = 4;
public const int MAX_LPC_SECTIONS = 32;
public unsafe static void window_welch(float* window, int L)
{
@@ -21,7 +22,7 @@ namespace CUETools.Codecs
}
}
public unsafe static void window_bartlett(float* window, int L)
public unsafe static void window_bartlett(float* window, int L)
{
int N = L - 1;
double N2 = (double)N / 2.0;
@@ -33,35 +34,65 @@ namespace CUETools.Codecs
}
}
public unsafe static void window_rectangle(float* window, int L)
public unsafe static void window_rectangle(float* window, int L)
{
for (int n = 0; n < L; n++)
window[n] = 1.0F;
}
public unsafe static void window_flattop(float* window, int L)
public unsafe static void window_flattop(float* window, int L)
{
int N = L - 1;
for (int n = 0; n < L; n++)
window[n] = (float)(1.0 - 1.93 * Math.Cos(2.0 * Math.PI * n / N) + 1.29 * Math.Cos(4.0 * Math.PI * n / N) - 0.388 * Math.Cos(6.0 * Math.PI * n / N) + 0.0322 * Math.Cos(8.0 * Math.PI * n / N));
}
public unsafe static void window_tukey(float* window, int L)
public unsafe static void window_tukey(float* window, int L)
{
window_rectangle(window, L);
double p = 0.5;
int Np = (int)(p / 2.0 * L) - 1;
int z = 0;
int Np = (int)(p / 2.0 * L) - z;
if (Np > 0)
{
for (int n = 0; n <= Np; n++)
{
window[n] = (float)(0.5 - 0.5 * Math.Cos(Math.PI * n / Np));
window[L - Np - 1 + n] = (float)(0.5 - 0.5 * Math.Cos(Math.PI * (n + Np) / Np));
}
for (int n = 0; n < z; n++)
window[n] = window[L - n - 1] = 0;
for (int n = 0; n < Np - 1; n++)
window[n + z] = window[L - n - 1 - z] = (float)(0.5 - 0.5 * Math.Cos(Math.PI * (n + 1) / Np));
for (int n = z + Np - 1; n < L - z - Np + 1; n++)
window[n] = 1.0F;
}
}
public unsafe static void window_hann(float* window, int L)
public unsafe static void window_punchout_tukey(float* window, int L, double p, double start, double end)
{
int start_n = (int)(start * L);
int end_n = (int)(end * L);
int Np = (int)(p / 2.0 * L);
int i, n = 0;
if (start_n != 0)
{
for (i = 1; n < Np; n++, i++)
window[n] = (float)(0.5 - 0.5 * Math.Cos(Math.PI * i / Np));
for (; n < start_n - Np; n++)
window[n] = 1.0f;
for (i = Np; n < start_n; n++, i--)
window[n] = (float)(0.5 - 0.5 * Math.Cos(Math.PI * i / Np));
}
for (; n < end_n; n++)
window[n] = 0.0f;
if (end_n != L)
{
for (i = 1; n < end_n + Np; n++, i++)
window[n] = (float)(0.5 - 0.5 * Math.Cos(Math.PI * i / Np));
for (; n < L - Np; n++)
window[n] = 1.0f;
for (i = Np; n < L; n++, i--)
window[n] = (float)(0.5 - 0.5 * Math.Cos(Math.PI * i / Np));
}
}
public unsafe static void window_hann(float* window, int L)
{
int N = L - 1;
for (int n = 0; n < L; n++)
@@ -73,6 +104,7 @@ namespace CUETools.Codecs
return (short)((val >> 31) + ((val - 1) >> 31) + 1);
}
#if XXX
static public unsafe void
compute_corr_int(/*const*/ short* data1, short* data2, int len, int min, int lag, int* autoc)
{
@@ -92,13 +124,14 @@ namespace CUETools.Codecs
autoc[i] = temp + temp2;
}
}
#endif
/**
* Calculates autocorrelation data from audio samples
* A window function is applied before calculation.
*/
static public unsafe void
compute_autocorr(/*const*/ int* data, int len, int min, int lag, double* autoc, float* window)
compute_autocorr(/*const*/ int* data, float* window, int len, int min, int lag, double* autoc, int prev, int next)
{
#if FPAC
short* data1 = stackalloc short[len + 1];
@@ -125,30 +158,224 @@ namespace CUETools.Codecs
for (int coeff = min; coeff <= lag; coeff++)
autoc[coeff] = (c1[coeff] * (double)(1 << 18) + (c2[coeff] + c3[coeff]) * (double)(1 << 9) + c4[coeff]);
#else
double* data1 = stackalloc double[(int)len + 16];
int i;
#if XXX
if (min == 0 && lag >= 4)
{
int* pdata = data;
float* pwindow = window;
for (i = 0; i < len; i++)
data1[i] = data[i] * window[i];
data1[len] = 0;
double temp0 = 1.0;
double temp1 = 1.0;
double temp2 = 1.0;
double temp3 = 1.0;
double temp4 = 1.0;
for (i = min; i <= lag; ++i)
double c0 = *(pdata++) * *(pwindow++);
float c1 = *(pdata++) * *(pwindow++);
float c2 = *(pdata++) * *(pwindow++);
float c3 = *(pdata++) * *(pwindow++);
float c4 = *(pdata++) * *(pwindow++);
int* finish = data + len;
while (pdata <= finish)
{
temp0 += c0 * c0;
temp1 += c0 * c1;
temp2 += c0 * c2;
temp3 += c0 * c3;
temp4 += c0 * c4;
c0 = c1;
c1 = c2;
c2 = c3;
c3 = c4;
c4 = *(pdata++) * *(pwindow++);
}
temp0 += c0 * c0;
temp1 += c0 * c1;
temp2 += c0 * c2;
temp3 += c0 * c3;
c0 = c1;
c1 = c2;
c2 = c3;
temp0 += c0 * c0;
temp1 += c0 * c1;
temp2 += c0 * c2;
c0 = c1;
c1 = c2;
temp0 += c0 * c0;
temp1 += c0 * c1;
c0 = c1;
temp0 += c0 * c0;
autoc[0] += temp0;
autoc[1] += temp1;
autoc[2] += temp2;
autoc[3] += temp3;
autoc[4] += temp4;
min = 5;
if (lag < min) return;
}
#endif
double* data1 = stackalloc double[lag + len + lag];
int i;
for (i = 0; i < lag; i++)
data1[i] = prev != 0 ? data[i - lag] : 0;
for (i = 0; i < len; i++)
data1[lag + i] = data[i] * window[i];
for (i = 0; i < lag; i++)
data1[lag + len + i] = next != 0 ? data[len + i] : 0;
for (i = min; i <= lag; ++i)
{
double temp = 1.0;
double temp2 = 1.0;
double* finish = data1 + len - i;
double temp = 0;
double temp2 = 0;
double* pdata = data1 + lag - i;
double* finish = data1 + lag + len - 1;
for (double* pdata = data1; pdata < finish; pdata += 2)
{
temp += pdata[i] * pdata[0];
temp2 += pdata[i + 1] * pdata[1];
}
autoc[i] = temp + temp2;
while (pdata < finish)
{
temp += pdata[i] * (*pdata++);
temp2 += pdata[i] * (*pdata++);
}
if (pdata <= finish)
temp += pdata[i] * (*pdata++);
autoc[i] += temp + temp2;
}
#endif
}
/**
static public unsafe void
compute_autocorr_windowless(/*const*/ int* data, int len, int min, int lag, double* autoc)
{
// if databits*2 + log2(len) <= 64
#if !XXX
#if XXX
if (min == 0 && lag >= 4)
{
long temp0 = 0;
long temp1 = 0;
long temp2 = 0;
long temp3 = 0;
long temp4 = 0;
int* pdata = data;
int* finish = data + len - 4;
while (pdata < finish)
{
long c0 = *(pdata++);
temp0 += c0 * c0;
temp1 += c0 * pdata[0];
temp2 += c0 * pdata[1];
temp3 += c0 * pdata[2];
temp4 += c0 * pdata[3];
}
{
long c0 = *(pdata++);
temp0 += c0 * c0;
temp1 += c0 * pdata[0];
temp2 += c0 * pdata[1];
temp3 += c0 * pdata[2];
}
{
long c0 = *(pdata++);
temp0 += c0 * c0;
temp1 += c0 * pdata[0];
temp2 += c0 * pdata[1];
}
{
long c0 = *(pdata++);
temp0 += c0 * c0;
temp1 += c0 * pdata[0];
}
{
long c0 = *(pdata++);
temp0 += c0 * c0;
}
autoc[0] += temp0;
autoc[1] += temp1;
autoc[2] += temp2;
autoc[3] += temp3;
autoc[4] += temp4;
min = 5;
if (lag < min) return;
}
#endif
for (int i = min; i <= lag; ++i)
{
long temp = 0;
long temp2 = 0;
int* pdata = data;
int* finish = data + len - i - 1;
while (pdata < finish)
{
temp += (long)pdata[i] * (*pdata++);
temp2 += (long)pdata[i] * (*pdata++);
}
if (pdata <= finish)
temp += (long)pdata[i] * (*pdata++);
autoc[i] += temp + temp2;
}
#else
for (int i = min; i <= lag; ++i)
{
double temp = 0;
double temp2 = 0;
int* pdata = data;
int* finish = data + len - i - 1;
while (pdata < finish)
{
temp += (double)pdata[i] * (double)(*pdata++);
temp2 += (double)pdata[i] * (double)(*pdata++);
}
if (pdata <= finish)
temp += (double)pdata[i] * (double)(*pdata++);
autoc[i] += temp + temp2;
}
#endif
}
static public unsafe void
compute_autocorr_windowless_large(/*const*/ int* data, int len, int min, int lag, double* autoc)
{
for (int i = min; i <= lag; ++i)
{
double temp = 0;
double temp2 = 0;
int* pdata = data;
int* finish = data + len - i - 1;
while (pdata < finish)
{
temp += (long)pdata[i] * (*pdata++);
temp2 += (long)pdata[i] * (*pdata++);
}
if (pdata <= finish)
temp += (long)pdata[i] * (*pdata++);
autoc[i] += temp + temp2;
}
}
static public unsafe void
compute_autocorr_glue(/*const*/ int* data, int min, int lag, double* autoc)
{
for (int i = min; i <= lag; ++i)
{
long temp = 0;
int* pdata = data - i;
int* finish = data;
while (pdata < finish)
temp += (long)pdata[i] * (*pdata++);
autoc[i] += temp;
}
}
/**
* Levinson-Durbin recursion.
* Produces LPC coefficients from autocorrelation data.
*/

View File

@@ -1,7 +1,173 @@
using System;
using System.Collections.Generic;
namespace CUETools.Codecs
{
unsafe public class LpcSubframeInfo
{
public LpcSubframeInfo()
{
autocorr_section_values = new double[lpc.MAX_LPC_SECTIONS, lpc.MAX_LPC_ORDER + 1];
autocorr_section_orders = new int[lpc.MAX_LPC_SECTIONS];
}
// public LpcContext[] lpc_ctx;
public double[,] autocorr_section_values;
public int[] autocorr_section_orders;
//public int obits;
public void Reset()
{
for (int sec = 0; sec < autocorr_section_orders.Length; sec++)
autocorr_section_orders[sec] = 0;
}
}
unsafe public struct LpcWindowSection
{
public enum SectionType
{
Zero,
One,
Data,
Glue
};
public int m_start;
public int m_end;
public SectionType m_type;
public int m_id;
public LpcWindowSection(int end)
{
m_id = -1;
m_start = 0;
m_end = end;
m_type = SectionType.Data;
}
public void setData(int start, int end)
{
m_id = -1;
m_start = start;
m_end = end;
m_type = SectionType.Data;
}
public void setOne(int start, int end)
{
m_id = -1;
m_start = start;
m_end = end;
m_type = SectionType.One;
}
public void setGlue(int start)
{
m_id = -1;
m_start = start;
m_end = start;
m_type = SectionType.Glue;
}
public void setZero(int start, int end)
{
m_id = -1;
m_start = start;
m_end = end;
m_type = SectionType.Zero;
}
unsafe public static void Detect(int _windowcount, float* window_segment, int stride, int sz, LpcWindowSection* sections)
{
int section_id = 0;
var boundaries = new List<int>();
var types = new LpcWindowSection.SectionType[_windowcount, lpc.MAX_LPC_SECTIONS * 2];
for (int x = 0; x < sz; x++)
{
for (int i = 0; i < _windowcount; i++)
{
float w = window_segment[i * stride + x];
types[i, boundaries.Count] =
boundaries.Count >= lpc.MAX_LPC_SECTIONS * 2 - 2 ?
LpcWindowSection.SectionType.Data : w == 0.0 ?
LpcWindowSection.SectionType.Zero : w == 1.0 ?
LpcWindowSection.SectionType.One :
LpcWindowSection.SectionType.Data;
}
bool isBoundary = false;
for (int i = 0; i < _windowcount; i++)
{
isBoundary |= boundaries.Count == 0 ||
types[i, boundaries.Count - 1] != types[i, boundaries.Count];
}
if (isBoundary) boundaries.Add(x);
}
boundaries.Add(sz);
var ones = new int[boundaries.Count - 1];
// Reconstruct segments list.
for (int i = 0; i < _windowcount; i++)
{
int secs = 0;
for (int j = 0; j < boundaries.Count - 1; j++)
{
if (types[i, j] == LpcWindowSection.SectionType.Zero)
{
if (secs > 0 && sections[i * lpc.MAX_LPC_SECTIONS + secs - 1].m_end == boundaries[j] && sections[i * lpc.MAX_LPC_SECTIONS + secs - 1].m_type == LpcWindowSection.SectionType.Zero)
{
sections[i * lpc.MAX_LPC_SECTIONS + secs - 1].m_end = boundaries[j + 1];
continue;
}
sections[i * lpc.MAX_LPC_SECTIONS + secs++].setZero(boundaries[j], boundaries[j + 1]);
continue;
}
if (types[i, j] == LpcWindowSection.SectionType.Data
|| secs + 1 >= lpc.MAX_LPC_SECTIONS
|| (boundaries[j + 1] - boundaries[j] < lpc.MAX_LPC_ORDER))
{
if (secs > 0 && sections[i * lpc.MAX_LPC_SECTIONS + secs - 1].m_end == boundaries[j] && sections[i * lpc.MAX_LPC_SECTIONS + secs - 1].m_type == LpcWindowSection.SectionType.Data)
{
sections[i * lpc.MAX_LPC_SECTIONS + secs - 1].m_end = boundaries[j + 1];
continue;
}
sections[i * lpc.MAX_LPC_SECTIONS + secs++].setData(boundaries[j], boundaries[j + 1]);
continue;
}
if (secs > 0 && sections[i * lpc.MAX_LPC_SECTIONS + secs - 1].m_end == boundaries[j] && sections[i * lpc.MAX_LPC_SECTIONS + secs - 1].m_type == LpcWindowSection.SectionType.One)
sections[i * lpc.MAX_LPC_SECTIONS + secs++].setGlue(boundaries[j]);
sections[i * lpc.MAX_LPC_SECTIONS + secs++].setOne(boundaries[j], boundaries[j + 1]);
ones[j] |= 1 << i;
}
while (secs < lpc.MAX_LPC_SECTIONS)
sections[i * lpc.MAX_LPC_SECTIONS + secs++].setZero(sz, sz);
}
for (int j = 0; j < boundaries.Count - 1; j++)
{
if (j > 0 && ones[j - 1] == ones[j])
{
for (int i = 0; i < _windowcount; i++)
{
for (int sec = 0; sec < lpc.MAX_LPC_SECTIONS; sec++)
if (sections[i * lpc.MAX_LPC_SECTIONS + sec].m_type == LpcWindowSection.SectionType.Glue &&
sections[i * lpc.MAX_LPC_SECTIONS + sec].m_start == boundaries[j])
{
sections[i * lpc.MAX_LPC_SECTIONS + sec - 1].m_end = sections[i * lpc.MAX_LPC_SECTIONS + sec + 1].m_end;
for (int sec1 = sec; sec1 + 2 < lpc.MAX_LPC_SECTIONS; sec1++)
sections[i * lpc.MAX_LPC_SECTIONS + sec1] = sections[i * lpc.MAX_LPC_SECTIONS + sec1 + 2];
}
}
continue;
}
if ((ones[j] & (ones[j] - 1)) != 0 && section_id < lpc.MAX_LPC_SECTIONS)
{
for (int i = 0; i < _windowcount; i++)
for (int sec = 0; sec < lpc.MAX_LPC_SECTIONS; sec++)
if (sections[i * lpc.MAX_LPC_SECTIONS + sec].m_type == LpcWindowSection.SectionType.One &&
sections[i * lpc.MAX_LPC_SECTIONS + sec].m_start == boundaries[j])
{
sections[i * lpc.MAX_LPC_SECTIONS + sec].m_id = section_id;
}
section_id++;
}
}
}
}
/// <summary>
/// Context for LPC coefficients calculation and order estimation
/// </summary>
@@ -36,18 +202,62 @@ namespace CUETools.Codecs
/// <param name="samples">Samples pointer</param>
/// <param name="blocksize">Block size</param>
/// <param name="window">Window function</param>
public void GetReflection(int order, int* samples, int blocksize, float* window)
public void GetReflection(LpcSubframeInfo subframe, int order, int* samples, float* window, LpcWindowSection* sections, bool large)
{
if (autocorr_order > order)
return;
fixed (double* reff = reflection_coeffs, autoc = autocorr_values, err = prediction_error)
{
lpc.compute_autocorr(samples, blocksize, autocorr_order, order, autoc, window);
for (int i = autocorr_order; i <= order; i++) autoc[i] = 0;
int prev = 0;
for (int section = 0; section < lpc.MAX_LPC_SECTIONS; section++)
{
if (sections[section].m_type == LpcWindowSection.SectionType.Zero)
{
prev = 0;
continue;
}
if (sections[section].m_type == LpcWindowSection.SectionType.Data)
{
int next = section + 1 < lpc.MAX_LPC_SECTIONS && sections[section + 1].m_type == LpcWindowSection.SectionType.One ? 1 : 0;
lpc.compute_autocorr(samples + sections[section].m_start, window + sections[section].m_start, sections[section].m_end - sections[section].m_start, autocorr_order, order, autoc, prev, next);
}
else if (sections[section].m_type == LpcWindowSection.SectionType.Glue)
lpc.compute_autocorr_glue(samples + sections[section].m_start, autocorr_order, order, autoc);
else if (sections[section].m_type == LpcWindowSection.SectionType.One)
{
if (sections[section].m_id >= 0)
{
if (subframe.autocorr_section_orders[sections[section].m_id] <= order)
{
fixed (double* autocsec = &subframe.autocorr_section_values[sections[section].m_id, 0])
{
for (int i = subframe.autocorr_section_orders[sections[section].m_id]; i <= order; i++) autocsec[i] = 0;
if (large)
lpc.compute_autocorr_windowless_large(samples + sections[section].m_start, sections[section].m_end - sections[section].m_start, subframe.autocorr_section_orders[sections[section].m_id], order, autocsec);
else
lpc.compute_autocorr_windowless(samples + sections[section].m_start, sections[section].m_end - sections[section].m_start, subframe.autocorr_section_orders[sections[section].m_id], order, autocsec);
}
subframe.autocorr_section_orders[sections[section].m_id] = order + 1;
}
for (int i = autocorr_order; i <= order; i++)
autoc[i] += subframe.autocorr_section_values[sections[section].m_id, i];
}
else
{
if (large)
lpc.compute_autocorr_windowless_large(samples + sections[section].m_start, sections[section].m_end - sections[section].m_start, autocorr_order, order, autoc);
else
lpc.compute_autocorr_windowless(samples + sections[section].m_start, sections[section].m_end - sections[section].m_start, autocorr_order, order, autoc);
}
prev = 1;
}
}
lpc.compute_schur_reflection(autoc, (uint)order, reff, err);
autocorr_order = order + 1;
}
}
#if XXX
public void GetReflection1(int order, int* samples, int blocksize, float* window)
{
if (autocorr_order > order)
@@ -83,11 +293,12 @@ namespace CUETools.Codecs
autocorr_order = order + 1;
}
}
#endif
public double Akaike(int blocksize, int order, double alpha, double beta)
{
//return (blocksize - order) * (Math.Log(prediction_error[order - 1]) - Math.Log(1.0)) + Math.Log(blocksize) * order * (alpha + beta * order);
return blocksize * Math.Log(prediction_error[order - 1]) + Math.Log(blocksize) * order * (alpha + beta * order);
//return blocksize * (Math.Log(prediction_error[order - 1]) - Math.Log(autocorr_values[0]) / 2) + Math.Log(blocksize) * order * (alpha + beta * order);
return blocksize * (Math.Log(prediction_error[order - 1])) + Math.Log(blocksize) * order * (alpha + beta * order);
}
/// <summary>

View File

@@ -475,25 +475,45 @@ namespace CUETools.FlakeExe
if (debug)
{
settings = flake.Settings as FlakeWriterSettings;
Console.SetOut(stdout);
Console.Out.WriteLine("{0}\t{1:0.00}\t{2}\t{3}\t{4}\t{5}\t{6}..{7}\t{8}..{9}\t{10}..{11}\t{12}..{13}\t{14}\t{15}\t{16}",
Console.Out.WriteLine("{17}\t{0}\t{1:0.00}\t{2}\t{3}\t{4}\t{5}\t{6}..{7}\t{8}..{9}\t{10}..{11}\t{12}..{13}\t{14}\t{15}\t{16}\t{18}",
flake.TotalSize,
flake.UserProcessorTime.TotalSeconds > 0 ? flake.UserProcessorTime.TotalSeconds : totalElapsed.TotalSeconds,
flake.PredictionType.ToString().PadRight(15),
flake.StereoMethod.ToString().PadRight(15),
(flake.OrderMethod.ToString() + "(" + flake.EstimationDepth.ToString() + ")").PadRight(15),
flake.WindowFunction,
(flake.Settings as FlakeWriterSettings).MinPartitionOrder,
(flake.Settings as FlakeWriterSettings).MaxPartitionOrder,
(flake.Settings as FlakeWriterSettings).MinLPCOrder,
(flake.Settings as FlakeWriterSettings).MaxLPCOrder,
(flake.Settings as FlakeWriterSettings).MinFixedOrder,
(flake.Settings as FlakeWriterSettings).MaxFixedOrder,
(flake.WindowMethod.ToString().PadRight(10) + "(" +
((flake.WindowFunction & WindowFunction.Tukey) != 0 ? "T" : " ") +
((flake.WindowFunction & WindowFunction.PartialTukey) != 0 ? "P" : " ") +
(((flake.WindowFunction & WindowFunction.PunchoutTukey) != 0 ? "O" : "") +
((flake.WindowFunction & WindowFunction.Welch) == 0 ? "" : "W") +
((flake.WindowFunction & WindowFunction.Hann) == 0 ? "" : "H") +
((flake.WindowFunction & WindowFunction.Flattop) == 0 ? "" : "F") +
((flake.WindowFunction & WindowFunction.Bartlett) == 0 ? "" : "B")).PadRight(1))
+")",
settings.MinPartitionOrder,
settings.MaxPartitionOrder,
settings.MinLPCOrder,
settings.MaxLPCOrder,
settings.MinFixedOrder,
settings.MaxFixedOrder,
flake.MinPrecisionSearch,
flake.MaxPrecisionSearch,
flake.Settings.BlockSize,
flake.VBRMode,
coeffs ?? ""
coeffs ?? "",
audioSource.Position * audioSource.PCM.BlockAlign,
(flake.PredictionType == PredictionType.Search && flake.StereoMethod == StereoMethod.Evaluate && flake.WindowMethod == WindowMethod.EvaluateN && flake.WindowFunction == (WindowFunction.PunchoutTukey | WindowFunction.PartialTukey | WindowFunction.Tukey) && flake.EstimationDepth == 2 && flake.MinPrecisionSearch == 0 && settings.HasDefaultValuesForMode(8)) ? "8" :
(flake.PredictionType == PredictionType.Levinson && flake.StereoMethod == StereoMethod.Evaluate && flake.WindowMethod == WindowMethod.EvaluateN && flake.WindowFunction == (WindowFunction.PunchoutTukey | WindowFunction.PartialTukey | WindowFunction.Tukey) && flake.EstimationDepth == 1 && flake.MinPrecisionSearch == 1 && settings.HasDefaultValuesForMode(7)) ? "7" :
(flake.PredictionType == PredictionType.Levinson && flake.StereoMethod == StereoMethod.Evaluate && flake.WindowMethod == WindowMethod.EstimateN && flake.WindowFunction == (WindowFunction.PunchoutTukey | WindowFunction.PartialTukey) && flake.EstimationDepth == 1 && flake.MinPrecisionSearch == 1 && settings.HasDefaultValuesForMode(6)) ? "6" :
(flake.PredictionType == PredictionType.Levinson && flake.StereoMethod == StereoMethod.Evaluate && flake.WindowMethod == WindowMethod.Estimate && flake.WindowFunction == WindowFunction.PunchoutTukey && flake.EstimationDepth == 1 && flake.MinPrecisionSearch == 1 && settings.HasDefaultValuesForMode(5)) ? "5" :
(flake.PredictionType == PredictionType.Levinson && flake.StereoMethod == StereoMethod.Estimate && flake.WindowMethod == WindowMethod.Estimate && flake.WindowFunction == WindowFunction.PunchoutTukey && flake.EstimationDepth == 2 && flake.MinPrecisionSearch == 1 && settings.HasDefaultValuesForMode(4)) ? "4" :
(flake.PredictionType == PredictionType.Levinson && flake.StereoMethod == StereoMethod.Estimate && flake.WindowMethod == WindowMethod.Estimate && flake.WindowFunction == WindowFunction.PunchoutTukey && flake.EstimationDepth == 1 && flake.MinPrecisionSearch == 1 && settings.HasDefaultValuesForMode(3)) ? "3" :
(flake.PredictionType == PredictionType.Levinson && flake.StereoMethod == StereoMethod.Estimate && flake.WindowMethod == WindowMethod.Estimate && flake.WindowFunction == WindowFunction.PartialTukey && flake.EstimationDepth == 1 && flake.MinPrecisionSearch == 1 && settings.HasDefaultValuesForMode(2)) ? "2" :
(flake.PredictionType == PredictionType.Fixed && flake.StereoMethod == StereoMethod.Independent && flake.WindowMethod == WindowMethod.Estimate && flake.MinPrecisionSearch == 1 && flake.WindowFunction == WindowFunction.Tukey && settings.HasDefaultValuesForMode(1)) ? "1" :
(flake.PredictionType == PredictionType.Fixed && flake.StereoMethod == StereoMethod.Independent && flake.WindowMethod == WindowMethod.Estimate && flake.MinPrecisionSearch == 1 && flake.WindowFunction == WindowFunction.Tukey && settings.HasDefaultValuesForMode(0)) ? "0" :
"?"
);
}
#endif