Copy.py
4.1 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
# coding=utf-8
#author: 4N
#createtime: 2021/6/11
#email: nheweijun@sina.com
import copy
from osgeo.ogr import *
from osgeo import gdal,ogr
import uuid
import time
import os
def get_info_from_sqlachemy_uri(uri):
parts = uri.split(":")
user = parts[1][2:]
password_list = parts[2].split("@")
if password_list.__len__() > 2:
password = "@".join(password_list[:-1])
else:
password = parts[2].split("@")[0]
host = parts[2].split("@")[-1]
port = parts[3].split("/")[0]
database = parts[3].split("/")[1]
return user, password, host, port, database
def open_pg_data_source(iswrite, uri):
"""
# 获取PostGIS数据源
:return:
"""
db_conn_tuple = get_info_from_sqlachemy_uri(uri)
fn = "PG: user=%s password=%s host=%s port=%s dbname=%s " % db_conn_tuple
driver = ogr.GetDriverByName("PostgreSQL")
if driver is None:
raise Exception("打开PostgreSQL驱动失败，可能是当前GDAL未支持PostgreSQL驱动！")
ds = driver.Open(fn, iswrite)
if ds is None:
raise Exception("打开数据源失败！")
return ds
def move(geo:Geometry,offx,offy):
g = geo.GetGeometryRef(0)
num = g.GetPointCount()
coor = []
try:
for j in range(num):
point = g.GetPoint(j)
x = point[0]
y = point[1]
x += offx
y += offy
coor.append([x, y])
if num==1:
point = ogr.Geometry(ogr.wkbPoint)
point.AddPoint(coor[0][0], coor[0][1])
return point
elif coor[0].__eq__(coor[-1]):
ring = ogr.Geometry(ogr.wkbLinearRing)
for co in coor:
ring.AddPoint(co[0], co[1])
poly = ogr.Geometry(ogr.wkbPolygon)
poly.AddGeometry(ring)
return poly
else :
line = ogr.Geometry(ogr.wkbLineString)
for co in coor:
line.AddPoint(co[0], co[1])
return line
except:
return None
# work_dir =os.path.dirname(os.path.abspath(__file__))
work_dir =r"E:\D2\桌面\FSSD_RES_PY_FWM"
data_path=r"E:\D2\桌面\FSSD_RES_PY_FWM\fs4326.shp"
base_name="fs900w"
driver: Driver = ogr.GetDriverByName("ESRI Shapefile")
ds: DataSource = driver.Open(data_path, 1)
if not ds:
raise Exception("打开数据失败！")
layer: Layer = ds.GetLayer(0)
schema = layer.schema
schema = [s for s in schema if not s.GetName().lower().__eq__("objectid")]
for s in schema :
if s.GetName().lower().__eq__("objectid") and not s.GetTypeName().__eq__("Interger"):
s.SetType(ogr.OFTInteger)
# pg_ds:DataSource= open_pg_data_source(1,"postgresql://postgres:chinadci@172.26.99.173:5433/postgres")
# pg_layer: Layer = pg_ds.CreateLayer(base_name, layer.GetSpatialRef(), layer.GetGeomType(),["GEOMETRY_NAME=g","OVERWRITE=YES","FID=fid"])
# pg_layer.CreateFields(schema)
gdb_driver:Driver = ogr.GetDriverByName("FileGDB")
gdb_ds: DataSource = gdb_driver.CreateDataSource(os.path.join(work_dir,"{}.gdb".format(base_name)))
gdb_layer:Layer = gdb_ds.CreateLayer(base_name, layer.GetSpatialRef(),layer.GetGeomType())
# gdb_layer.CreateFields(schema)
print(gdb_layer.GetFIDColumn())
extent = layer.GetExtent()
xxrange=extent[1]-extent[0]
yrange = extent[3]-extent[2]
feature_defn: FeatureDefn = layer.GetLayerDefn()
begin = time.time()
count=0
work_dir =os.path.dirname(os.path.abspath(__file__))
for f in layer:
if count%10000==0:
print(count)
with open(os.path.join(work_dir, "copy.txt"), "w") as fi:
fi.write("已完成{}".format(count))
g:Geometry=f.GetGeometryRef()
new_f:Feature = copy.copy(f)
for xt in range(3):
for yt in range(3):
out_g = move(g,xxrange*xt,yrange*yt)
new_f.SetGeometry(out_g)
new_f.SetFID(count)
dd:FieldDefn=new_f.GetField("OBJECTID")
# new_f.UnsetField("OBJECTID")
gdb_layer.CreateFeature(new_f)
count += 1
print(time.time()-begin)
gdb_ds.Destroy()