C++ 和 OpenCV 实现卷积神经网络并加载 Keras 训练好的参数进行预测

2023-11-02

一. 背景

我们会经常在 Python 的环境下进行神经网络的搭建和训练, 特别是 Keras 这样的高级库用起来很是方便, 分分钟可以弄一个神经网络出来. 一些小的网络训练也是分分钟的事情. But, 很多应用是用 C++ 之类的语言开发的, 而且 OpenCV 对 C++ 的支持也是极好的. 怎样在 C++ 环境下利用 Keras 训练好的模型和参数进行预测或者分类着实成了一些同学的难题, 今天给大家讲一个用 C++ 实现一个 卷积神经网络 模型, 然后再加载训练好的参数的例子. 这样就不用通过 Python 或者 TensorFlow 的 C++ 接口来调用模型, 而是直接从 HDF5 文件中读取模型参数. 看完这个例子之后, 你应该能够写一个 C++ 类直接加载模型和参数进行预测了

二. Keras 定义神经网络结构

MNIST 作为一个 Hello word 级的入门例子, 作为本文的例子也是极好的. 具体实现不是本文要讲的主题, 所以这里就不讲怎么实现了, 直接上一个卷积神经网络, 实现的模型如下. 不用理会为什么这样设计结构, 或者正确率是多少, 这个只是例子, 能说明问题就行

project_name = "mnist"
num_classes = 10

# 因为OpenCV的图像是数据是channels_last所以这里的input_size通道数放到最后
model = keras.Sequential(name = project_name)

model.add(keras.layers.Conv2D(32, kernel_size = (3, 3), activation = "relu",
                              input_shape = input_shape, name = "conv_1"))
model.add(keras.layers.MaxPool2D(pool_size = (2, 2), strides = (2, 2),
                                 name = "max_pool_1"))

model.add(keras.layers.Conv2D(64, kernel_size = (3, 3), activation = "relu",
                              name = "conv_2"))
model.add(keras.layers.MaxPool2D(pool_size = (2, 2), strides = (1, 1),
                                 name = "max_pool_2"))

model.add(keras.layers.Flatten(name = "flatten"))
model.add(keras.layers.Dense(256, activation = "relu", name = "dense_1"))
model.add(keras.layers.Dense(num_classes, activation = "softmax", name = "output"))

model.compile(optimizer = "adam",
              loss = "sparse_categorical_crossentropy",
              metrics = ["accuracy"])
  
model.summary()

_________________________________________________________________
Model: "mnist"
_________________________________________________________________
Layer (type)                 Output Shape              Param #   
=================================================================
conv_1 (Conv2D)              (None, 26, 26, 32)        320       
_________________________________________________________________
max_pool_1 (MaxPooling2D)    (None, 13, 13, 32)        0         
_________________________________________________________________
conv_2 (Conv2D)              (None, 11, 11, 64)        18496     
_________________________________________________________________
max_pool_2 (MaxPooling2D)    (None, 10, 10, 64)        0         
_________________________________________________________________
flatten (Flatten)            (None, 6400)              0         
_________________________________________________________________
dense_1 (Dense)              (None, 256)               1638656   
_________________________________________________________________
output (Dense)               (None, 10)                2570      
=================================================================
Total params: 1,660,042
Trainable params: 1,660,042
Non-trainable params: 0
_________________________________________________________________

你发现没有, 在上面的每一层我都填了 name 这个参数, 这个参数是可以省略的, 如果你有强迫症的话, 就把它填了, 要不然每次 model.summary() 出来的名称可能都不一样

channels_first 与 channels_last

所谓 channels_last 就是通道维度在各个 shape 的最后一维, 比如上面的 (None, 26, 26, 32), 32 个 Filter 也就是 32 个通道的意思, depth 就是 32

channels_last 的数据排列特点是相邻两个数据先按 channel 变化. 比如 OpenCV 的 Mat, 三通道的 RGB 图像数据在内存中的排列是 BGRBGRBGR… …, 而不是先把一张图像的第一个通道全部数据排完之后再排后面的通道, 后者的排列方式是 BBB…BBB|GGG…GGG|RRR…RRR, 而这种方式叫 channels_first

channels_first 与 channels_last 转换

只需要改变一下维度信息就可以了

channels_first 转换成 channels_last

if bkd.image_data_format() == 'channels_first':
    x_train = x_train.reshape(x_train.shape[0], 1, img_rows, img_cols)
    x_test = x_test.reshape(x_test.shape[0], 1, img_rows, img_cols)
    
    // 这里是转换代码, 将第 1 个维度移动到 第 3 个维度位置(0 开始)
    x_train = np.moveaxis(x_train, 1, 3)
    x_test = np.moveaxis(x_test, 1, 3)
    
    input_shape = (img_rows, img_cols, 1)

channels_last 转换成 channels_first

if bkd.image_data_format() == 'channels_last':
    x_train = x_train.reshape(x_train.shape[0], img_rows, img_cols, 1)
    x_test = x_test.reshape(x_test.shape[0], img_rows, img_cols, 1)
    
    // 这里是转换代码, 将第 3 个维度移动到 第 1 个维度位置(0 开始)
    x_train = np.moveaxis(x_train, 3, 1)
    x_test = np.moveaxis(x_test, 3, 1)
    
    input_shape = (1, img_rows, img_cols)

三. 用 C++ 和 OpenCV 实现网络结构

如果你想用 C++ 实现神经网络结构并进行训练的话, 这个难度还是有一点大的, 但是, 如果你只是实现结构并进行预测的话, 难度一下子就降低到地面了

输入图像处理

虽然输入是单通道图像, 但是我还是按三通道来写, 这样方便扩展

const int img_rows = 28;
const int img_cols = 28;

Mat img_src; // 输入图像, 这个是你传递进来的需要识别的图像, 假设是 CV_U8 的类型,

// 如果输入图像的尺寸不是 28 * 28, 则需要改变其尺寸
if (img_rows != img_src.rows || img_cols != img_src.cols)
{
	resize(img_src, img_src, Size2i(img_cols, img_rows), 0, 0, CV_INTER_CUBIC);
}

// 因为是以三通道来写, 所以要将每个通道分离出来才能卷积
vector<Mat> img_input;

if (3 == img_src.channels())
{
	// 将图像转换成 0~1 的范围
	img_src.convertTo(img_src, CV_32FC3, 1.0 / 255.0, 0);

	split(img_src, img_input);
}
else
{
	// 将图像转换成 0~1 的范围
	img_src.convertTo(img_src, CV_32FC1, 1.0 / 255.0, 0);

	img_input.push_back(img_src);
}

卷积操作要点

Filter 的通道数与输入图像的通道数是一样的. 卷积时, 输入图像每个通道用 Filter 对应的通道来卷积, 比如 RGB 三通道的图像卷积时, R 通道用 Filter 的 R 通道来卷积, 其他两个通道也一样对应. 就算是100个通道也是相同的操作

最最重要的一点是每个通道分别卷积完成后, 要把这些卷积结果图像相加变成一个通道

煮个粟子, 一个 Filter 有 32 个通道, 卷积完成后便有 32 张结果图像, 将这 32 张结果图像相加变成一张图像. 如果有 64 个 Filer, 这样的操作完成后就有 64 张图像, 所以输出图像通道数就是 64, 这就是第三层网络的输出 Output Shape = (None, 11, 11, 64) 的原因

第一层 (conv_1)

这一层就是卷积操作, 有 32 个 3 * 3 大小的 Filter , 因为输入图像是单通道的, 所以每个 Filter 是单通道的. 上面已经讲过了, 还是以三通道的情况来写, 这样的写法也适用于单通道输入, 代码中的 const int ker_channels_1 = img_input.size() 就可以指明通道数
想一下, OpenCV 里面已经有卷积操作的函数了, 所以可以直接调用

// 后缀 _1 表示第 1 层的参数, 依次类推
const int kernels_1 = 32;
const int ker_channels_1 = img_input.size(); // 用于指明输入图像的通道数

// 卷积核, 3 * 3大小的矩阵, 因为神经网络的定义里面就是 3 * 3 的大小
// 这里 mKernels_1 没有初始化, 先定义在这里说明问题, 要到后面加载训练好的参数时初始化, 也就是训练好的权重(Weight)
// 这里定义成 vector<vector<Mat>> 类型, 可以将输入扩展为 3 通道的情况
vector<vector<Mat>> mKernels_1(kernels_1);

// 卷积后的结果图像
vector<Mat> imgFilter2d_1(kernels_1);

for (int i = 0; i < kernels_1; i++)
{
	for (int j = 0; j < ker_channels_1; j++)
	{
		Mat img_tmp;

		filter2D(img_input[j], img_tmp, mKernels_1[i][j].depth(), mKernels_1[i][j]);

		// 多通道的图像卷积后, 同一个 Filter(多通道的 Filter 也算一个 Filer) 的卷积结果要相加成变成一张图像
		// 不然会造成输出与下一层输入不匹配
		// 这里也适用于单通道图像, 因为单通道图像 ker_channels_1 == 1, j 只能取 0
		if (0 == j)
		{
			imgFilter2d_1[i] = img_tmp;
		}
		else
		{
			imgFilter2d_1[i] += img_tmp;
		}
	}
}

经过上面的步骤之后, 就得到了 kernels_1(==32) 张卷积后的图像(为什么是 32 张?).

但是还没有完成. 在定义网络时有一个东东没有定义, 可是在网络里面默认又是有的, 那就是 Bias, 这个东东是用来抑制或者增强某一个神经元或者 Filter 的, 每个一个 Filter 对应了一个 Scaler(标量)的 bias, 32 个 Filter 就对应了一个 Bias 的向量, 所以还要在每个卷积后的图像上面加上这个 Bias . 上面的代码要修改成如下

const int kernels_1 = 32;
const int ker_channels_1 = img_input.size();

vector<vector<Mat>> mKernels_1(kernels_1);
vector<Mat> imgFilter2d_1(kernels_1);

for (int i = 0; i < kernels_1; i++)
{
	for (int j = 0; j < ker_channels_1; j++)
	{
		Mat img_tmp;

		filter2D(img_input[j], img_tmp, mKernels_1[i][j].depth(), mKernels_1[i][j]);

		if (0 == j)
		{
			imgFilter2d_1[i] = img_tmp;
		}
		else
		{
			imgFilter2d_1[i] += img_tmp;
		}
	}

	imgFilter2d_1[i] += mBias_1.at<float>(i, 0);
}

是不是很简单
从上面的 model.summary() 中可以看出, 第一层 Output Shape 是 (None, 26, 26, 32), 怎么解释?
None: 表示样本数量暂时不知道
26, 26: 表示卷积之后的图像大小, 输入图像是28 * 28, 因为卷积核是 3 * 3, 所以要去掉边缘的一个像素变成 26 * 26
32: 输出图像的通道数, 也就是 Filter 的数量
所以卷积后有 32 张图像(32 个通道, channels_last), 每张大小为 26 * 26, 也说明了卷积操作没有进行图像边缘的填充, 因为在定义网络时没有指明要填充

还要解释一点, 后面的 Param #, 就是参数数量, 第一层是 320, 怎么得来的?
计算式: Param # = 32 * 3 * 3 + 32(Bias) = 320, 这个不用解释吧

不要以为第一层就这么定义完成了, 还有一点就是激活函数, 这里我选择的是 relu, 用 OpenCV 二值化函数就可以搞定

const int kernels_1 = 32;
const int ker_channels_1 = img_input.size();

vector<vector<Mat>> mKernels_1(kernels_1);
vector<Mat> imgFilter2d_1(kernels_1);

for (int i = 0; i < kernels_1; i++)
{
	for (int j = 0; j < ker_channels_1; j++)
	{
		Mat img_tmp;

		filter2D(img_input[j], img_tmp, mKernels_1[i][j].depth(), mKernels_1[i][j]);
		
		if (0 == j)
		{
			imgFilter2d_1[i] = img_tmp;
		}
		else
		{
			imgFilter2d_1[i] += img_tmp;
		}
	}
	
	imgFilter2d_1[i] += mBias_1.at<float>(i, 0);
	threshold(imgFilter2d_1[i], imgFilter2d_1[i], 0, 1, CV_THRESH_TOZERO);
}

现在才算定义完成了
如果你想可视化结果的话, 只要把卷积结果当成图像显示出来就可以了, 只是要归一化为 0~255 再转换成 CV_8U 才能正常显示, 其他层的输出也一样

Mat img_show;

for (int i = 0; i < kernels_1; i++)
{
	normalize(imgFilter2d_1[i], img_show, 0, 255, CV_MINMAX);
	img_show.convertTo(img_show, CV_8UC1);

	char show_name[128] = {0};
	printf_s("img_show-%d", i + 1);
		
	namedWindow(show_name, WINDOW_NORMAL);
	imshow(show_name, img_show);
}

第一层输出的前8张图像
visualize

第二层 (max_pool_1)

第一层的输出是 32 通道的图像, 作为第二层的输入, 第二层是一个 MaxPool2D 层, 滑动窗口大小为 2 * 2, stride 为 2, 意思就是在 2 * 2 大小的范围内取最大值, 滑动步长为 2, 这个实现起来也很简单, 但是有一个坑要小心, 因为在第一层里边缘没有填充, 所以图像尺寸有所减小, 但是 OpenCV filter2D 会对边缘进行填充, 所以图像尺寸并没有减小, So, for 循环的起始和结束要改一下

// 滑动步长
const int stride_2 = 2;

// 输出图像, 通道数不变, 只是 Size 不一样
vector<Mat> imgPool2d_2(kernels_1);

for (int i = 0; i < kernels_1; i++) // 第一重循环, 通道
{
	// 第二层不需要加载参数, 所以直接创建图像
	// 行和列都减 2 是因为前一层没有填充, 并且和 Filter 的尺寸有关, 右移一位是因为 stride_2 == 2
	imgPool2d_2[i].create((imgFilter2d_1[i].rows - 2) >> 1, (imgFilter2d_1[i].cols - 2) >> 1, CV_32FC1);
	
	// 起始从 1 开始, 结束为 imgFilter2d_1[i].rows - 1
	for (int r = 1; r < imgFilter2d_1[i].rows - 1; r += stride_2) // 第二重循环, 图像每一行
	{
		for (int c = 1; c < imgFilter2d_1[i].cols - 1; c += stride_2) // 第三重循环, 图像每一列
		{
			const Mat m = imgFilter2d_1[i](cv::Rect(c, r, 2, 2));

			double max_val = 0;
			minMaxLoc(m, nullptr, &max_val);

			// 因为输出图像尺寸减小一半, 所以坐标要除以2
			imgPool2d_2[i].at<float>((r - 1) >> 1, (c - 1) >> 1) = (float)max_val;
		}
	}
}

Pooling 之后, 图像的大小缩小一半, 变成 13 * 13 * 32, 和上面的 model.summary() 也能对应起来, 这一层是没有参数的, 因为只是一个取最大值操作, 第二层也就完成了

第三层 (conv_2)

第三层和第一层原理差不多, 有 64 个 Filter , 大小还是一样, 但是每一个 Filter 的通道数变成了 32, 因为这些 Filter 都要在输入图像上做卷积操作, 而输入图像有 32 个通道, 前面已经解释过了. 一个 Filter 卷积完成后, 结果会有 32 个通道, 要把这 32 个通道相加变成一张图像, 这样才能保证 64 个 Filter 输出就有 64 个通道, 第一层虽然是单通道 , 但是也是这个逻辑. 上面的参数个数计算式为: 64 * 3 * 3 * 32 + 64 = 18496, 这个怎么解释? 你自己想一下吧

对比第一层, 第三层的代码几乎是一样的, 只是把数量改一下

const int kernels_3 = 64;
const int ker_channels_3 = kernels_1;

Mat mBias_3;
vector<vector<Mat>> mKernels_3(kernels_3);
vector<Mat> imgFilter2d_3(kernels_3);

for (int i = 0; i < kernels_3; i++)				// 第一重循环, Filte r的数量
{
	for (int j = 0; j < ker_channels_3; j++)	// 第二重循环, Filter 的通道, 也就是输入图像的通道数
	{
		Mat img_tmp;

		// imgPool2d_2[j] 为上一层的输出
		filter2D(imgPool2d_2[j], img_tmp, mKernels_3[i][j].depth(), mKernels_3[i][j]);

		// 这里就是把每个通道相加
		if (0 == j)
		{
			imgFilter2d_3[i] = img_tmp;
		}
		else
		{
			imgFilter2d_3[i] += img_tmp;
		}
	}
	
	imgFilter2d_3[i] += mBias_3.at<float>(i, 0);

	threshold(imgFilter2d_3[i], imgFilter2d_3[i], 0, 1, CV_THRESH_TOZERO);
}

就这样, 第三层就完成了

第四层 (max_pool_2)

第四层和第二层也差不多, 只是通道数多了一点

const int stride_4 = 1;
// 输出图像
vector<Mat> imgPool2d_4(kernels_3);

for (int i = 0; i < kernels_3; i++)									  // 第一重循环, 通道
{
	imgPool2d_4[i].create(imgFilter2d_3[i].rows - 3, imgFilter2d_3[i].cols - 3, CV_32FC1);
	
	// 这里 r < imgFilter2d_3[i].rows - 2 是因为前一层的卷积把图像尺寸减小了 2
	// 这一层的滑动步长是 1, 所以还要减 1, 一共减少了 3, 所以起始为 1, 结束为 imgFilter2d_3[i].rows - 2
	for (int r = 1; r < imgFilter2d_3[i].rows - 2; r += stride_4)     // 第二重循环, 图像每一行
	{
		for (int c = 1; c < imgFilter2d_3[i].cols - 2; c += stride_4) // 第三重循环, 图像每一列
		{
			const Mat m = imgFilter2d_3[i](cv::Rect(c, r, 2, 2));

			double max_val = 0;
			minMaxLoc(m, nullptr, &max_val);

			// 因为输出图像尺寸只减小了 1, 就不除以 2 了
			imgPool2d_4[i].at<float>(r - 1, c - 1) = (float)max_val;
		}
	}
}

第四层结束

第五层 (flatten)

这一层只是把前一层的输出展开成一个列向量, 前一层的输出有 6400 个, 希望你能算出来为什么是 6400. 展开的时候只是把 64 个通道的图像依次展开就可以了, 展开后不需要通过激活函数

这一层有一个比较大的坑, 就是要注意区分是 channels_first 还是 channels_last. 如果反了, 或者展开错了, 你会死在这层上, 找几天都找不到问题, 不要问我怎么知道的

下面的代码写两个版本, 一个是 channels_first, 另一个是 channels_last, 不同点就循环里的数据处理

// channels_first
const int data_size = imgPool2d_4[0].rows * imgPool2d_4[0].cols * kernels_3;
float *flat_data = new float[data_size];

int k = 0;

for (int i = 0; i < kernels_3; i++)
{
	for (int r = 0; r < imgPool2d_4[i].rows; r++)
	{
		const float *pData = imgPool2d_4[i].ptr<float>(r);

		for (int c = 0; c < imgPool2d_4[i].cols; c++)
		{
			flat_data[k++] = pData[c];
		}
	}
}

// 这个就是列向量
Mat mFlat(data_size, 1, CV_32FC1, flat_data);

delete []flat_data;
flat_data = nullptr;
// channels_last, 以下两种方法选一种
// 方法 1:--------------------------------------------
Mat multi_channel;
merge(imgPool2d_4, multi_channel);

Mat mFlat(multi_channel.rows * multi_channel.cols * multi_channel.channels(), 1,
		CV_32FC1, multi_channel.data);

// 方法 2:--------------------------------------------
const int data_size = imgPool2d_4[0].rows * imgPool2d_4[0].cols * kernels_3;
float *flat_data = new float[data_size];

for (int i = 0; i < kernels_3; i++)
{
	for (int r = 0; r < imgPool2d_4[i].rows; r++)
	{
		for (int c = 0; c < imgPool2d_4[i].cols; c++)
		{
			const int k =
					c * kernels_3 +
					r * kernels_3 * imgPool2d_4[i].cols +
					i;
				flat_data[k] = imgPool2d_4[i].at<float>(r, c);
		}
	}
}

// 这个就是列向量
Mat mFlat(data_size, 1, CV_32FC1, flat_data);

delete []flat_data;
flat_data = nullptr;

展开之后就是第六层的输入了

第六层 (dense_1)

这一层是全连接层, 就是矩阵乘法, 输出为一个矩阵乘以一个列向量再加 Bias, 公式如下, 如果不明白的话, 可能看要回炉一下线性代数的矩阵乘法了
o u t p u t 6 = r e l u ( M 6 ∗ v 5 + B i a s 6 ) output_6= relu(M_6 * v_5 + Bias_6) output6=relu(M6v5+Bias6)
M 6 M_6 M6: 第五层到第六层的 Weight
v 5 v_5 v5: 第五层展开的那个向量
B i a s 6 Bias_6 Bias6: 第六层的偏置
r e l u relu relu: 激活函数

用代码表示如下

// 这里也没有初始化, 希望你能知道
// mWeight_6 的 rows == 256(第六层的神精元数量), cols == data_size
// mBias_6 的 rows == 256(第六层的神精元数量), cols == 1
Mat mWeight_6;
Mat mBias_6;

Mat output_6 = mWeight_6 * mFlat + mBias_6;
threshold(output_6, output_6, 0, 1, CV_THRESH_TOZERO);

第七层 (output)

这一层和第六层不一样的是激活函数 softmax, softmax 使输出变得复杂一点, 公式如下
y i = e o i ∑ e o i + b i y_i = \frac{e^{o_i}}{\sum e ^{o_i + b_i}} yi=eoi+bieoi
o i o_i oi: 第 i i i 个输出(此时还没有经过激活函数)
b i b_i bi: 第 i i i 个偏置
y i y_i yi: 第 i i i 个输出值

const int num_classes = 10;

// mWeight_7 的 rows == num_classes, cols == mBias6.rows
// mBias_7 的 rows == num_classes, cols == 1

Mat mWeight_7;
Mat mBias_7;

Mat output_7 = mWeight_7 * output_6 + mBias_7;

float sum = 0;
vector<float> tmp(output_7.rows);

for (int i = 0; i < output_7.rows; i++)
{
	tmp[i] = exp(output_7.at<float>(i, 0));
	sum += tmp[i];
}

for (int i = 0; i < output_7.rows; i++)
{
	output_7.at<float>(i, 0) = tmp[i] / sum;
}

// max_pt.y就是最终的类别序号了
Point2i max_pt;
minMaxLoc(output_7, nullptr, nullptr, nullptr, &max_pt);

到这里, 神经网络的结构已经完成, 现在只差加载 Weight 了

其实你都可以根据这个思路来实现一个类一层一层增加神经网络层, 就和 Keras 差不多, 因为不需要训练, 只是一层一层地计算得到输出而已

四. 加载 Keras 训练完成的参数

加载参数的原理可以参考 C++ 从 HDF5 文件读取 Keras 神经网络模型和参数

加载第一层

第一层需要加载的是卷积 Filter 和 Bias 的参数, 其定义如下

// kernels_1 = 32;
Mat mBias_1;
vector<Mat> mKernels_1(kernels_1); 

打开文档并加载数据

H5File file("文件路径, 你需要自己修改, 文件名包括.扩展名", H5F_ACC_RDONLY);

// 这里的路径可能和你的不一样, 你可以先用 HDF5 View 打开看一下
// 后面的函数的实现方式会自动识别路径, 这里只是讲原理
DataSet dset = file.openDataSet("/conv_1/conv_1/kernel:0");

// Dataset中数据的数据类型
H5::DataType dt = dset.getDataType();
H5T_class_t t = dt.getClass();

// 数据在内存中的字节数除以数据类型得到 buf 的大小
hsize_t data_Len = dset.getInMemDataSize() / sizeof(float);
float *buf = new float[data_Len];

// 将数据读到 buf
dset.read(buf, dt);

// 这里需要把数据拆分到每一个 mKernels_1
// ----拆分代码----

// ----拆分代码----

// 不用之后要关闭
dt.close();
dset.close();

delete []buf;
buf = nullptr;

上面将数据读到 buf 之后, 需要把数据拆分到 mKernels_1 每一个元素, 我把拆分写在下面

划重点, buf 中的数据是按以下方式排列的

第 1 个 Filter 第 1 行第 1 列第 1 通道 | 第 2 个 Filter 第 1 行第 1 列第 1 通道… | 第 m 个 Filter 第 1 行第 1 列第 1 通道
第 1 个 Filter 第 1 行第 1 列第 2 通道 | 第 2 个 Filter 第 1 行第 1 列第 2 通道… | 第 m 个 Filter 第 1 行第 1 列第 2 通道

第 1 个 Filter 第 1 行第 1 列第 n 通道 | 第 2 个 Filter 第 1 行第 1 列第 n 通道… | 第 m 个 Filter 第 1 行第 1 列第 n 通道

第 1 个 Filter 第 1 行第 2 列第 1 通道 | 第 2 个 Filter 第 1 行第 2 列第 1 通道… | 第 m 个 Filter 第 1 行第 2 列第 1 通道
第 1 个 Filter 第 1 行第 2 列第 1 通道 | 第 2 个 Filter 第 1 行第 2 列第 1 通道… | 第 m 个 Filter 第 1 行第 2 列第 2 通道

第 1 个 Filter 第 1 行第 2 列第 n 通道 | 第 2 个 Filter 第 1 行第 2 列第 n 通道… | 第 m 个 Filter 第 1 行第 2 列第 n 通道

从以上可以看出, 的确是 channels_last 的排列方式, 最先变化的是 Filter 序号, 然后是通道, 再然后是列, 最后是行

那要怎样把这些数据放到 Mat 中去, 只需要 4 重循环就可以办到. 这里没有区分 channels_first 和 channels_last 因为 Keras 都会保存为 channels_last

// 这里的代码是上面代码的拆分部分
for (int i = 0; i < kernels_1; i++)	// 第一重: Filter 序号
{
	mKernels_1[i].resize(ker_channels_1);

	for (int j = 0; j < ker_channels_1; j++)	// 第二重: channel 序号
	{
		mKernels_1[i][j].create(3, 3, CV_32FC1);
			
		for (int r = 0; r < mKernels_1[i][j].rows; r++)	// 第三重: Filter 行序号
		{
			for (int c = 0; c < mKernels_1[i][j].cols; c++)	// 第四重: Filter 列序号
			{
				const int k = 
					c * kernels_1 * ker_channels_1 +
					r * kernels_1 * ker_channels_1 * mKernels_1[i][j].cols +
					j * kernels_1 +
					i;

				mKernels_1[i][j].at<float>(r, c) = buf[k];
			}
		}
	}
}

kernel:0 加载完成后, 再加载 bias:0, 因为上面已经有打开文件的代码了, 这里就直接放码过来了

dset = file.openDataSet("/conv_1/conv_1/bias:0");
dt = dset.getDataType();
t = dt.getClass();

// 数据在内存中的字节数除以数据类型得到 buf 的大小
data_Len = dset.getInMemDataSize() / sizeof(float);
buf = new float[data_Len];

dset.read(buf, dt);

// Bias 只是一个向量, 所以比较简单
mBias_1.create(kernels_1, 1, CV_32FC1);
memcpy(mBias_1.data, buf, dset.getInMemDataSize());

dt.close();
dset.close();

delete []buf;
buf = nullptr;

就这样, 第一层的参数就加载完成了, 其他层的参数加载也差不多, 只是需要注意 Dataset 的路径和名字, 矩阵的尺寸和通道数, 第三层和第一层的方式一样, 这里再演示第六层怎么加载

加载第六层

第六层是一个矩阵和一个向量(也是矩阵)

// data_size 是第五层展开的长度, 这里的维度顺序要反一下, 变成, data_size, 256, 加载参数后再转置, 不然数据对不上
// mWeight_6 和 mBias_6 维度顺序都反了, 这样是为了加载的数据是正确的, 上面的加载第一层没有反是因为 1 维的向量不会影响
// 到他用函数时, 都统一反一下
Mat mWeight_6(data_size, 256, CV_32FC1);
Mat mBias_6(1, 256, CV_32FC1);

dset = file.openDataSet("/dense_1/dense_1/kernel:0");
dt = dset.getDataType();
t = dt.getClass();

data_Len = dset.getInMemDataSize() / sizeof(float);
buf = new float[data_Len];

dset.read(buf, dt);
memcpy(mWeight_6.data, buf, dset.getInMemDataSize());

// 转置交换行和列
transpose(mWeight_6, mWeight_6);

dt.close();
dset.close();

delete []buf;
buf = nullptr;

dset = file.openDataSet("/dense_1/dense_1/bias:0");
dt = dset.getDataType();
t = dt.getClass();

data_Len = dset.getInMemDataSize() / sizeof(float);
buf = new float[data_Len];

dset.read(buf, dt);
memcpy(mBias_6.data, buf, dset.getInMemDataSize());

// 转置交换行和列
transpose(mBias_6, mBias_6);

dt.close();
dset.close();

delete []buf;
buf = nullptr;

最后用完记得 rg.close() 和 file.close();
其他的层加载参数也是一样的原理, 加载参数要在模型定义之前, 模型定义完成就可以用于预测了

五. 加载函数

上面只是讲了步骤, 你可以定义函数或者类简化过程, 完整过程可参考 C++ 和 OpenCV 实现 Keras Sequential 网络. 这里做两个函数示例, 代码放到了 六. 完整过程 里面了. 里面自动计算 Dataset 路径的方法最好是根据保存的网络模型(model)中的路径去读取, 这样才能自适应其他网络结构, 读取网络结构可定义参考 C++ 从 HDF5 文件读取 Keras 神经网络模型和参数

两个加载函数的使用示例

LoadWeight(file, "conv_1", mKernels_1, mBias_1);
LoadWeight(file, "output", mWeight_7, mBias_7);

六. 完整过程

工程的配置可参考 C\C++ 写 HDF5 文件示例C++ 从 HDF5 文件读取 Keras 神经网络模型和参数

// cnn_test.cpp : 定义控制台应用程序的入口点。
//

#include "stdafx.h"

#include <stdint.h>
#include <hdf5.h>
#include <H5Cpp.h>
#include <iostream>

#include <opencv2\opencv.hpp>
#include <opencv2\features2d\features2d.hpp>
#include <opencv2\nonfree\features2d.hpp>

using namespace H5;
using namespace cv;
using namespace std;

#ifdef _DEBUG
#pragma comment(lib, "hdf5_D.lib")
#pragma comment(lib, "hdf5_cpp_D.lib")
#else
#pragma comment(lib, "hdf5.lib")
#pragma comment(lib, "hdf5_cpp.lib")
#endif

bool LoadWeight(const H5File & file, const char * layer_name, vector<vector<Mat>> & kernel, Mat & bias);
bool LoadWeight(const H5File & file, const char * layer_name, Mat & weight, Mat & bias);

int _tmain(int argc, _TCHAR* argv[])
{
	const int img_rows = 28;
	const int img_cols = 28;

	// 图像路径你自己写, 包括.扩展名称
	Mat img_src = imread("F:\\Tmp\\num.bmp", CV_LOAD_IMAGE_UNCHANGED);

	// 如果输入图像的尺寸不是 28 * 28, 则需要改变其尺寸
	if (img_rows != img_src.rows || img_cols != img_src.cols)
	{
		resize(img_src, img_src, Size2i(img_cols, img_rows), 0, 0, CV_INTER_CUBIC);
	}

	// 因为是以三通道来写, 所以要将每个通道分离出来才能卷积
	vector<Mat> img_input;

	if (3 == img_src.channels())
	{
		// 将图像转换成 0~1 的范围
		img_src.convertTo(img_src, CV_32FC3, 1.0 / 255.0, 0);

		split(img_src, img_input);
	}
	else
	{
		// 将图像转换成 0~1 的范围
		img_src.convertTo(img_src, CV_32FC1, 1.0 / 255.0, 0);

		img_input.push_back(img_src);
	}

	const int kernels_1 = 32;
	const int ker_channels_1 = img_input.size(); // 用于指明输入图像的通道数

	Mat mBias_1;
	vector<vector<Mat>> mKernels_1;
	vector<Mat> imgFilter2d_1(kernels_1);

	// 图像路径你自己写, 包括.扩展名称
	H5File file("F:\\Tmp\\mnist.h5", H5F_ACC_RDONLY);

	LoadWeight(file, "conv_1", mKernels_1, mBias_1);

	for (int i = 0; i < kernels_1; i++)
	{
		for (int j = 0; j < ker_channels_1; j++)
		{
			Mat img_tmp;

			filter2D(img_input[j], img_tmp, mKernels_1[i][j].depth(), mKernels_1[i][j]);

			if (0 == j)
			{
				imgFilter2d_1[i] = img_tmp;
			}
			else
			{
				imgFilter2d_1[i] += img_tmp;
			}
		}

		imgFilter2d_1[i] += mBias_1.at<float>(i, 0);

		threshold(imgFilter2d_1[i], imgFilter2d_1[i], 0, 1, CV_THRESH_TOZERO);
	}

	const int stride_2 = 2;
	vector<Mat> imgPool2d_2(kernels_1);

	for (int i = 0; i < kernels_1; i++)
	{
		imgPool2d_2[i].create((imgFilter2d_1[i].rows - 2) >> 1, (imgFilter2d_1[i].cols - 2) >> 1, CV_32FC1);

		// 起始从 1 开始, 结束为 ImgFilter2d_1[i].rows - 1
		for (int r = 1; r < imgFilter2d_1[i].rows - 1; r += stride_2)
		{
			for (int c = 1; c < imgFilter2d_1[i].cols - 1; c += stride_2)
			{
				const Mat m = imgFilter2d_1[i](cv::Rect(c, r, 2, 2));

				double max_val = 0;
				minMaxLoc(m, nullptr, &max_val);

				// 因为输出图像尺寸减小一半, 所以坐标要除以2
				imgPool2d_2[i].at<float>((r - 1) >> 1, (c - 1) >> 1) = (float)max_val;
			}
		}
	}

	const int kernels_3 = 64;
	const int ker_channels_3 = kernels_1;

	Mat mBias_3;
	vector<vector<Mat>> mKernels_3;
	vector<Mat> imgFilter2d_3(kernels_3);

	LoadWeight(file, "conv_2", mKernels_3, mBias_3);

	for (int i = 0; i < kernels_3; i++)					// 第二重循环, Filter的数量
	{
		for (int j = 0; j < ker_channels_3; j++)		// 第二重循环, Filter的通道, 也就是输入图像的通道数
		{
			Mat ImgTmp;

			// ImgPool2d_2[j] 为上一层的输出
			filter2D(imgPool2d_2[j], ImgTmp, mKernels_3[i][j].depth(), mKernels_3[i][j]);

			if (0 == j)
			{
				imgFilter2d_3[i] = ImgTmp;
			}
			else
			{
				imgFilter2d_3[i] += ImgTmp;
			}
		}

		imgFilter2d_3[i] += mBias_3.at<float>(i, 0);

		threshold(imgFilter2d_3[i], imgFilter2d_3[i], 0, 1, CV_THRESH_TOZERO);
	}

	const int stride_4 = 1;
	vector<Mat> imgPool2d_4(kernels_3);

	for (int i = 0; i < kernels_3; i++)									  // 第一重循环, 通道
	{
		imgPool2d_4[i].create(imgFilter2d_3[i].rows - 3, imgFilter2d_3[i].cols - 3, CV_32FC1);

		for (int r = 1; r < imgFilter2d_3[i].rows - 2; r += stride_4)     // 第二重循环, 图像每一行
		{
			for (int c = 1; c < imgFilter2d_3[i].cols - 2; c += stride_4) // 第三重循环, 图像每一列
			{
				const Mat m = imgFilter2d_3[i](cv::Rect(c, r, 2, 2));

				double max_val = 0;
				minMaxLoc(m, nullptr, &max_val);

				// 因为输出图像尺寸只减小了 1, 就不除以 2 了
				imgPool2d_4[i].at<float>(r - 1, c - 1) = (float)max_val;
			}
		}
	}

	// 这里使用方法 1 展开
	// 如果用方法 2, 千万不要在 output_6 计算完成之前删除 delete []flat_data,
	// 要不然你会找不到问题在哪里

	Mat multi_channel;
	merge(imgPool2d_4, multi_channel);

	Mat mFlat(multi_channel.rows * multi_channel.cols * multi_channel.channels(), 1,
		CV_32FC1, multi_channel.data);

	Mat mWeight_6;
	Mat mBias_6;

	LoadWeight(file, "dense_1", mWeight_6, mBias_6);

	Mat output_6 = mWeight_6 * mFlat + mBias_6;

	threshold(output_6, output_6, 0, 1, CV_THRESH_TOZERO);

	Mat mWeight_7;
	Mat mBias_7;

	LoadWeight(file, "output", mWeight_7, mBias_7);

	file.close();

	Mat output_7 = mWeight_7 * output_6 + mBias_7;

	float sum = 0;
	vector<float> tmp(output_7.rows);

	for (int i = 0; i < output_7.rows; i++)
	{
		tmp[i] = exp(output_7.at<float>(i, 0));
		sum += tmp[i];
	}

	for (int i = 0; i < output_7.rows; i++)
	{
		output_7.at<float>(i, 0) = tmp[i] / sum;

		cout << "output: " << output_7.at<float>(i, 0) << endl;
	}

	Point2i max_pt;
	minMaxLoc(output_7, nullptr, nullptr, nullptr, &max_pt);

	cout << "预测结果是: " << max_pt.y << endl << endl;

	system("pause");

	return 0;
}


/*================================================================
功能: 加载神经网络卷积Filter 和 bias 参数
传入参数:	
	1. file: 已经打开的 H5File 变量
	2. layer_name: 神经网络层名称
	3. kernel: 卷积 Filter
	4. bias: 偏置
返回值: bool 类型, true: 成功, flase: 失败
================================================================*/
bool LoadWeight(const H5File & file, const char * layer_name, vector<vector<Mat>> & kernel, Mat & bias)
{
	H5std_string path = "/";

	Group rg(file.getObjId(path));

	bool find = false;

	hsize_t objs = rg.getNumObjs();

	for (hsize_t i = 0; i < objs; i++)
	{
		const H5std_string name = rg.getObjnameByIdx(i);

		if (layer_name == name)
		{
			find = true;

			Group sub_group(rg.getObjId(name));

			path = name + "/";
			path.append(sub_group.getObjnameByIdx(0));

			sub_group.close();

			break;
		}
	}

	if (!find)
	{
		rg.close();
		
		return false;
	}

	DataSet ds_kernel(rg.getObjId(path + "/kernel:0"));

	DataSpace dsp = ds_kernel.getSpace();
	H5::DataType dt = ds_kernel.getDataType();

	int rank = dsp.getSimpleExtentNdims();

	hsize_t *dims = new hsize_t[rank];
	int ndims = dsp.getSimpleExtentDims(dims);

	hsize_t data_Len = ds_kernel.getInMemDataSize() / sizeof(float);
	float *buf = new float[data_Len];

	ds_kernel.read(buf, dt);

	kernel.resize(dims[rank - 1]);

	for (int i = 0; i < dims[rank - 1]; i++)	// 第一重: Filter 序号
	{
		// Filter 的通道数存储在 dims[rank - 2] 中
		kernel[i].resize(dims[rank - 2]);

		for (int j = 0; j < dims[rank - 2]; j++)	// 第二重: channel 序号
		{
			// dims[0], dims[1] 中分别是 Filter 的行数和列数
			kernel[i][j].create(dims[0], dims[1], CV_32FC1);

			for (int r = 0; r < kernel[i][j].rows; r++)	// 第三重: Filter 行序号
			{
				for (int c = 0; c < kernel[i][j].cols; c++)	// 第四重: Filter 列序号
				{
					const int k =
						c * dims[rank - 1] * dims[rank - 2] +
						r * dims[rank - 1] * dims[rank - 2] * kernel[i][j].cols +
						j * dims[rank - 1] +
						i;

					kernel[i][j].at<float>(r, c) = buf[k];
				}
			}
		}
	}

	dt.close();
	dsp.close();
	ds_kernel.close();

	delete []dims;
	dims = nullptr;

	delete []buf;
	buf = nullptr;

	//-------------------------------------------
	DataSet ds_bias(rg.getObjId(path + "/bias:0"));

	dsp = ds_bias.getSpace();
	dt = ds_bias.getDataType();

	rank = dsp.getSimpleExtentNdims();

	dims = new hsize_t[rank];
	ndims = dsp.getSimpleExtentDims(dims);

	if (1 == rank)
	{
		// 这里这样就不用转置了
		bias.create((int)dims[0], 1, CV_32FC1);
	}
	else
	{
		rg.close();
		
		return false;
	}

	ds_bias.read(bias.data, dt);

	dt.close();
	dsp.close();
	ds_bias.close();

	rg.close();

	delete []dims;
	dims = nullptr;

	return true;
}

/*================================================================
功能: 加载神经网络全连接层 weight 和 bias 参数
传入参数:	
	1. file: 已经打开的 H5File 变量
	2. layer_name: 神经网络层名称
	3. weight: 全连接层权重
	4. bias: 偏置
返回值: bool 类型, true: 成功, flase: 失败
================================================================*/

bool LoadWeight(const H5File & file, const char * layer_name, Mat & weight, Mat & bias)
{
	H5std_string path = "/";

	Group rg(file.getObjId(path));

	bool find = false;

	hsize_t objs = rg.getNumObjs();

	for (hsize_t i = 0; i < objs; i++)
	{
		const H5std_string name = rg.getObjnameByIdx(i);

		if (layer_name == name)
		{
			find = true;

			Group sub_group(rg.getObjId(name));

			path = name + "/";
			path.append(sub_group.getObjnameByIdx(0));

			sub_group.close();

			break;
		}
	}

	if (!find)
	{
		rg.close();
		
		return false;
	}

	DataSet ds_kernel(rg.getObjId(path + "/kernel:0"));

	DataSpace dsp = ds_kernel.getSpace();
	H5::DataType dt = ds_kernel.getDataType();

	int rank = dsp.getSimpleExtentNdims();

	hsize_t *dims = new hsize_t[rank];
	int ndims = dsp.getSimpleExtentDims(dims);

	if (1 == rank)
	{
		weight.create(1, (int)dims[0], CV_32FC1);
	}
	else if (2 == rank)
	{
		weight.create((int)dims[0], (int)dims[1], CV_32FC1);
	}
	else
	{
		rg.close();
		
		return false;
	}

	ds_kernel.read(weight.data, dt);

	transpose(weight, weight);

	dt.close();
	dsp.close();
	ds_kernel.close();

	delete []dims;
	dims = nullptr;

	//-------------------------------------------
	DataSet ds_bias(rg.getObjId(path + "/bias:0"));

	dsp = ds_bias.getSpace();
	dt = ds_bias.getDataType();

	rank = dsp.getSimpleExtentNdims();

	dims = new hsize_t[rank];
	ndims = dsp.getSimpleExtentDims(dims);

	if (1 == rank)
	{
		// 这里这样就不用转置了
		bias.create((int)dims[0], 1, CV_32FC1);
	}
	else
	{
		rg.close();
		
		return false;
	}

	ds_bias.read(bias.data, dt);

	dt.close();
	dsp.close();
	ds_bias.close();

	rg.close();

	delete []dims;
	dims = nullptr;

	return true;
}

七. 测试

测试图像为数字 0, 如下
zero
测试结果
test_result
你看到的值和我的肯定不一样, 因为你训练的参数和我的不一样, 但是也应该把图像识别为 0

八. 代码下载

完整的代码可下载 VS2015 x64 代码示例

本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

C++ 和 OpenCV 实现卷积神经网络并加载 Keras 训练好的参数进行预测 的相关文章

  • 没有强命名的代码签名是否会让您的应用程序容易被滥用?

    尝试了解authenticode代码签名和强命名 我是否正确地认为 如果我对引用一些 dll 非强命名 的 exe 进行代码签名 恶意用户就可以替换我的 DLL 并以看似由我签名但正在运行的方式分发应用程序他们的代码 假设这是真的 那么您似
  • 以文化中立的方式将字符串拆分为单词

    我提出了下面的方法 旨在将可变长度的文本拆分为单词数组 以进行进一步的全文索引处理 删除停止词 然后进行词干分析 结果似乎不错 但我想听听关于这种实现对于不同语言的文本的可靠性的意见 您会建议使用正则表达式来代替吗 请注意 我选择不使用 S
  • Web 客户端和 Expect100Continue

    使用 WebClient C NET 时设置 Expect100Continue 的最佳方法是什么 我有下面的代码 我仍然在标题中看到 100 continue 愚蠢的 apache 仍然抱怨 505 错误 string url http
  • 为什么两个不同的 Base64 字符串的转换会返回相等的字节数组?

    我想知道为什么从 base64 字符串转换会为不同的字符串返回相同的字节数组 const string s1 dg const string s2 dq byte a1 Convert FromBase64String s1 byte a2
  • 秒表有最长运行时间吗?

    多久可以Stopwatch在 NET 中运行 如果达到该限制 它会回绕到负数还是从 0 重新开始 Stopwatch Elapsed返回一个TimeSpan From MSDN https learn microsoft com en us
  • 嵌套接口:将 IDictionary> 转换为 IDictionary>?

    我认为投射一个相当简单IDictionary
  • 如何使用 ICU 解析汉字数字字符?

    我正在编写一个使用 ICU 来解析由汉字数字字符组成的 Unicode 字符串的函数 并希望返回该字符串的整数值 五 gt 5 三十一 gt 31 五千九百七十二 gt 5972 我将区域设置设置为 Locale getJapan 并使用
  • HTTPWebResponse 响应字符串被截断

    应用程序正在与 REST 服务通信 Fiddler 显示作为 Apps 响应传入的完整良好 XML 响应 该应用程序位于法属波利尼西亚 在新西兰也有一个相同的副本 因此主要嫌疑人似乎在编码 但我们已经检查过 但空手而归 查看流读取器的输出字
  • OleDbDataAdapter 未填充所有行

    嘿 我正在使用 DataAdapter 读取 Excel 文件并用该数据填充数据表 这是我的查询和连接字符串 private string Query SELECT FROM Sheet1 private string ConnectStr
  • C#中如何移动PictureBox?

    我已经使用此代码来移动图片框pictureBox MouseMove event pictureBox Location new System Drawing Point e Location 但是当我尝试执行时 图片框闪烁并且无法识别确切
  • 创建链表而不将节点声明为指针

    我已经在谷歌和一些教科书上搜索了很长一段时间 我似乎无法理解为什么在构建链表时 节点需要是指针 例如 如果我有一个节点定义为 typedef struct Node int value struct Node next Node 为什么为了
  • 将多个表映射到实体框架中的单个实体类

    我正在开发一个旧数据库 该数据库有 2 个具有 1 1 关系的表 目前 我为每个定义的表定义了一种类型 1Test 1Result 我想将这些特定的表合并到一个类中 当前的类型如下所示 public class Result public
  • 显示UnityWebRequest的进度

    我正在尝试使用下载 assetbundle统一网络请求 https docs unity3d com ScriptReference Networking UnityWebRequest GetAssetBundle html并显示进度 根
  • 使用 Bearer Token 访问 IdentityServer4 上受保护的 API

    我试图寻找此问题的解决方案 但尚未找到正确的搜索文本 我的问题是 如何配置我的 IdentityServer 以便它也可以接受 授权带有 BearerTokens 的 Api 请求 我已经配置并运行了 IdentityServer4 我还在
  • 什么时候虚拟继承是一个好的设计? [复制]

    这个问题在这里已经有答案了 EDIT3 请务必在回答之前清楚地了解我要问的内容 有 EDIT2 和很多评论 有 或曾经 有很多答案清楚地表明了对问题的误解 我知道这也是我的错 对此感到抱歉 嗨 我查看了有关虚拟继承的问题 class B p
  • 通过指向其基址的指针删除 POD 对象是否安全?

    事实上 我正在考虑那些微不足道的可破坏物体 而不仅仅是POD http en wikipedia org wiki Plain old data structure 我不确定 POD 是否可以有基类 当我读到这个解释时is triviall
  • 基于 OpenCV 边缘的物体检测 C++

    我有一个应用程序 我必须检测场景中某些项目的存在 这些项目可以旋转并稍微缩放 更大或更小 我尝试过使用关键点检测器 但它们不够快且不够准确 因此 我决定首先使用 Canny 或更快的边缘检测算法 检测模板和搜索区域中的边缘 然后匹配边缘以查
  • 混合 ExecutionContext.SuppressFlow 和任务时 AsyncLocal.Value 出现意外值

    在应用程序中 由于 AsyncLocal 的错误 意外值 我遇到了奇怪的行为 尽管我抑制了执行上下文的流程 但 AsyncLocal Value 属性有时不会在新生成的任务的执行范围内重置 下面我创建了一个最小的可重现示例来演示该问题 pr
  • 如何防止用户控件表单在 C# 中处理键盘输入(箭头键)

    我的用户控件包含其他可以选择的控件 我想实现使用箭头键导航子控件的方法 问题是家长控制拦截箭头键并使用它来滚动其视图什么是我想避免的事情 我想自己解决控制内容的导航问题 我如何控制由箭头键引起的标准行为 提前致谢 MTH 这通常是通过重写
  • 使用.NET技术录制屏幕视频[关闭]

    Closed 这个问题正在寻求书籍 工具 软件库等的推荐 不满足堆栈溢出指南 help closed questions 目前不接受答案 有没有一种方法可以使用 NET 技术来录制屏幕 无论是桌面还是窗口 我的目标是免费的 我喜欢小型 低

随机推荐