Google Earth Engine（GEE）深度学习入门教程-Python数据读入篇

齐硕君

已于 2024-05-23 14:24:51 修改

阅读量954

点赞数 6

分类专栏： GEE 深度学习文章标签：深度学习 python 人工智能

于 2024-05-23 14:23:48 首次发布

本文链接：https://blog.csdn.net/weixin_41099712/article/details/139146766

版权

GEE 同时被 2 个专栏收录

6 篇文章 1 订阅

订阅专栏

深度学习

6 篇文章 2 订阅

订阅专栏

Python数据读入篇

前置条件：

GEE预处理影像导出保存为tfrecord的数据包，并下载到本地
tensorflow的深度学习环境

本篇文章的目的主要是把Tfrecord格式的数据加载为tf可使用的数据集格式

设定超参数

首先需要设定导出时的波段名称和数据格式，用于解析数据。

变量名：

KERNEL_SIZE ：单个patch的大小
FEATURES：波段的名称
dtype=tf.float32 ：数据格式

# 所有波段的名称
BANDS = ["B2","B3","B4","B5","B6","B7","B8","B8A","B11","B12"]
TIMES = 12 
ALLBANDS = [ "%d_"%(i)+s for i in range(TIMES) for s in BANDS ]
print('所有的波段名称为:',ALLBANDS)
RESPONSE = 'soya'
FEATURES = ALLBANDS + [RESPONSE]
#输出：
#所有的波段名称为: ['0_B2', '0_B3', '0_B4', '0_B5', '0_B6', '0_B7', '0_B8', '0_B8A', '0_B11', '0_B12', '1_B2', '1_B3', '1_B4', '1_B5', '1_B6', '1_B7', '1_B8', '1_B8A', '1_B11', '1_B12', '2_B2', '2_B3', '2_B4', '2_B5', '2_B6', '2_B7', '2_B8', '2_B8A', '2_B11', '2_B12', '3_B2', '3_B3', '3_B4', '3_B5', '3_B6', '3_B7', '3_B8', '3_B8A', '3_B11', '3_B12', '4_B2', '4_B3', '4_B4', '4_B5', '4_B6', '4_B7', '4_B8', '4_B8A', '4_B11', '4_B12', '5_B2', '5_B3', '5_B4', '5_B5', '5_B6', '5_B7', '5_B8', '5_B8A', '5_B11', '5_B12', '6_B2', '6_B3', '6_B4', '6_B5', '6_B6', '6_B7', '6_B8', '6_B8A', '6_B11', '6_B12', '7_B2', '7_B3', '7_B4', '7_B5', '7_B6', '7_B7', '7_B8', '7_B8A', '7_B11', '7_B12', '8_B2', '8_B3', '8_B4', '8_B5', '8_B6', '8_B7', '8_B8', '8_B8A', '8_B11', '8_B12', '9_B2', '9_B3', '9_B4', '9_B5', '9_B6', '9_B7', '9_B8', '9_B8A', '9_B11', '9_B12', '10_B2', '10_B3', '10_B4', '10_B5', '10_B6', '10_B7', '10_B8', '10_B8A', '10_B11', '10_B12', '11_B2', '11_B3', '11_B4', '11_B5', '11_B6', '11_B7', '11_B8', '11_B8A', '11_B11', '11_B12','soya']

# 单个patch的大小和数据类型
KERNEL_SIZE = 64
KERNEL_SHAPE = [KERNEL_SIZE, KERNEL_SIZE]
COLUMNS = [
  tf.io.FixedLenFeature(shape=KERNEL_SHAPE, dtype=tf.float32) for k in FEATURES
]
FEATURES_DICT = dict(zip(FEATURES, COLUMNS))

# 单个数据包中样本的数量
BUFFER_SIZE = 160

解析数据

编写解析数据包的函数

parse_tfrecord( )和get_dataset( )函数是基本固定的

to_HWTCtuple(inputs)函数可根据自己的数据格式进行修改。数据的归一化和标准化、拓展特征等处理操作都可以放在此函数里。本文此处仅仅将折叠120维的数据展开为12(时相)×10(波段)，数据集单个样本（HWTC格式）的大小为此时64×64×12×10。

def parse_tfrecord(example_proto):
  """The parsing function.
  Read a serialized example into the structure defined by FEATURES_DICT.
  Args:
    example_proto: a serialized Example.
  Returns:
    A dictionary of tensors, keyed by feature name.
  """
  return tf.io.parse_single_example(example_proto, FEATURES_DICT)
def to_HWTCtuple(inputs):
    """Function to convert a dictionary of tensors to a tuple of (inputs, outputs).
    Turn the tensors returned by parse_tfrecord into a stack in HWCT shape.
    Args:
    inputs: A dictionary of tensors, keyed by feature name.
    Returns:
    A tuple of (inputs, outputs).
    """
    # 读取出多时相多波段影像
    inputsList = [inputs.get(key) for key in ALLBANDS]
    label = inputs.get(RESPONSE)
    stacked = tf.stack(inputsList, axis=0)
    # 从通道维展开出时间维,TIMES在前，len(BANDS)在后，不可以反
    #stacked = tf.reshape(stacked, [12,10,KERNEL_SIZE, KERNEL_SIZE])
    stacked = tf.reshape(stacked, [TIMES,len(BANDS),KERNEL_SIZE, KERNEL_SIZE])
    # Convert from TCHW to HWTC
    stacked = tf.transpose(stacked, [2,3,0,1])
    
    return stacked, label
def get_dataset(namelist):
    """Function to read, parse and format to tuple a set of input tfrecord files.
    Get all the files matching the pattern, parse and convert to tuple.
    Args:
    namelist: A file list in a Cloud Storage bucket.
    Returns:
    A tf.data.Dataset
    """
    dataset = tf.data.TFRecordDataset(namelist, compression_type='GZIP')
    dataset = dataset.map(parse_tfrecord, num_parallel_calls=5)
    dataset = dataset.map(to_HWTCtuple, num_parallel_calls=5)
    return dataset

编写训练集和验证机的加载函数，根据自己的需求加载指定路径的数据包。BATCH_SIZE可根据自己GPU显存设置，此处为8

trainPath = 'data/guoyang/training_patches_g%d.tfrecord.gz'
trainIndex = [14,20,19,25,13,31,24,15,18,7,30,26,21,29,16,8,28,27,12,11,17,33,4,32,9,10,22]
training_patches_size = 27
BATCH_SIZE = 8
def get_training_dataset():
	"""Get the preprocessed training dataset
  Returns: 
    A tf.data.Dataset of training data.
  """
	namelist = [trainPath%(i) for i in trainIndex]
	dataset = get_dataset(namelist)
	dataset = dataset.shuffle(BUFFER_SIZE).batch(BATCH_SIZE,drop_remainder=True).repeat()
	return dataset

evalPath = 'data/guoyang/eval_patches_g%d.tfrecord.gz'
evalIndex = [9,12,13,3,5,1,6,4,2]
eval_patches_size = 9
EVAL_BATCH_SIZE = 8
def get_eval_dataset():
	"""Get the preprocessed evaluation dataset
  Returns: 
    A tf.data.Dataset of evaluation data.
  """
	namelist = [evalPath%(i) for i in evalIndex]
	dataset = get_dataset(namelist)
	dataset = dataset.batch(EVAL_BATCH_SIZE,drop_remainder=True).repeat()
	return dataset

evaluation = get_eval_dataset()

training = get_training_dataset()

读取单个样本，测试数据集是否加载成功。

data,label = iter(training.take(1)).next()
i = 3
print(data.shape)
#(8, 64, 64, 12, 10)

可视化

首先编写并测试显示所需要的函数。

#测试显示函数
import tensorflow as tf
import matplotlib.pyplot as plt
import numpy as np
from PIL import Image
from matplotlib.colors import ListedColormap

def show_image(tensor,index ,savename =None):
    # 从HWTC到THWC
    tensor = np.transpose(tensor, (2, 0,1,3))
    # 获取 tensor 中每个通道的数据
    r_channel = tensor[:, :, :, index[0]-2]  # 红色通道
    g_channel = tensor[:, :, :, index[1]-2]  # 绿色通道
    b_channel = tensor[:, :, :, index[2]-2]  # 蓝色通道
    # 将通道的数据拼接成一个 RGB 图像
    image = np.stack([r_channel, g_channel, b_channel], axis=3)
    
    # 定义一个变量来保存图片数量
    num_images = image.shape[0]
    # 定义每行显示的图片数量
    num_images_per_row = 4
    # 计算行数和列数
    num_rows = int(np.ceil(num_images / num_images_per_row))
    num_columns = num_images_per_row
    
    # 定义画布的大小
    fig, axes = plt.subplots(num_rows, num_columns, figsize=(15, 15))
    
    # 展示每张图片
    for i, ax in enumerate(axes.flat):
        # 如果图片数量不足，则跳过
        if i >= num_images:
            break
        # 分别计算每个波段的%2的像素值
        img = image[i]
        vmin_b1, vmax_b1 = np.percentile(img[:, :, 0], (2, 98))
        vmin_b2, vmax_b2 = np.percentile(img[:, :, 1], (2, 98))
        vmin_b3, vmax_b3 = np.percentile(img[:, :, 2], (2, 98))
        # 线性拉伸每个波段并合并为RGB彩色图像
        image_r = np.interp(img[:, :, 0], (vmin_b1, vmax_b1), (0, 1))
        image_g = np.interp(img[:, :, 1], (vmin_b2, vmax_b2), (0, 1))
        image_b = np.interp(img[:, :, 2], (vmin_b3, vmax_b3), (0, 1))
        image_rgb = np.dstack((image_r, image_g, image_b))
        # 显示图片
        ax.imshow(image_rgb)
        if savename:
            Image.fromarray(np.uint8(image_rgb*255)).save(f"img/{savename}_{i}.png")  #.resize((512,512),resample = Image.NEAREST)
        ax.set_title('TimeSeries Image%d'%(i))
        # 关闭刻度线
        ax.axis('off')
    # 展示所有图片
    plt.show()
def show_label(tensor , savename =None):
    # 创建一个浅绿色的 colormap
    cmap = ListedColormap(['white', 'palegreen'])

    # 显示二值化图像
    fig, ax = plt.subplots()
    ax.imshow(tensor, cmap=cmap)
    if savename:
        Image.fromarray(np.uint8(tensor*255)).save(f"img/{savename}.png")   #.resize((512,512),resample = Image.NEAREST)
        # 隐藏坐标轴
    ax.axis('off')
    ax.set_title('Label Image')
    plt.show()

# 随机生成一个大小为 [2, 128, 128, 10] 的张量
tf.random.set_seed(0)
tensor = tf.random.normal([128, 128,4,10])

# 调用 show_image 函数进行展示
show_image(tensor.numpy(),[2,3,4])

# 随机生成一个 8x8 的二值化图像
image = np.random.randint(2, size=(8, 8))
show_label(image)

编写区块拼接函数，测试区块通常是通过滑窗的方式裁切出多个样本，因此在导出最终结果的时候需要将整个区块拼接起来并保存。

此函数部分需要根据实际情况进行具体调整，因为GEE导出时边缘似乎是不太固定的（未解决的问题），所以导出的时候需要增加缓冲区，在此处再将缓冲区裁剪掉。num_rows 和 num_cols代表区块内patch的行列数，可以从GEE导出的json格式文件中获取，num_rows和num_cols不一致是因为GEE导出坐标系为默认EPSG：4326（地理坐标系），作者后来更改为EPSG：32650（投影坐标系），num_rows和num_cols就相同了（区块为矩形）。

def mosicImage(data,core_size = 32,num_rows = 15,num_cols = 18):
    buffer = int(core_size/2)
    # 创建大图的空数组
    big_image = np.zeros(((num_rows+1) * core_size, (num_cols+1) * core_size, 12, 10))
    # 遍历原始数据
    index = 0
    for i in range(num_rows):
        for j in range(num_cols):
            # 计算当前核心区域的索引范围
            start_row = i * core_size  
            end_row = (i + 1) * core_size 
            start_col = j * core_size
            end_col = (j + 1) * core_size
            
            # 提取当前核心区域的数据
            core_data = data[index,buffer:buffer+core_size, buffer:buffer+core_size,:,:]

            # 将核心区域的数据放入大图的对应位置
            big_image[start_row:end_row, start_col:end_col,:,:] = core_data
            
            index = index +1 
    index = index - num_cols
    #拼接完（480*608）还要拼接最后一列
    #拼接最后一边缘行    拼接完（16*608）
    for j in range(num_cols):
        # 计算当前核心区域的索引范围
        start_row = num_rows * core_size 
        end_row = (num_rows) * core_size +  buffer
        start_col = j * core_size
        end_col = (j + 1) * core_size

        # 提取当前核心区域的数据
        core_data = data[index,buffer+core_size:, buffer:buffer+core_size,:,:]

        # 将核心区域的数据放入大图的对应位置
        big_image[start_row:end_row, start_col:end_col,:,:] = core_data

        index = index +1 
    index = 1
    for j in range(num_rows):
        # 计算当前核心区域的索引范围
        start_col = num_cols * core_size 
        end_col = (num_cols) * core_size +  buffer
        start_row = j * core_size
        end_row = (j + 1) * core_size

        # 提取当前核心区域的数据
        core_data = data[index*num_cols-1,buffer:buffer+core_size,buffer+core_size:,:,:]

        # 将核心区域的数据放入大图的对应位置
        big_image[start_row:end_row, start_col:end_col,:,:] = core_data

        index = index +1 
    #拼接完（496*608）
    #根据实际情况调整
    big_image = big_image[8:8+500,8:8+500]
    big_image = tf.constant(big_image)
    return big_image
def mosicLabel(data,core_size = 32,num_rows = 15,num_cols = 18):
    buffer = int(core_size/2)
    # 创建大图的空数组
    big_image = np.zeros(((num_rows+1) * core_size, (num_cols+1) * core_size))
    # 遍历原始数据
    index = 0
    for i in range(num_rows):
        for j in range(num_cols):
            # 计算当前核心区域的索引范围
            start_row = i * core_size  
            end_row = (i + 1) * core_size 
            start_col = j * core_size
            end_col = (j + 1) * core_size

            # 提取当前核心区域的数据
            core_data = data[index,buffer:buffer+core_size, buffer:buffer+core_size]

            # 将核心区域的数据放入大图的对应位置
            big_image[start_row:end_row, start_col:end_col] = core_data
            
            index = index +1 
    index = index - num_cols
    #拼接完（480*608）还要拼接最后一列
    #拼接最后一边缘行    拼接完（16*608）
    for j in range(num_cols):
        # 计算当前核心区域的索引范围
        start_row = num_rows * core_size 
        end_row = (num_rows) * core_size +  buffer
        start_col = j * core_size
        end_col = (j + 1) * core_size

        # 提取当前核心区域的数据
        core_data = data[index,buffer+core_size:, buffer:buffer+core_size]

        # 将核心区域的数据放入大图的对应位置
        big_image[start_row:end_row, start_col:end_col] = core_data

        index = index +1 
    #拼接完（496*608）
    index = 1
    for j in range(num_rows):
        # 计算当前核心区域的索引范围
        start_col = num_cols * core_size 
        end_col = (num_cols) * core_size +  buffer
        start_row = j * core_size
        end_row = (j + 1) * core_size

        # 提取当前核心区域的数据
        core_data = data[index*num_cols-1,buffer:buffer+core_size,buffer+core_size:]

        # 将核心区域的数据放入大图的对应位置
        big_image[start_row:end_row, start_col:end_col] = core_data

        index = index +1 
    #根据实际情况调整
    big_image = big_image[8:8+500,8:8+500]
    big_image = tf.constant(big_image)
    return big_image

设置数据包路径，拼接、显示并保存区块。

directory = 'data/yangfang/'  # 替换为你要列出文件的目录路径
files = ['2022training_patches_g2.tfrecord.gz']
# 遍历多个样方
for file in files:
    data = []
    labels = []
    dataset = get_dataset([directory+file])
    for example in dataset:
        data.append(example[0])
        labels.append(example[1])   
    data = tf.stack(data)
    labels = tf.stack(labels)
    print(data.shape)     

    # 拼接数据 
    big_image = mosicImage(data)
    big_label = mosicLabel(labels)
       
    #显示 并保存
    savename =  file.split('.')[0]
    savePath = f"yangfang/{savename}/"
    if not os.path.exists('img/'+savePath):
        os.makedirs('img/'+savePath)
        print(f'{savePath}文件夹已创建!')
        
    show_image(big_image.numpy(),[8,4,3],savePath+f"{savename}")
    show_label(big_label.numpy(),savePath+f"{savename}-label")

文件夹内保存的文件：

其他

数据增广

数据增广是一种有用的技术，可以增加训练数据的数量和多样性。首先将原始瓷砖图像顺时针随机旋转 90°、180° 或 270° ;然后，将旋转的图像乘以 0.9 和 1.1 之间的随机值（波段间相同，时相间不同），以增加样本数据时间序列曲线的波动。此外，每个原始平铺图像被垂直或水平翻转，随后将翻转的图像乘以0.9和1.1之间的随机数。训练数据集增广为原来的4倍。

trainPath = 'data/guoyang/training_patches_g%d.tfrecord.gz'
trainIndex = [14,20,19,25,13,31,24,15,18,7,30,26,21,29,16,8,28,27,12,11,17,33,4,32,9,10,22]
training_patches_size = 27
BATCH_SIZE = 8
def to_HWTCtuple(inputs):
    # 读取出多时相多波段影像
    inputsList = [inputs.get(key) for key in ALLBANDS]
    label = inputs.get(RESPONSE)
    stacked = tf.stack(inputsList, axis=0)
    # 从通道维展开出时间维,TIMES在前，len(BANDS)在后，不可以反
    #stacked = tf.reshape(stacked, [12,10,KERNEL_SIZE, KERNEL_SIZE])
    stacked = tf.reshape(stacked, [TIMES,len(BANDS),KERNEL_SIZE, KERNEL_SIZE])
    # Convert from TCHW to HWTC
    stacked = tf.transpose(stacked, [2,3,0,1])   
    return stacked, label
def generate_augment():
    #生成增强函数
    #method = 1      2      3       4       5 
    #对应    旋转90°、 180°、 270°、水平翻转、竖直翻转 
    #用法：dataset = dataset.map(generate_augment(), num_parallel_calls=5)

    def augment_data(time_series,label):
        label = tf.expand_dims(label, -1)        #(H,W)-->(H,W,1)
        # 旋转 90°、180°、 270°、水平翻转、数值翻转 
        time_series = tf.transpose(time_series, perm=[2, 0, 1, 3])   #(T ,H ,W,C)
        
        # 生成1到5之间的整数 [1,6)
        method = np.random.randint(1,6)
        if method <=3:
            time_series = tf.image.rot90(time_series, k=method)
            label = tf.image.rot90(label, k=method)
        elif method == 4:
            time_series = tf.image.flip_left_right(time_series)
            label = tf.image.flip_left_right(label)
        elif method == 5:
            time_series = tf.image.flip_up_down(time_series)
            label = tf.image.flip_up_down(label)
        else:
            raise ValueError("没有此类数据增广方法")

        time_series = tf.transpose(time_series, perm=[1, 2, 0, 3])   #(H ,W, T,C) 

        # 随机选择时间序列缩放因子
        scale_factor = tf.random.uniform(shape=(12,), minval=0.9, maxval=1.1,seed = 0)
        scale_factor = tf.expand_dims(scale_factor, 0)     #(1,12)
        scale_factor = tf.expand_dims(scale_factor, 0)     #(1,1,12)
        scale_factor = tf.expand_dims(scale_factor, -1)     #(1,1,12,1)
        # 时间序列缩放
        time_series *= scale_factor
        
        label = tf.squeeze(label)
        return time_series,label
    return augment_data
def get_augment_training_dataset():
	namelist = [trainPath%(i) for i in trainIndex]
	dataset = get_dataset(namelist)
	augment_dataset1 = dataset.map(generate_augment(), num_parallel_calls=5)
	augment_dataset2 = dataset.map(generate_augment(), num_parallel_calls=5)
	augment_dataset3 = dataset.map(generate_augment(), num_parallel_calls=5)   
	dataset = dataset.concatenate(augment_dataset1).concatenate(augment_dataset2).concatenate(augment_dataset3)
	dataset = dataset.shuffle(BUFFER_SIZE).batch(BATCH_SIZE,drop_remainder=True).repeat()
	return dataset

齐硕君

关注

6
点赞
踩
7

收藏

觉得还不错? 一键收藏
0
评论
Google Earth Engine（GEE）深度学习入门教程-Python数据读入篇

num_rows 和 num_cols代表区块内patch的行列数，可以从GEE导出的json格式文件中获取，num_rows和num_cols不一致是因为GEE导出坐标系为默认EPSG：4326（地理坐标系），作者后来更改为EPSG：32650（投影坐标系），num_rows和num_cols就相同了（区块为矩形）。然后，将旋转的图像乘以 0.9 和 1.1 之间的随机值（波段间相同，时相间不同），以增加样本数据时间序列曲线的波动。编写训练集和验证机的加载函数，根据自己的需求加载指定路径的数据包。
复制链接

扫一扫