中文亚洲精品无码_熟女乱子伦免费_人人超碰人人爱国产_亚洲熟妇女综合网

當(dāng)前位置: 首頁(yè) > news >正文

上海頻道做網(wǎng)站怎么樣希愛(ài)力雙效片副作用

上海頻道做網(wǎng)站怎么樣,希愛(ài)力雙效片副作用,如何建立營(yíng)銷型網(wǎng)站,茌平網(wǎng)站建設(shè)菜譜制作完整代碼見(jiàn) Github:https://github.com/ChaoQiezi/read_fy3d_tshs,由于代碼中注釋較為詳細(xì),因此博客中部分操作一筆帶過(guò)。 01 FY3D的HDF轉(zhuǎn)TIFF 1.1 數(shù)據(jù)集說(shuō)明 FY3D TSHS數(shù)據(jù)集是二級(jí)產(chǎn)品(TSHS即MWTS/MWHS 融合大氣溫濕度廓線/穩(wěn)定度指數(shù)/…

完整代碼見(jiàn) Github:https://github.com/ChaoQiezi/read_fy3d_tshs,由于代碼中注釋較為詳細(xì),因此博客中部分操作一筆帶過(guò)。

01 FY3D的HDF轉(zhuǎn)TIFF

1.1 數(shù)據(jù)集說(shuō)明

FY3D TSHS數(shù)據(jù)集是二級(jí)產(chǎn)品(TSHS即MWTS/MWHS 融合大氣溫濕度廓線/穩(wěn)定度指數(shù)/位勢(shì)高度產(chǎn)品);
文件格式為HDF5
空間分辨率為33KM(星下點(diǎn))
范圍為全球區(qū)域;(FY3D為極軌衛(wèi)星,因此對(duì)于得到的單個(gè)數(shù)據(jù)集并沒(méi)有完全覆蓋全球區(qū)域),類似下方(其余地區(qū)均未掃描到):

在這里插入圖片描述

提取數(shù)據(jù)涉及的多個(gè)數(shù)據(jù)集:

在這里插入圖片描述

位于DATA/MWHS_Ch_BT和DATA/MWTS_Ch_BT 路徑下的兩個(gè)數(shù)據(jù)集分別為MWTS 通道亮溫MWHS 通道亮溫,這是需要提取的數(shù)據(jù)集;
另外位于GEO/Latitude 和GEO/Longitude 路徑下的兩個(gè)數(shù)據(jù)集分別經(jīng)緯度數(shù)據(jù)集,用于確定像元的地理位置以便后續(xù)進(jìn)行重投影;

1.1.1 MWTS 通道亮溫和MWHS 通道亮溫?cái)?shù)據(jù)集的說(shuō)明

這是關(guān)于數(shù)據(jù)集的基本屬性說(shuō)明:

在這里插入圖片描述

在這里插入圖片描述

在產(chǎn)品說(shuō)明書中關(guān)于數(shù)據(jù)集的說(shuō)明如下:

在這里插入圖片描述

可以發(fā)現(xiàn),其共有三個(gè)維度,譬如MWHS_Ch_BT數(shù)據(jù)集的Shape為(15, 1212, 90),其表示極軌衛(wèi)星的傳感器15個(gè)通道(即波段數(shù))在每條掃描線總共 90 個(gè)觀測(cè)像元的MWHS亮溫值,此處掃描線共有1212個(gè)。對(duì)于MWTS數(shù)據(jù)集類似格式。(這意味著三個(gè)維度中并沒(méi)有空間上的關(guān)系)

1.1.2 經(jīng)緯度數(shù)據(jù)集的說(shuō)明

在這里插入圖片描述
在這里插入圖片描述

這是官方產(chǎn)品說(shuō)明對(duì)于經(jīng)緯度數(shù)據(jù)集的介紹:

在這里插入圖片描述

可以發(fā)現(xiàn),經(jīng)緯度數(shù)據(jù)集的Shape均為(1212, 90),這正好對(duì)應(yīng)前文提及的兩個(gè)數(shù)據(jù)集的所有像元(除去波段數(shù)),其表示每條1212條掃描線上的90個(gè)觀測(cè)像元的經(jīng)緯度。

1.2 讀取HDF5文件數(shù)據(jù)集

def read_h5(hdf_path, ds_path, scale=True):"""讀取指定數(shù)據(jù)集并進(jìn)行預(yù)處理:param hdf_path: 數(shù)據(jù)集所在HDF文件的絕對(duì)路徑:param ds_path: 數(shù)據(jù)集在HDF文件內(nèi)的路徑:return: 返回處理好的數(shù)據(jù)集"""with h5py.File(hdf_path) as f:# 讀取目標(biāo)數(shù)據(jù)集屬矩陣和相關(guān)屬性ds = f[ds_path]ds_values = np.float32(ds[:])  # 獲取數(shù)據(jù)集valid_range = ds.attrs['valid_range']  # 獲取有效范圍slope = ds.attrs['Slope'][0]  # 獲取斜率(類似scale_factor)intercept = ds.attrs['Intercept'][0]  # 獲取截距(類似add_offset)""""原始數(shù)據(jù)集可能存在縮放(可能是為存儲(chǔ)空間全部存為整數(shù)(需要通過(guò)斜率和截距再還原為原來(lái)的值,或者是需要進(jìn)行單位換算甚至物理量的計(jì)算例如最常見(jiàn)的DN值轉(zhuǎn)大氣層表觀反射率等(這多出現(xiàn)于一級(jí)產(chǎn)品的輻射定標(biāo), 但二級(jí)產(chǎn)品可能因?yàn)閱挝粨Q算等也可能出現(xiàn)));如果原始數(shù)據(jù)集不存在縮放, 那么Slope屬性和Intercept屬性分別為1和0;這里為了確保所有迭代的HDF文件的數(shù)據(jù)集均正常得到, 這里依舊進(jìn)行slope和intercept的讀取和計(jì)算(盡管可能冗余)"""# 目標(biāo)數(shù)據(jù)集的預(yù)處理(無(wú)效值, 范圍限定等)ds_values[(ds_values < valid_range[0]) | (ds_values > valid_range[1])] = np.nanif scale:ds_values = ds_values * slope + intercept  # 還原縮放"""Note: 這里之所以選擇是否進(jìn)行縮放, 原因?yàn)榻?jīng)緯度數(shù)據(jù)集中的slope為1, intercept也為1, 但是進(jìn)行縮放后超出地理范圍1°即出現(xiàn)了90.928對(duì)于緯度。其它類似, 因此認(rèn)為這里可能存在問(wèn)題如果進(jìn)行縮放, 所以對(duì)于經(jīng)緯度數(shù)據(jù)集這里不進(jìn)行縮放"""return ds_values

上述代碼用于讀取指定HDF5文件的指定數(shù)據(jù)集的數(shù)組/矩陣,scale參數(shù)用于是否對(duì)數(shù)據(jù)集進(jìn)行slope和intercept線性轉(zhuǎn)換。

1.3 重組和重投影

這一部分是整個(gè)數(shù)據(jù)處理的核心。

1.3.1 重組

def reform_ds(ds, lon, lat, reform_range=None):"""重組數(shù)組:param ds: 目標(biāo)數(shù)據(jù)集(三維):param lon: 對(duì)應(yīng)目標(biāo)數(shù)據(jù)集的經(jīng)度數(shù)據(jù)集():param lat: 對(duì)應(yīng)目標(biāo)數(shù)據(jù)集的緯度數(shù)據(jù)集(二維):param reform_range: 重組范圍, (lon_min, lat_max, lon_max, lat_min), 若無(wú)則使用全部數(shù)據(jù):return: 以元組形式依次返回: 重組好的目標(biāo)數(shù)據(jù)集, 經(jīng)度數(shù)據(jù)集, 緯度數(shù)據(jù)集(均為二維數(shù)組)"""# 裁選指定地理范圍的數(shù)據(jù)集if reform_range:lon_min, lat_max, lon_max, lat_min = reform_rangex, y = np.where((lon > lon_min) & (lon < lon_max) & (lat > lat_min) & (lat < lat_max))ds = ds[:, x, y]lon = lon[x, y]lat = lat[x, y]else:ds = ds.reshape(ds.shape[0], -1)lon = lon.flatten()lat = lat.flatten()# 無(wú)效值去除(去除地理位置為無(wú)效值的元素)valid_pos = ~np.isnan(lat) & ~np.isnan(lon)ds = ds[:, valid_pos]lon = lon[valid_pos]lat = lat[valid_pos]# 重組數(shù)組的初始化bands = []for band in ds:reform_ds_size = np.int32(np.sqrt(band.size))  # int向下取整band = band[:reform_ds_size ** 2].reshape(reform_ds_size, reform_ds_size)bands.append(band)else:lon = lon[:reform_ds_size ** 2].reshape(reform_ds_size, reform_ds_size)lat = lat[:reform_ds_size ** 2].reshape(reform_ds_size, reform_ds_size)bands = np.array(bands)return bands, lon, lat

這部分是對(duì)原始數(shù)據(jù)集(此處理中待重組數(shù)據(jù)集shape為(15, 1212, 90))進(jìn)行重組。
重組原理即基于經(jīng)緯度數(shù)據(jù)集,依據(jù)裁剪范圍將滿足地理范圍內(nèi)的所有有效像元一維化,然后重新reshape為二維數(shù)組,數(shù)組行列數(shù)均為原維度的平方根。

1.3.2 重投影

def data_glt(out_path, src_ds, src_x, src_y, out_res, zoom_scale=6, glt_range=None, windows_size=9):"""基于經(jīng)緯度數(shù)據(jù)集對(duì)目標(biāo)數(shù)據(jù)集進(jìn)行GLT校正/重投影(WGS84), 并輸出為TIFF文件:param out_path: 輸出tiff文件的路徑:param src_ds: 目標(biāo)數(shù)據(jù)集:param src_x: 對(duì)應(yīng)的橫軸坐標(biāo)系(對(duì)應(yīng)地理坐標(biāo)系的經(jīng)度數(shù)據(jù)集):param src_y: 對(duì)應(yīng)的縱軸坐標(biāo)系(對(duì)應(yīng)地理坐標(biāo)系的緯度數(shù)據(jù)集):param out_res: 輸出分辨率(單位: 度/°):param zoom_scale::return: None"""if glt_range:# lon_min, lat_max, lon_max, lat_min = -180.0, 90.0, 180.0, -90.0lon_min, lat_max, lon_max, lat_min = glt_rangeelse:lon_min, lat_max, lon_max, lat_min = np.nanmin(src_x), np.nanmax(src_y), \np.nanmax(src_x), np.nanmin(src_y)zoom_lon = zoom(src_x, (zoom_scale, zoom_scale), order=0)  # 0為最近鄰插值zoom_lat = zoom(src_y, (zoom_scale, zoom_scale), order=0)# # 確保插值結(jié)果正常# zoom_lon[(zoom_lon < -180) | (zoom_lon > 180)] = np.nan# zoom_lat[(zoom_lat < -90) | (zoom_lat > 90)] = np.nanglt_cols = np.ceil((lon_max - lon_min) / out_res).astype(int)glt_rows = np.ceil((lat_max - lat_min) / out_res).astype(int)deal_bands = []for src_ds_band in src_ds:glt_ds = np.full((glt_rows, glt_cols), np.nan)glt_lon = np.full((glt_rows, glt_cols), np.nan)glt_lat = np.full((glt_rows, glt_cols), np.nan)geo_x_ix, geo_y_ix = np.floor((zoom_lon - lon_min) / out_res).astype(int), \np.floor((lat_max - zoom_lat) / out_res).astype(int)glt_lon[geo_y_ix, geo_x_ix] = zoom_longlt_lat[geo_y_ix, geo_x_ix] = zoom_latglt_x_ix, glt_y_ix = np.floor((src_x - lon_min) / out_res).astype(int), \np.floor((lat_max - src_y) / out_res).astype(int)glt_ds[glt_y_ix, glt_x_ix] = src_ds_band# write_tiff('H:\\Datasets\\Objects\\ReadFY3D\\Output\\py_lon.tiff', [glt_lon],#            [lon_min, out_res, 0, lat_max, 0, -out_res])# write_tiff('H:\\Datasets\\Objects\\ReadFY3D\\Output\\py_lat.tiff', [glt_lat],#            [lon_min, out_res, 0, lat_max, 0, -out_res])# # 插值# interpolation_ds = np.full_like(glt_ds, fill_value=np.nan)# jump_size = windows_size // 2# for row_ix in range(jump_size, glt_rows - jump_size):#     for col_ix in range(jump_size, glt_cols - jump_size):#         if ~np.isnan(glt_ds[row_ix, col_ix]):#             interpolation_ds[row_ix, col_ix] = glt_ds[row_ix, col_ix]#             continue#         # 定義當(dāng)前窗口的邊界#         row_start = row_ix - jump_size#         row_end = row_ix + jump_size + 1  # +1 因?yàn)榍衅话Y(jié)束索引#         col_start = col_ix - jump_size#         col_end = col_ix + jump_size + 1#         rows, cols = np.ogrid[row_start:row_end, col_start:col_end]#         distances = np.sqrt((rows - row_ix) ** 2 + (cols - col_ix) ** 2)#         window_ds = glt_ds[(row_ix - jump_size):(row_ix + jump_size + 1),#                     (col_ix - jump_size):(col_ix + jump_size + 1)]#         if np.sum(~np.isnan(window_ds)) == 0:#             continue#         distances_sort_pos = np.argsort(distances.flatten())#         window_ds_sort = window_ds[np.unravel_index(distances_sort_pos, distances.shape)]#         interpolation_ds[row_ix, col_ix] = window_ds_sort[~np.isnan(window_ds_sort)][0]deal_bands.append(glt_ds)# print('處理波段: {}'.format(len(deal_bands)))# if len(deal_bands) == 6:#     breakwrite_tiff(out_path, deal_bands, [lon_min, out_res, 0, lat_max, 0, -out_res])# write_tiff('H:\\Datasets\\Objects\\ReadFY3D\\Output\\py_underlying.tiff', [interpolation_ds], [lon_min, out_res, 0, lat_max, 0, -out_res])# write_tiff('H:\\Datasets\\Objects\\ReadFY3D\\Output\\py_lon.tiff', [glt_lon], [x_min, out_res, 0, y_max, 0, -out_res])# write_tiff('H:\\Datasets\\Objects\\ReadFY3D\\Output\\py_lat.tiff', [glt_lat], [x_min, out_res, 0, y_max, 0, -out_res])

這是對(duì)重組后的數(shù)組進(jìn)行重投影,其基本思路就是對(duì)經(jīng)緯度數(shù)據(jù)集進(jìn)行zoom (congrid),將其重采樣放大譬如此處為原行列數(shù)的六倍。
接著再將zoom后的經(jīng)緯度數(shù)據(jù)集按照角點(diǎn)信息套入到輸出glt數(shù)組中,而對(duì)重組后的目標(biāo)數(shù)組直接套入,無(wú)需進(jìn)行zoom操作。
接著對(duì)套入的目標(biāo)數(shù)組進(jìn)行最近鄰插值,如果沒(méi)有進(jìn)行插值,情況如下:

在這里插入圖片描述

進(jìn)行最近鄰插值后(9*9窗口內(nèi)的最近有效像元為填充值),結(jié)果如下:

在這里插入圖片描述

但是由于數(shù)據(jù)集巨大,關(guān)于最近鄰插值此處進(jìn)行了注釋,將該操作移至后面鑲嵌操作之后,數(shù)據(jù)量大大減少,所花費(fèi)時(shí)間成本也極大降低。

2 鑲嵌和最近鄰插值

2.1 鑲嵌

def img_mosaic(mosaic_paths: list, return_transform: bool = True, mode: str = 'last'):"""該函數(shù)用于對(duì)列表中的所有TIFF文件進(jìn)行鑲嵌:param mosaic_paths: 多個(gè)TIFF文件路徑組成的字符串列表:param return_transform: 是否一同返回仿射變換:param mode: 鑲嵌模式, 默認(rèn)是Last(即如果有存在像元重疊, mosaic_paths中靠后影像的像元將覆蓋其),可選: last, mean, max, min:return:"""# 獲取鑲嵌范圍x_mins, x_maxs, y_mins, y_maxs = [], [], [], []for mosaic_path in mosaic_paths:ds = gdal.Open(mosaic_path)x_min, x_res, _, y_max, _, y_res_negative = ds.GetGeoTransform()x_size, y_size = ds.RasterXSize, ds.RasterYSizex_mins.append(x_min)x_maxs.append(x_min + x_size * x_res)y_mins.append(y_max+ y_size * y_res_negative)y_maxs.append(y_max)else:y_res = -y_res_negativeband_count = ds.RasterCountds = Nonex_min, x_max, y_min, y_max = min(x_mins), max(x_maxs), min(y_mins), max(y_maxs)# 鑲嵌col = ceil((x_max - x_min) / x_res).astype(int)row = ceil((y_max - y_min) / y_res).astype(int)mosaic_imgs = []  # 用于存儲(chǔ)各個(gè)影像for ix, mosaic_path in enumerate(mosaic_paths):mosaic_img = np.full((band_count, row, col), np.nan)  # 初始化ds = gdal.Open(mosaic_path)ds_bands = ds.ReadAsArray()# 計(jì)算當(dāng)前鑲嵌范圍start_row = floor((y_max - (y_maxs[ix] - x_res / 2)) / y_res).astype(int)start_col = floor(((x_mins[ix] + x_res / 2) - x_min) / x_res).astype(int)end_row = (start_row + ds_bands.shape[1]).astype(int)end_col = (start_col + ds_bands.shape[2]).astype(int)mosaic_img[:, start_row:end_row, start_col:end_col] = ds_bandsmosaic_imgs.append(mosaic_img)# 判斷鑲嵌模式if mode == 'last':mosaic_img = mosaic_imgs[0].copy()for img in mosaic_imgs:mask = ~np.isnan(img)mosaic_img[mask] = img[mask]elif mode == 'mean':mosaic_imgs = np.asarray(mosaic_imgs)mask = np.isnan(mosaic_imgs)mosaic_img = np.ma.array(mosaic_imgs, mask=mask).mean(axis=0).filled(np.nan)elif mode == 'max':mosaic_imgs = np.asarray(mosaic_imgs)mask = np.isnan(mosaic_imgs)mosaic_img = np.ma.array(mosaic_imgs, mask=mask).max(axis=0).filled(np.nan)elif mode == 'min':mosaic_imgs = np.asarray(mosaic_imgs)mask = np.isnan(mosaic_imgs)mosaic_img = np.ma.array(mosaic_imgs, mask=mask).min(axis=0).filled(np.nan)else:raise ValueError('不支持的鑲嵌模式: {}'.format(mode))if return_transform:return mosaic_img, [x_min, x_res, 0, y_max, 0, -y_res]return mosaic_img

這里的鑲嵌不僅僅可以解決相同地理位置的鑲嵌,也可以解決不同地理位置的拼接,拼接方式支持最大最小值計(jì)算、均值計(jì)算、Last等模式。
思路非常簡(jiǎn)單,就是將輸入的每個(gè)數(shù)據(jù)集重新裝入統(tǒng)一地理范圍的箱子(數(shù)組)中,使所有數(shù)組的角點(diǎn)信息一致,然后基于數(shù)據(jù)集個(gè)數(shù)這一維度進(jìn)行Mean、Max、Min等計(jì)算,得到鑲嵌數(shù)組。

2.2 最近鄰插值

def window_interp(arr, distances):if np.sum(~np.isnan(arr)) == 0:return np.nan# 距離最近的有效像元arr = arr.flatten()arr_sort = arr[np.argsort(distances)]if np.sum(~np.isnan(arr_sort)) == 0:return np.nanelse:return arr_sort[~np.isnan(arr_sort)][0]

思路與之前一致,不過(guò)重構(gòu)為函數(shù)了。

具體代碼見(jiàn)項(xiàng)目:https://github.com/ChaoQiezi/read_fy3d_tshs

這是原始數(shù)據(jù)集:

在這里插入圖片描述

這是目標(biāo)結(jié)果:

在這里插入圖片描述

在這里插入圖片描述

在這里插入圖片描述

http://www.risenshineclean.com/news/54959.html

相關(guān)文章:

  • 商務(wù)網(wǎng)站建設(shè)需要多少錢seo常用的優(yōu)化工具
  • PHP做的彩票網(wǎng)站好用嗎電腦培訓(xùn)中心
  • 研究政府網(wǎng)站建設(shè)的意義優(yōu)秀的網(wǎng)絡(luò)搜索引擎營(yíng)銷案例
  • win2003做網(wǎng)站百度關(guān)鍵詞怎么優(yōu)化
  • 教育類型網(wǎng)站如何優(yōu)化網(wǎng)絡(luò)
  • 定做微信小程序搜索引擎優(yōu)化的技巧
  • 樂(lè)清 做網(wǎng)站 多少錢中國(guó)經(jīng)濟(jì)網(wǎng)人事
  • 日本做電子賀卡網(wǎng)站軟件推廣怎么賺錢
  • 蘇州做網(wǎng)站設(shè)計(jì)的公司網(wǎng)站推廣的渠道有
  • ui培訓(xùn)的課程都有哪些seo建站收費(fèi)地震
  • 東莞做網(wǎng)站一年費(fèi)用百度指數(shù)明星搜索排名
  • 阿里云做網(wǎng)站流程網(wǎng)絡(luò)營(yíng)銷戰(zhàn)略的內(nèi)容
  • 網(wǎng)站電子地圖怎么做百度保障中心人工電話
  • 網(wǎng)站建設(shè)站長(zhǎng)免費(fèi)網(wǎng)絡(luò)營(yíng)銷軟件
  • 繪本館網(wǎng)站建設(shè)百度服務(wù)
  • dw網(wǎng)站模板下載西安seo顧問(wèn)培訓(xùn)
  • 網(wǎng)站換ip對(duì)優(yōu)化有影響嗎武漢百度推廣公司
  • 活在永久免費(fèi)服務(wù)器西安seo優(yōu)化顧問(wèn)
  • 網(wǎng)站上怎樣做輪播圖競(jìng)價(jià)推廣專員
  • 邯鄲建設(shè)網(wǎng)站的公司廣告引流推廣平臺(tái)
  • 網(wǎng)站開(kāi)發(fā)與維護(hù)能做什么職業(yè)互聯(lián)網(wǎng)營(yíng)銷師證書含金量
  • 做網(wǎng)站必須要備案嗎免費(fèi)行情軟件網(wǎng)站大全
  • 怎么識(shí)別網(wǎng)站是用什么語(yǔ)言做的網(wǎng)絡(luò)推廣崗位職責(zé)和任職要求
  • 統(tǒng)計(jì)網(wǎng)站建設(shè)青島設(shè)計(jì)優(yōu)化公司
  • 電子商務(wù)網(wǎng)站項(xiàng)目預(yù)算選擇寧波seo優(yōu)化公司
  • 旅游網(wǎng)站制作過(guò)程教育培訓(xùn)機(jī)構(gòu)招生方案
  • div css做網(wǎng)站微信社群營(yíng)銷
  • ui設(shè)計(jì)和網(wǎng)站開(kāi)發(fā)官方網(wǎng)站營(yíng)銷
  • 泉州開(kāi)發(fā)網(wǎng)站的公司有哪些網(wǎng)站推廣費(fèi)用
  • 如何找網(wǎng)站互聯(lián)網(wǎng)公司