Cách thực hiện tính toán GIS trên cơ sở dữ liệu không phải GIS

Giới thiệu:

Đôi khi, bạn có thể gặp phải nhu cầu thực hiện các chức năng không gian địa lý trong ứng dụng của mình, chẳng hạn như ánh xạ vị trí người dùng hoặc phân tích dữ liệu địa lý. Có rất nhiều thư viện dành riêng cho ngôn ngữ có sẵn cho các tác vụ này, chẳng hạn như GDAL, Shapely và Geopandas cho Python.

Ngoài ra, chức năng không gian địa lý có thể được thực hiện thông qua cơ sở dữ liệu; ví dụ: bạn có thể sử dụng tiện ích mở rộng PostGIS với cơ sở dữ liệu quan hệ như PostgreSQL hoặc tận dụng hỗ trợ riêng cho các loại dữ liệu không gian trong cơ sở dữ liệu phân tán như Azure CosmosDB.

Tuy nhiên, nếu bộ lưu trữ phụ trợ của bạn, chẳng hạn như Redis hoặc Google Spanner, về cơ bản không hỗ trợ các truy vấn không gian và bạn cần xử lý các truy vấn không gian địa lý quy mô lớn thì bài viết này được thiết kế riêng cho bạn.

Những lựa chọn của tôi là gì?

Bạn luôn có thể xây dựng một vi dịch vụ khác để xử lý dữ liệu không gian, nhưng tùy chọn này thường đòi hỏi chi phí duy trì một ứng dụng bổ sung. Một cách tiếp cận khác là sử dụng các thư viện lập chỉ mục địa lý như S2 và H3. S2, do Google phát triển, dựa trên đường cong Hilbert, trong khi H3, do Uber phát triển, dựa trên hệ thống lưới toàn cầu rời rạc trắc địa. S2 và H3 có nhiều điểm tương đồng: cả hai đều chia một vùng nhất định thành các ô và sử dụng số nguyên 64 bit để lập chỉ mục cho các ô này.

Tuy nhiên, sự khác biệt chính nằm ở hình dạng của tế bào; S2 sử dụng các ô hình vuông, trong khi H3 sử dụng các ô hình lục giác. Đối với một số ứng dụng, H3 có thể mang lại hiệu suất tốt hơn. Tuy nhiên, nhìn chung, một trong hai thư viện là đủ. Trong bài viết này, chúng tôi sẽ sử dụng S2, nhưng bạn có thể thực hiện các chức năng tương tự bằng H3.

Các khái niệm cơ bản về Thư viện Google S2

• Ô: S2 chia hình cầu thành các ô, mỗi ô có một mã định danh 64 bit duy nhất.

• Cấp độ ô: Hệ thống phân cấp cho phép các cấp độ chi tiết khác nhau, từ các vùng lớn đến các vùng chính xác nhỏ. Mỗi cấp độ đại diện cho một độ phân giải khác nhau:

  
  

  • Cấp độ 0: Các tế bào lớn nhất, bao phủ một phần đáng kể bề mặt Trái đất.

    
    

  • Cấp độ cao hơn: Các ô được chia dần thành các góc phần tư nhỏ hơn. Ví dụ: mỗi ô Cấp 1 được chia thành bốn ô Cấp 2, v.v.

    
    

  • Độ phân giải và Diện tích: Mức cao hơn tương ứng với độ phân giải tốt hơn và diện tích ô nhỏ hơn. Hệ thống phân cấp này cho phép lập chỉ mục và truy vấn chính xác ở các mức độ chi tiết khác nhau.

    
    

Bảng dưới đây hiển thị các cấp độ ô khác nhau cùng với các khu vực tương ứng của chúng.

Từ bảng được cung cấp, rõ ràng là bạn có thể đạt được độ chính xác ánh xạ xuống tới 0,44 cm^2 bằng S2. Trong mỗi ô vuông của ô S2, tồn tại một ô con có chung ô cha, biểu thị cấu trúc phân cấp. Cấp độ của ô có thể là giá trị tĩnh (cùng mức áp dụng cho tất cả các ô) hoặc có thể động trong đó S2 quyết định độ phân giải nào hoạt động tốt nhất.

Tính toán hàng xóm gần nhất

Hãy bắt đầu với một ví dụ. Hãy xem xét việc chúng tôi đang viết một ứng dụng cung cấp các tính năng giống như dịch vụ lân cận cho khu vực Seattle. Chúng tôi muốn trả về danh sách các quán cà phê trong vùng lân cận nhất định. Để thực hiện các thao tác này, chúng ta sẽ chia nhiệm vụ này thành 4 nhiệm vụ phụ:

• Đang tải bản đồ Seattle
• Trực quan hóa các ô S2 trên bản đồ Seattle
• Lưu trữ một số địa điểm quán cà phê trong cơ sở dữ liệu
• Truy vấn quán cà phê gần nhất

Đang tải bản đồ Seattle

Để tải bản đồ Google, chúng tôi sẽ sử dụng thư viện gmplot. Thư viện này yêu cầu khóa API Google Maps để tải. Để tạo khóa API, hãy làm theo hướng dẫn tại đây .

 import gmplot import const # plot seattle with zoom level 13 gmap = gmplot.GoogleMapPlotter(47.6061, -122.3328, 13, apikey=const.API_KEY) # Draw the map to an HTML file: gmap.draw('map.html')

Đoạn mã trên tạo ra một tệp map.html như dưới đây:

Trực quan hóa các ô S2 trên Bản đồ Seattle

Bây giờ chúng ta đã có bản đồ, hãy vẽ một số ô S2 cho bản đồ:

 from s2 import * import gmplot # plot seattle with zoom level 13 gmap = gmplot.GoogleMapPlotter(47.6061, -122.3328, 13, apikey=const.API_KEY) areatobeplotted = [ (47.64395531736767,-122.43597221319135), (47.51369277846956,-122.43597221319135), (47.51369277846956,-122.24156866779164), (47.64395531736767,-122.24156866779164), (47.64395531736767,-122.43597221319135) ] region_rect = S2LatLngRect( S2LatLng.FromDegrees(47.51369277846956,-122.43597221319135), S2LatLng.FromDegrees(47.64395531736767, -122.24156866779164)) coverer = S2RegionCoverer() coverer.set_min_level(8) coverer.set_max_level(15) covering = coverer.GetCovering(region_rect) geoms = 0 for cellid in covering: new_cell = S2Cell(cellid) vertices = [] for i in range(0, 4): vertex = new_cell.GetVertex(i) latlng = S2LatLng(vertex) vertices.append((latlng.lat().degrees(), latlng.lng().degrees())) gmap.polygon(*zip(*vertices), face_color='pink', edge_color='cornflowerblue', edge_width=5) geoms+=1 gmap.polygon(*zip(*areatobeplotted), face_color='red', edge_color='green', edge_width=5) print(f"Total Geometries: {geoms}") gmap.draw('/tmp/map.html')

 Output: Total Geometries: 273

Trong đoạn mã trên, trước tiên chúng tôi căn giữa máy vẽ Google Map xung quanh khu vực Seattle. Trong S2RegionCoverer , chúng tôi khởi tạo trình bao phủ vùng để có các mức động trong khoảng từ mức tối thiểu là 8 đến mức tối đa là 15. Điều này cho phép S2 tự động điều chỉnh tất cả các ô vào các kích thước ô cụ thể để phù hợp nhất. Phương thức GetCovering trả về lớp phủ cho một hình chữ nhật xung quanh khu vực Seattle.

Sau đó, chúng tôi lặp lại từng ô, tính toán các đỉnh cho các ô và vẽ chúng trên bản đồ. Chúng tôi giữ số lượng ô được tạo ở khoảng 273. Cuối cùng, chúng tôi vẽ hình chữ nhật đầu vào bằng màu đỏ. Mã này sẽ vẽ các ô S2 trên bản đồ Seattle tại /tmp/map.html , như hiển thị bên dưới:

Lưu trữ một vài địa điểm quán cà phê trong cơ sở dữ liệu

Hãy tạo cơ sở dữ liệu về các quán cà phê cùng với mã định danh ô S2 của họ. Bạn có thể lưu trữ các ô này trong bất kỳ cơ sở dữ liệu nào bạn chọn. Đối với hướng dẫn này, chúng tôi sẽ sử dụng cơ sở dữ liệu dữ liệu SQLite. Trong mẫu mã bên dưới, chúng tôi kết nối với cơ sở dữ liệu SQLite để tạo bảng CoffeeShops với 3 trường Id , name và cell_id .

Tương tự như ví dụ trước, chúng tôi sử dụng S2RegionCoverer để tính toán các ô nhưng lần này, chúng tôi sử dụng một mức cố định để vẽ các điểm. Cuối cùng, ID được tính toán sẽ được chuyển đổi thành chuỗi và được lưu trữ trong cơ sở dữ liệu.

 import sqlite3 from s2 import S2CellId,S2LatLng,S2RegionCoverer # Connect to SQLite database conn = sqlite3.connect('/tmp/sqlite_cells.db') cursor = conn.cursor() # Create a table to store cell IDs cursor.execute('''CREATE TABLE IF NOT EXISTS CoffeeShops ( id INTEGER PRIMARY KEY AUTOINCREMENT, name TEXT, cell_id TEXT )''') coverer = S2RegionCoverer() # Function to generate S2 cell ID for a given latitude and longitude def generate_cell_id(latitude, longitude, level=16): cell=S2CellId(S2LatLng.FromDegrees(latitude, longitude)) return str(cell.parent(level)) # Function to insert cell IDs into the database def insert_cell_ids(name,lat,lng): cell_id = generate_cell_id(lat, lng) cursor.execute("INSERT INTO CoffeeShops (name, cell_id) VALUES (?, ?)", (name, cell_id)) conn.commit() # Insert cell IDs into the database insert_cell_ids("Overcast Coffee", 47.616656277302155, -122.31156460382837) insert_cell_ids("Seattle Sunshine", 47.67366852914391, -122.29051997415843) insert_cell_ids("Sip House", 47.6682364706238, -122.31328618043693) insert_cell_ids("Victoria Coffee",47.624408595334536, -122.3117362652041) # Close connection conn.close()

Tại thời điểm này, chúng tôi có cơ sở dữ liệu lưu trữ các quán cà phê cùng với ID ô của họ, được xác định bởi độ phân giải đã chọn cho cấp ô.

Truy vấn các quán cà phê gần nhất

Cuối cùng, hãy truy vấn các quán cà phê trong khu vực University District.

 import sqlite3 from s2 import S2RegionCoverer,S2LatLngRect, S2LatLng # Connect to SQLite database conn = sqlite3.connect('/tmp/sqlite_cells.db') cursor = conn.cursor() # Function to query database for cells intersecting with the given polygon def query_intersecting_cells(start_x,start_y,end_x,end_y): # Create S2RegionCoverer region_rect = S2LatLngRect( S2LatLng.FromDegrees(start_x,start_y), S2LatLng.FromDegrees(end_x,end_y)) coverer = S2RegionCoverer() coverer.set_min_level(8) coverer.set_max_level(15) covering = coverer.GetCovering(region_rect) # Query for intersecting cells intersecting_cells = set() for cell_id in covering: cursor.execute("SELECT name FROM CoffeeShops WHERE cell_id >= ? and cell_id<=?", (str(cell_id.range_min()),str(cell_id.range_max()),)) intersecting_cells.update(cursor.fetchall()) return intersecting_cells # Query for intersecting cells intersecting_cells = query_intersecting_cells(47.6527847,-122.3286438,47.6782181, -122.2797203) # Print intersecting cells print("Intersecting cells:") for cell_id in intersecting_cells: print(cell_id[0]) # Close connection conn.close()

 Output: Intersecting cells: Sip House Seattle Sunshine

Dưới đây là hình ảnh trực quan của các ô. Để duy trì sự ngắn gọn, mã trực quan bên dưới không được thêm vào.

Vì tất cả các ô con và ô cha đều có chung một tiền tố nên chúng ta có thể truy vấn các phạm vi ô trong khoảng từ tối thiểu đến tối đa để lấy tất cả các ô nằm giữa hai giá trị đó. Trong ví dụ của chúng tôi, chúng tôi sử dụng nguyên tắc tương tự để truy vấn quán cà phê

Phần kết luận:

Trong bài viết này, chúng tôi đã trình bày cách sử dụng lập chỉ mục địa lý để lưu trữ và truy vấn dữ liệu không gian địa lý trong cơ sở dữ liệu không hỗ trợ truy vấn không gian địa lý. Điều này có thể được mở rộng hơn nữa cho nhiều trường hợp sử dụng như tính toán định tuyến giữa 2 điểm hoặc tìm các điểm lân cận gần nhất.

Thông thường, để truy vấn cơ sở dữ liệu được lập chỉ mục địa lý, bạn sẽ phải thực hiện một số xử lý hậu kỳ bổ sung trên dữ liệu. Cần phải xem xét cẩn thận logic truy vấn và xử lý hậu kỳ để đảm bảo chúng tôi không áp đảo nút.

Người giới thiệu:

• Máy vẽ bản đồ Google - https://github.com/gmplot/gmplot/wiki/GoogleMapPlotter
• Hướng dẫn dành cho nhà phát triển S2 - http://s2geometry.io/devguide/
• Sqlite - https://docs.python.org/3/library/sqlite3.html
• Không gian địa lý Azure Cosmos DB - https://learn.microsoft.com/en-us/azure/cosmos-db/nosql/query/geospatial?tabs=javascript
• GDAL - https://gdal.org/index.html
• Shapley- https://shapely.readthedocs.io/en/stable/manual.html
• Geopandas - https://geopandas.org
• Posgis - https://postgis.net/
• Cờ lê Google - https://cloud.google.com/spanner?hl=vi
• Redis - https://redis.io/