一段写坏掉的快速DCT实现

最新推荐文章于 2021-01-18 19:33:55 发布
最新推荐文章于 2021-01-18 19:33:55 发布
阅读量2.4k 收藏
点赞数
分类专栏： Algorithm and Arithmetic Code Antique Gallery Code Scrolls 音视频编解码与处理 编码偶记 .NET C#
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/quanben/article/details/8719919
版权

想当然了，用递归实现DCT，没想到DCT有4个分支需要递归下去，这样的规模非但无法快速实现，反而由于本身时间复杂度没有多大减少加上递归开销等等比慢速实现往往还慢。

这个代码片段将由于清洁需要从QSharp中删除而保留在这里，对其分析将在代码之后有空时进行。过两天想想是不是能用动态规划或备忘录来改进这个算法。

/// <summary>
///  Type-IV DCT implemented using recursive method
///  Warning: This method is mathematically crippled even slower than the normal one and thereby not suitable for fast caluculation
/// </summary>
/// <param name="input">The input sequence</param>
/// <param name="output">The output sequence</param>
public static void FctIvRecursive(this IList<double> input, IList<double> output)
{
	var n = new DftMethods.FactorAidedInteger(input.Count);
	input.FctIvRecursiveCore(n, 0, 1, output);
}

private static void DctIvFCore(this IList<double> input, int start, int skip, double a, IList<double> output)
{
	var ni = input.Count;
	var n = output.Count;
	for (var k = 0; k < n; k++)
	{
		var d = 0.0;
		var j = 0;
		for (var i = start; i < ni; i += skip, j++)
		{
			var v = input[i] * Math.Cos(a * (k + 0.5) * (j + 0.5));
			if (j % 2 != 0)
			{
				v = -v;
			}
			d += v;
		}
		output[k] = d;
	}
}

private static void DstIvCore(this IList<double> input, int start, int skip, double a, IList<double> output)
{
	var ni = input.Count;
	var n = output.Count;
	for (var k = 0; k < n; k++)
	{
		var d = 0.0;
		var j = 0;
		for (var i = start; i < ni; i += skip, j++)
		{
			d += input[i]*Math.Sin(a*(k + 0.5)*(j + 0.5));
		}
		output[k] = d;
	}
}

private static void DstIvFCore(this IList<double> input, int start, int skip, double a, IList<double> output)
{
	var ni = input.Count;
	var n = output.Count;
	for (var k = 0; k < n; k++)
	{
		var d = 0.0;
		var j = 0;
		for (var i = start; i < ni; i += skip, j++)
		{
			var v = input[i] * Math.Sin(a * (k + 0.5) * (j + 0.5));
			if (j % 2 != 0)
			{
				v = -v;
			}
			d += v;
		}
		output[k] = d;
	}
}


private static void FctIvRecursiveCore(this IList<double> input,
									   DftMethods.FactorAidedInteger n, int start, int skip,
									   IList<double> output)
{
	var n0 = n.N;
	var n1 = n.Decompose();
	var n2 = n.N;
	var a = Math.PI/n0;

	if (n1 > 1 && n2 > 1)
	{
		for (var k = 0; k < n0; k++)
		{
			output[k] = 0;
		}

		var b = 0.25*a;
		for (var r = 0; r < n1; r++)
		{
			int newStart, newSkip;

			DftMethods.GetSubsequenceIndexGenerator(start, skip, r, n1, out newStart, out newSkip);

			// N2-point DCT and DST
			var tmpDct = new double[n2];
			var tmpDst = new double[n2];
			var tmpDctF = new double[n2];
			var tmpDstF = new double[n2];
			input.FctIvRecursiveCore(n, newStart, newSkip, tmpDct);
			input.FstIvRecursiveCore(n, newStart, newSkip, tmpDst);
			input.FctIvFRecursiveCore(n, newStart, newSkip, tmpDctF);
			input.FstIvFRecursiveCore(n, newStart, newSkip, tmpDstF);

			var c = b*(2*r + 1 - n1);
			for (var k = 0; k < n0; k++)
			{
				var k2 = k%n2;
				var l = k/n2;
				var phi = c*(2*k + 1);
				var lMod4 = l%4;
				switch (lMod4)
				{
					case 0:
						output[k] += tmpDct[k2] * Math.Cos(phi) - tmpDst[k2] * Math.Sin(phi);
						break;
					case 2:
						output[k] += -tmpDct[k2] * Math.Cos(phi) + tmpDst[k2] * Math.Sin(phi);
						break;
					case 1:
						output[k] += -tmpDstF[k2] * Math.Cos(phi) - tmpDctF[k2] * Math.Sin(phi);
						break;
					case 3:
						output[k] += tmpDstF[k2] * Math.Cos(phi) + tmpDctF[k2] * Math.Sin(phi);
						break;
				}
			}
		}
	}
	else
	{
		input.DctIvCore(start, skip, a, output);
	}

	n.Restore(n1);
}

private static void FstIvRecursiveCore(this IList<double> input,
									   DftMethods.FactorAidedInteger n, int start, int skip,
									   IList<double> output)
{
	var n0 = n.N;
	var n1 = n.Decompose();
	var n2 = n.N;
	var a = Math.PI/n0;

	if (n1 > 1 && n2 > 1)
	{
		for (var k = 0; k < n0; k++)
		{
			output[k] = 0;
		}

		var b = 0.25*a;
		for (var r = 0; r < n1; r++)
		{
			int newStart, newSkip;

			DftMethods.GetSubsequenceIndexGenerator(start, skip, r, n1, out newStart, out newSkip);

			// N2-point DCT and DST
			var tmpDct = new double[n2];
			var tmpDst = new double[n2];
			var tmpDctF = new double[n2];
			var tmpDstF = new double[n2];
			input.FctIvRecursiveCore(n, newStart, newSkip, tmpDct);
			input.FstIvRecursiveCore(n, newStart, newSkip, tmpDst);
			input.FctIvFRecursiveCore(n, newStart, newSkip, tmpDctF);
			input.FstIvFRecursiveCore(n, newStart, newSkip, tmpDstF);

			var c = b*(2*r + 1 - n1);
			for (var k = 0; k < n0; k++)
			{
				var k2 = k%n2;
				var l = k/n2;
				var phi = c*(2*k + 1);
				var lMod4 = l % 4;
				switch (lMod4)
				{
					case 0:
						output[k] += tmpDst[k2] * Math.Cos(phi) + tmpDct[k2] * Math.Sin(phi);
						break;
					case 2:
						output[k] += -tmpDst[k2] * Math.Cos(phi) - tmpDct[k2] * Math.Sin(phi);
						break;
					case 1:
						output[k] += tmpDctF[k2] * Math.Cos(phi) - tmpDstF[k2] * Math.Sin(phi);
						break;
					case 3:
						output[k] += -tmpDctF[k2] * Math.Cos(phi) + tmpDstF[k2] * Math.Sin(phi);
						break;
				}
			}
		}
	}
	else
	{
		input.DstIvCore(start, skip, a, output);
	}

	n.Restore(n1);
}

private static void FctIvFRecursiveCore(this IList<double> input,
										DftMethods.FactorAidedInteger n, int start, int skip,
										IList<double> output)
{
	var n0 = n.N;
	var n1 = n.Decompose();
	var n2 = n.N;
	var a = Math.PI / n0;

	if (n1 > 1 && n2 > 1)
	{
		for (var k = 0; k < n0; k++)
		{
			output[k] = 0;
		}

		var b = 0.25 * a;
		for (var r = 0; r < n1; r++)
		{
			int newStart, newSkip;

			DftMethods.GetSubsequenceIndexGenerator(start, skip, r, n1, out newStart, out newSkip);

			// N2-point DCT and DST
			var tmpDct1 = new double[n2];
			var tmpDst1 = new double[n2];
			var tmpDct2 = new double[n2];
			var tmpDst2 = new double[n2];
			if (n1%2 == 0)
			{
				input.FctIvRecursiveCore(n, newStart, newSkip, tmpDct1);
				input.FstIvRecursiveCore(n, newStart, newSkip, tmpDst1);
				input.FctIvFRecursiveCore(n, newStart, newSkip, tmpDct2);
				input.FstIvFRecursiveCore(n, newStart, newSkip, tmpDst2);
			}
			else
			{
				input.FctIvRecursiveCore(n, newStart, newSkip, tmpDct2);
				input.FstIvRecursiveCore(n, newStart, newSkip, tmpDst2);
				input.FctIvFRecursiveCore(n, newStart, newSkip, tmpDct1);
				input.FstIvFRecursiveCore(n, newStart, newSkip, tmpDst1);
			}

			var c = b*(2*r + 1 - n1);
			var rIsOdd = r%2 != 0;
			for (var k = 0; k < n0; k++)
			{
				var k2 = k % n2;
				var l = k / n2;
				var phi = c * (2 * k + 1);
				var lMod4 = l % 4;
				var d = 0.0;
				switch (lMod4)
				{
					case 0:
						d = tmpDct1[k2] * Math.Cos(phi) - tmpDst1[k2] * Math.Sin(phi);
						break;
					case 2:
						d = -tmpDct1[k2] * Math.Cos(phi) + tmpDst1[k2] * Math.Sin(phi);
						break;
					case 1:
						d = -tmpDst2[k2] * Math.Cos(phi) - tmpDct2[k2] * Math.Sin(phi);
						break;
					case 3:
						d = tmpDst2[k2] * Math.Cos(phi) + tmpDct2[k2] * Math.Sin(phi);
						break;
				}
				if (rIsOdd)
				{
					d = -d;
				}
				output[k] += d;
			}
		}
	}
	else
	{
		input.DctIvFCore(start, skip, a, output);
	}

	n.Restore(n1);
}

private static void FstIvFRecursiveCore(this IList<double> input,
										DftMethods.FactorAidedInteger n, int start, int skip,
										IList<double> output)
{
	var n0 = n.N;
	var n1 = n.Decompose();
	var n2 = n.N;
	var a = Math.PI / n0;

	if (n1 > 1 && n2 > 1)
	{
		for (var k = 0; k < n0; k++)
		{
			output[k] = 0;
		}

		var b = 0.25 * a;
		for (var r = 0; r < n1; r++)
		{
			int newStart, newSkip;

			DftMethods.GetSubsequenceIndexGenerator(start, skip, r, n1, out newStart, out newSkip);

			// N2-point DCT and DST
			var tmpDct1 = new double[n2];
			var tmpDst1 = new double[n2];
			var tmpDct2 = new double[n2];
			var tmpDst2 = new double[n2];
			if (n1%2 == 0)
			{
				input.FctIvRecursiveCore(n, newStart, newSkip, tmpDct1);
				input.FstIvRecursiveCore(n, newStart, newSkip, tmpDst1);
				input.FctIvFRecursiveCore(n, newStart, newSkip, tmpDct2);
				input.FstIvFRecursiveCore(n, newStart, newSkip, tmpDst2);
			}
			else
			{
				input.FctIvRecursiveCore(n, newStart, newSkip, tmpDct2);
				input.FstIvRecursiveCore(n, newStart, newSkip, tmpDst2);
				input.FctIvFRecursiveCore(n, newStart, newSkip, tmpDct1);
				input.FstIvFRecursiveCore(n, newStart, newSkip, tmpDst1);
			}

			var c = b * (2 * r + 1 - n1);
			var rIsOdd = r%2 != 0;
			for (var k = 0; k < n0; k++)
			{
				var k2 = k % n2;
				var l = k / n2;
				var phi = c * (2 * k + 1);
				var lMod4 = l % 4;
				var d = 0.0;
				switch (lMod4)
				{
					case 0:
						d = tmpDst1[k2] * Math.Cos(phi) + tmpDct1[k2] * Math.Sin(phi);
						break;
					case 2:
						d = -tmpDst1[k2] * Math.Cos(phi) - tmpDct1[k2] * Math.Sin(phi);
						break;
					case 1:
						d = tmpDct2[k2] * Math.Cos(phi) - tmpDst2[k2] * Math.Sin(phi);
						break;
					case 3:
						d = tmpDct2[k2] * Math.Cos(phi) - tmpDst2[k2] * Math.Sin(phi);
						break;
				}
				if (rIsOdd)
				{
					d = -d;
				}
				output[k] += d;
			}
		}
	}
	else
	{
		input.DstIvFCore(start, skip, a, output);
	}

	n.Restore(n1);
}


 

 

quanben
关注
专栏目录
【语音识别】基于mfcc特征模板匹配算法实现语音识别matlab源码含GUI
qq_59747472的博客
09-01 44
在任意一个Automatic speech recognition 系统中，第一步就是提取特征。换句话说，我们需要把音频信号中具有辨识性的成分提取出来，然后把其他的乱七八糟的信息扔掉，例如背景噪声啊，情绪啊等等。 搞清语音是怎么产生的对于我们理解语音有很大帮助。人通过声道产生声音，声道的shape（形状？）决定了发出怎样的声音。声道的shape包括舌头，牙齿等。如果我们可以准确的知道这个形状，那么我们就可以对产生的音素phoneme进行准确的描述。声道的形状在语音短时功率谱的包络中显示出...
【STM32项目】基于STM32的音乐播放器（完整工程资料源码）
最新发布
妄北y
05-17 1275
设计一个基于STM32F103的多功能MP3，需要实现MP3播放、收音机、图片播放、温度计、彩色台灯、功率放大、万年历、闹钟、游戏和电子书等功能。整个系统由音频解码部分、收音机部分、音效处理部分、SD卡部分、功率放大、电源部分、人机交互等几个部分组成，通过一块MCU（STM32F103RBT6）控制运行，各部分协调运作，使用LCD和触摸屏作为人机交互工具，通过彩色液晶屏作为输出，使用触摸屏作为输入。这个系统的设计具有很多挑战性，因为涉及到FAT文件系统解析、JPEG/BMP数据格式解析和触摸屏控制等技
从“零”开始学习一下DCT
golfbears的博客
01-18 704
“汝果欲学诗，工夫在诗外”，这两句诗来形容对DCT的学习最贴切不过了，应为这个技术从DTFT技术衍生出来，具体实现又大量借鉴了FFT的蝶形运算，所以想搞懂DCT，先要拾起久违的傅里叶频率分析理论。
DCT(离散余弦变换)C#源代码
11-11
DCT C# 源代码 对媒体技术 离散余弦变换
MDCT/MDST, IMDCT与DCT-IV的关系以及快速FFT实现
xiaoyaoxiaoshenxian的博客
09-11 3734
本文主要是MDCT/MDST/IMDCT与DCT-IV的关系，然后是采用N/4点快速FFT傅里叶变换来实现N点的MDCT/IMDCT的公式推导以及重叠相加的基本原理，并根据MDCT/IMDCT的特点推导满足完全重建的窗函数的特性。 MDCT/MDST/IMDCT与DCT-IV的关系 如下是DCT-IV，MDCT/MDST/IMDCT/IMDST的公式 MDCT/MDST变
DCT 变换（几个简单的MATLAB的例子）
热门推荐
sxjk1987的专栏
09-08 2万+
原文地址：http://hi.baidu.com/huguosheng/item/06393d5744ee9f3795eb05f2 example 1(use real image): A=imread('class_f.png'); imshow(A) %A is unit8(0,255) C=dct2(A); %进行余弦变换 figure; B=log(abs(C));
DCT快速变换
dd167的专栏
02-15 7383
DCT快速变换作者：陈祖尚 (blog:darnshong)下载源代码一、引言DCT变换是数字图像处理中重要的变换，很多重要的图像算法、图像应用都是基于DCT变换的，如JPEG图像编码方式。对于大尺寸的二维数值矩阵，倘若采用普通的DCT变换来进行，其所花费的时间将是让人难以忍受甚至无法达到实用。而要克服这一难点，DCT变换的快速算法无非是非常吸引人的。就目前而言，DCT
关于操作系统设计的基本原理和设计原则
tugouxp的专栏
03-09 6860
操作系统设计的精妙之处就在于，在底层硬件之上创造了新的抽象， 对于系统初始化来说， 它呈现了执行流到线程转化这一概念，该概念远远比它的实现细节更加重要，处理器以“取指-执行”为周期开始串行执行指令，而初始化代码将自身转化为一个并发处理系统， 这里的关键之处在于，初始化代码并不是创建一个独立的并发系统，然后跳转到新的系统。
50、流式音视频之一（应用层）
跃寒的博客
03-20 1万+
引言 Web 应用和移动Web 并不是网络应用领域中唯一令人振奋的发展。对许多人来说，音频和视频才是网络的圣杯。当提及“多媒体”这个词时，无论技术人员和商人都会不约而 同地垂涎三尺。前者看到的是为每一台计算机提供E 语音和视频点播所隐含的巨大技术挑 战，而后者则看到它同样隐含的丰厚的利润。 虽然至少从20 世纪70 年代开始就有了通过Internet 发送音频和视频的想法，但大约 自2000 开...
AE自动变速器拆解.pptx
12-24
5. **双离合器变速器(DCT/DSG)**：DCT是一种结合了手动变速器的灵活性和自动变速器的舒适性的变速器，通过两套离合器分别控制奇数挡和偶数挡，实现快速换挡，无动力中断。 **丰田A340E自动变速器详解** 丰田A340E...
离散余弦变换（DCT）的C++实现
BigDream123的博客
02-19 5066
DCT变换具体概念性质参考博客https://blog.csdn.net/BigDream123/article/details/101426393 这里仅仅描述DCT-2变换的C++实现。 DCT-2正变换的公式如下： DCT-2反变换公式如下： #include<iostream> #include<stdio.h> #include<math...
DCT算法C++实现源码
厦门智创科技
08-07 1051
声明：此代码并非原创，参照网络上一位前辈的代码，但是测试发现他的代码有BUG，修改后贴出已做记录; #include<stdlib.h> #include<iostream> #include<math.h> #include<vector> #define PI acos(-1)//定义圆周率 using namespace std; v...
DCT离散余弦变换C版本
Calcular的博客
07-08 2308
#include #include using namespace std; #define PI 3.14159265354 #define DCT_SIZE 8 float DCT_Mat[DCT_SIZE][DCT_SIZE]; float DctMap[DCT_SIZE][DCT_SIZE]= { 130,130,130,129,134,133,129,130, 130
c++实现矩阵直接做dct变换
zyckhuntoria的博客
08-06 1926
#include &lt;iostream&gt; #include &lt;string.h&gt; #include &lt;opencv2/core/core.hpp&gt; #include &lt;opencv2/highgui/highgui.hpp&gt; #include &lt;opencv2\core\mat.hpp&gt; #include&lt;opencv2...
某阿里巴巴笔试题集第23题
quanben
10-30 1万+
废话少说了，原题集出处：http://blog.csdn.net/hackbuteer1/article/details/6902917 题目： 一个骰子，6面，1个面是 1， 2个面是2， 3个面是3， 问平均掷多少次能使1、2、3都至少出现一次。 方法： 面对面试概率题几乎屡试不爽的分叉树递归列方程法。 这是一个求数学期望的问题，最终是求1，2，3出现至少一次的最短长度的期望。 这
ADA语言Hello World
quanben
11-11 1万+
为什么要写关于这么古老的语言（1980年左右）的文章呢？因为我想既然mainframe用的cobol这么恶劣的语言都写了那么多了，为何不写一点我比较喜欢的语言呢？虽然除了高安全级控制和仿真领域之外ADA是不怎么用的。 首先说两句我个人对ADA语言的态度： 从编写大型程序的OO语言而言（我基本主要用的也就是这类语言），ADA在我评价中的排名是和C#和C++一个级别上的。ADA很大程度上是Java
给我写一段用DCT实现图像压缩的matlab代码
06-10
以下是一个基于DCT的图像压缩Matlab代码示例： ```matlab % 读入需要压缩的图像 img = imread('lena.png'); [m, n] = size(img); % 设置压缩比 ratio = 0.5; % 计算需要保留的系数数 num_coeffs = round(m * n * ratio); % 对图像进行DCT变换 dct_img = dct2(img); % 将DCT系数按大小排序，并仅保留前num_coeffs个系数 sorted_coeffs = sort(abs(dct_img(:)), 'descend'); thresh = sorted_coeffs(num_coeffs); dct_img(abs(dct_img) < thresh) = 0; % 对DCT系数进行逆变换，得到压缩后的图像 compressed_img = idct2(dct_img); % 显示压缩前后的图像 subplot(1,2,1); imshow(img); title('Original Image'); subplot(1,2,2); imshow(compressed_img); title(['Compressed Image (', num2str(round(ratio*100)), '%)']); ``` 这段代码会读入一个名为"lena.png"的图像，对其进行DCT变换，然后按照设定的压缩比例保留前num_coeffs个系数，并将其余系数置为0。最后对变换后的系数进行逆变换，得到压缩后的图像。
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值