303 lines
10 KiB
Plaintext
303 lines
10 KiB
Plaintext
{
|
||
"cells": [
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 1,
|
||
"id": "72eae2b0-1402-4fa4-a0ff-d859cf8e7bcf",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"import os\n",
|
||
"from osgeo import gdal\n",
|
||
"import numpy as np"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 2,
|
||
"id": "1fb1c1cb-6e43-442c-a8ea-e8b61a4108bc",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"def read_img(filename):\n",
|
||
" data = gdal.Open(filename) # 打开文件\n",
|
||
" im_width = data.RasterXSize # 读取图像行数\n",
|
||
" im_height = data.RasterYSize # 读取图像列数\n",
|
||
"\n",
|
||
" im_geotrans = data.GetGeoTransform()\n",
|
||
" # 仿射矩阵，左上角像素的大地坐标和像素分辨率。\n",
|
||
" # 共有六个参数，分表代表左上角x坐标；东西方向上图像的分辨率；如果北边朝上，地图的旋转角度，0表示图像的行与x轴平行；左上角y坐标；\n",
|
||
" # 如果北边朝上，地图的旋转角度，0表示图像的列与y轴平行；南北方向上地图的分辨率。\n",
|
||
" im_proj = data.GetProjection() # 地图投影信息\n",
|
||
" im_data = data.ReadAsArray(0, 0, im_width, im_height) # 此处读取整张图像\n",
|
||
" # ReadAsArray(<xoff>, <yoff>, <xsize>, <ysize>)\n",
|
||
" # 读出从(xoff,yoff)开始，大小为(xsize,ysize)的矩阵。\n",
|
||
" del data # 释放内存，如果不释放，在arcgis，envi中打开该图像时会显示文件被占用\n",
|
||
"\n",
|
||
" return im_proj, im_geotrans, im_data"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 3,
|
||
"id": "24c38f93-9d68-46b7-a7f8-0124d98384e2",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"def readtxt(file_path):\n",
|
||
" print(\"read file path:\\t\",file_path)\n",
|
||
" with open(file_path,'r',encoding='utf-8') as fp:\n",
|
||
" line=fp.readlines()\n",
|
||
" print(line)\n",
|
||
" return line"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 4,
|
||
"id": "04e4be06-3561-4e74-975a-5e8687cb1385",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"def writetiff(outname,data,nl,ns,im_proj,im_geotrans,nband=2):\n",
|
||
" # 创建.tif文件\n",
|
||
" driver = gdal.GetDriverByName(\"GTiff\")\n",
|
||
" out_tif = driver.Create(outname, ns, nl, nband, gdal.GDT_Float32)\n",
|
||
" #geotransform = (LonMin, Lon_Res, 0, LatMax, 0, -Lat_Res)\n",
|
||
" out_tif.SetGeoTransform(im_geotrans) \n",
|
||
" # 获取地理坐标系统信息，用于选取需要的地理坐标系统 \n",
|
||
" #srs = osr.SpatialReference()\n",
|
||
" #proj_type = 'GEOGCS[\"WGS 84\",DATUM[\"WGS_1984\",SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]],AUTHORITY[\"EPSG\",\"6326\"]],PRIMEM[\"Greenwich\",0],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AXIS[\"Latitude\",NORTH],AXIS[\"Longitude\",EAST],AUTHORITY[\"EPSG\",\"4326\"]]'\n",
|
||
" out_tif.SetProjection(im_proj) # 给新建图层赋予投影信息 \n",
|
||
" # 数据写出 \n",
|
||
" for i in range(nband):\n",
|
||
" out_tif.GetRasterBand(i+1).WriteArray(data[i]) # 将数据写入内存，此时没有写入硬盘\n",
|
||
" out_tif.FlushCache() # 将数据写入硬盘\n",
|
||
" out_tif = None # 注意必须关闭tif文件"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 5,
|
||
"id": "919de6ec-8441-4fa1-bbc3-4c9a7f1ec0cb",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"def writetiff_NoProj(outname,data,nl,ns,nband=2):\n",
|
||
" # 创建.tif文件\n",
|
||
" driver = gdal.GetDriverByName(\"GTiff\")\n",
|
||
" out_tif = driver.Create(outname,nl,ns,Lon_Res,Lat_Res,LonMin, LatMax, LonMax, nband, gdal.GDT_Float32)\n",
|
||
" #geotransform = (LonMin, Lon_Res, 0, LatMax, 0, -Lat_Res)\n",
|
||
" #out_tif.SetGeoTransform(geotransform) \n",
|
||
" # 获取地理坐标系统信息，用于选取需要的地理坐标系统 \n",
|
||
" #srs = osr.SpatialReference()\n",
|
||
" #proj_type = 'GEOGCS[\"WGS 84\",DATUM[\"WGS_1984\",SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]],AUTHORITY[\"EPSG\",\"6326\"]],PRIMEM[\"Greenwich\",0],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AXIS[\"Latitude\",NORTH],AXIS[\"Longitude\",EAST],AUTHORITY[\"EPSG\",\"4326\"]]'\n",
|
||
" #out_tif.SetProjection(proj_type) # 给新建图层赋予投影信息 \n",
|
||
" # 数据写出 \n",
|
||
" for i in range(nband):\n",
|
||
" out_tif.GetRasterBand(i+1).WriteArray(data[i]) # 将数据写入内存，此时没有写入硬盘\n",
|
||
" out_tif.FlushCache() # 将数据写入硬盘\n",
|
||
" out_tif = None # 注意必须关闭tif文件"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 6,
|
||
"id": "af5414eb-f97a-4811-9f32-a4b533b513eb",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"master_path=r\"D:\\backsigma\\microproduct\\vegetationHeight-L-SAR-V1.0\\vegetationPhenology\\track_master\"\n",
|
||
"slave_path=r\"D:\\backsigma\\microproduct\\vegetationHeight-L-SAR-V1.0\\vegetationPhenology\\track_slave\"\n",
|
||
"file_list=[\"s11.bin\",\"s12.bin\",\"s21.bin\",\"s22.bin\"]\n",
|
||
"config_path=r\"config.txt\"\n",
|
||
"mapinfo_path=r\"config_mapinfo.txt\""
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 7,
|
||
"id": "4100846c-faab-48e8-af6c-f86e02686ad1",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"read file path:\t D:\\backsigma\\microproduct\\vegetationHeight-L-SAR-V1.0\\vegetationPhenology\\track_master\\config.txt\n",
|
||
"['Nrow\\n', '105\\n', '---------\\n', 'Ncol\\n', '141\\n', '---------\\n', 'PolarCase\\n', 'monostatic\\n', '---------\\n', 'PolarType\\n', 'full\\n']\n",
|
||
"read file path:\t D:\\backsigma\\microproduct\\vegetationHeight-L-SAR-V1.0\\vegetationPhenology\\track_master\\config_mapinfo.txt\n",
|
||
"['Sensor\\n', 'NONE\\n', '---------\\n', 'MapInfo\\n', 'map info = {UTM,1,1,0.0,0.0,1.0,1.0,30,North}\\n', '---------\\n', 'MapProj\\n', 'UTM\\n', '1\\n', '1\\n', '1.0\\n', '1.0\\n']\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"config_text=readtxt(os.path.join(master_path,config_path))\n",
|
||
"mapinfo_text=readtxt(os.path.join(master_path,mapinfo_path))"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 9,
|
||
"id": "801b6bad-4fa5-422c-8894-97d3cc995874",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"for file in file_list:\n",
|
||
" filepath=os.path.join(master_path,file)\n",
|
||
" im_proj, im_geotrans, im_data=read_img(filepath)\n",
|
||
" data=[im_data.real,im_data.imag]\n",
|
||
" writetiff_NoProj(os.path.join(master_path,file.replace(\".bin\",\".tiff\")),data,105,141,nband=2)\n",
|
||
" "
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 10,
|
||
"id": "2c3be0df-cc37-43ff-8cd7-704a903438a7",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"(105, 141)"
|
||
]
|
||
},
|
||
"execution_count": 10,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"im_data.shape"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 11,
|
||
"id": "d2d4b4fb-89f2-4eb9-b93a-7b1ae7f29682",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"for file in file_list:\n",
|
||
" filepath=os.path.join(slave_path,file)\n",
|
||
" im_proj, im_geotrans, im_data=read_img(filepath)\n",
|
||
" data=[im_data.real,im_data.imag]\n",
|
||
" writetiff_NoProj(os.path.join(slave_path,file.replace(\".bin\",\".tiff\")),data,105,141,nband=2)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 25,
|
||
"id": "88faf53c-49e9-453b-b345-728701b82ede",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"c=np.ones((105,141))\n",
|
||
"c=c*np.array(list(range(141)))\n",
|
||
"r=np.ones((105,141)).transpose()\n",
|
||
"r=(r*np.array(list(range(105)))).transpose()"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 28,
|
||
"id": "7631cf12-160f-4be6-af62-240df8b70ff3",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"array([[ 0., 1., 2., ..., 138., 139., 140.],\n",
|
||
" [ 0., 1., 2., ..., 138., 139., 140.],\n",
|
||
" [ 0., 1., 2., ..., 138., 139., 140.],\n",
|
||
" ...,\n",
|
||
" [ 0., 1., 2., ..., 138., 139., 140.],\n",
|
||
" [ 0., 1., 2., ..., 138., 139., 140.],\n",
|
||
" [ 0., 1., 2., ..., 138., 139., 140.]])"
|
||
]
|
||
},
|
||
"execution_count": 28,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"c"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 26,
|
||
"id": "cc788f4a-96f3-4723-842c-7973ee998e1c",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"(105, 141) (105, 141)\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"print(c.shape,r.shape)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 29,
|
||
"id": "1027bafc-892d-4d03-99f2-0b030d433e9b",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"data=[r,c]\n",
|
||
"writetiff(os.path.join(master_path,\"ori_sim.tif\"),data,105,141,im_proj,im_geotrans,2)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 30,
|
||
"id": "84aab2a4-3b48-4aac-a3a5-b804e8f431ed",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"data=[np.ones((105,141))*45.0]\n",
|
||
"writetiff(os.path.join(master_path,\"out_orth_sar_incidence.tif\"),data,105,141,im_proj,im_geotrans,1)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "63191526-29c4-4e7c-bb63-30bc9c5c6788",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"data=[np.ones((105,141))*0]\n",
|
||
"writetiff(os.path.join(master_path,\"ori_sim.tif\"),data,105,141,im_proj,im_geotrans,1)"
|
||
]
|
||
}
|
||
],
|
||
"metadata": {
|
||
"kernelspec": {
|
||
"display_name": "Python 3 (ipykernel)",
|
||
"language": "python",
|
||
"name": "python3"
|
||
},
|
||
"language_info": {
|
||
"codemirror_mode": {
|
||
"name": "ipython",
|
||
"version": 3
|
||
},
|
||
"file_extension": ".py",
|
||
"mimetype": "text/x-python",
|
||
"name": "python",
|
||
"nbconvert_exporter": "python",
|
||
"pygments_lexer": "ipython3",
|
||
"version": "3.8.13"
|
||
}
|
||
},
|
||
"nbformat": 4,
|
||
"nbformat_minor": 5
|
||
}
|