feat(sf): geocode geom для ekburg_construction_permits из МСК-66 (#1088) #1124

Merged
bot-reviewer merged 2 commits from feat/ekburg-permits-geocode into main 2026-06-07 11:16:46 +00:00
3 changed files with 402 additions and 8 deletions

View file

@ -31,9 +31,12 @@ Format: Form 3 (РНС) + Form 4 (РВЭ) — стандартные Росст
15 rve_number Номер разрешения на ввод
16 rve_date Дата разрешения на ввод
Координаты raw_coord_x / raw_coord_y хранятся как строки CRS неизвестна (предположительно
СНСК-66 / СКФ-66 Свердловская обл., local offset system). Геокодирование через construction_address
планируется в Phase 3 (отдельный PR).
Координаты raw_coord_x / raw_coord_y СКФ-66 (МСК-66, Свердловская обл., зона 1).
Параметры проекции: tmerc, lon_0=60°E, k=1, false_easting=1496944, false_northing=-5911077,
эллипсоид Красовского, сдвиг datum towgs84=23.92,-141.27,-80.9,0,0.35,0.82,-0.12.
Параметры откалиброваны по cad_parcels ground-truth (738 точек): avg_err14м, median7м.
Ячейки могут содержать несколько точек через пробел берём первую.
Геокодирование выполняется при парсинге (geom, geocoded_lat/lon).
"""
from __future__ import annotations
@ -48,9 +51,89 @@ from typing import Any
import httpx
from openpyxl import load_workbook
from pyproj import CRS, Transformer
logger = logging.getLogger(__name__)
# ── МСК-66 → WGS84 трансформер ───────────────────────────────────────────────
# СКФ-66 (местная СК Свердловской обл., зона 1):
# tmerc, осевой меридиан 60°E, масштаб 1.0
# x_0=1496944 (false_easting), y_0=-5911077 (false_northing)
# Откалиброваны по cad_parcels ground-truth (738 точек permit'ов с кад.номером):
# avg_err≈14м, median≈7м, p95≈58м (без 8 записей с ошибочными Y в исходных данных).
# Предыдущее x_0=1500000 давало систематический сдвиг ~3050м западнее по долготе.
# эллипсоид Красовского + сдвиг Пулково→WGS84
_MSK66_PROJ4 = (
"+proj=tmerc +lat_0=0 +lon_0=60 +k=1"
" +x_0=1496944 +y_0=-5911077"
" +ellps=krass +towgs84=23.92,-141.27,-80.9,0,0.35,0.82,-0.12"
" +units=m +no_defs"
)
_MSK66_TO_WGS84: Transformer = Transformer.from_crs(
CRS.from_proj4(_MSK66_PROJ4),
CRS.from_epsg(4326),
always_xy=True,
)
# Валидный bbox ЕКБ и ближайших районов (WGS84)
_EKB_LON_MIN = 59.5
_EKB_LON_MAX = 61.5
_EKB_LAT_MIN = 56.3
_EKB_LAT_MAX = 57.5
# Паттерн для разбора первого числа из строки координат
# (ячейки могут содержать несколько точек через пробелы или запятую-десятичный разделитель)
_COORD_FIRST_RE = re.compile(r"[\d]+[.,][\d]+")
def msk66_to_wgs84(raw_x: str | None, raw_y: str | None) -> tuple[float, float] | None:
"""Преобразовать первую точку МСК-66 (raw_coord_x, raw_coord_y) в (lon, lat) WGS84.
Берёт первое числовое значение из каждой ячейки (в xlsx могут быть несколько
точек через пробел). Возвращает None если координаты отсутствуют, не парсятся
или результат выходит за пределы ЕКБ-региона.
В МСК-66 Свердловской обл.:
raw_coord_x 1 500 000 + смещение (откалиброванный false_easting = 1 496 944)
raw_coord_y northing + смещение (откалиброванный false_northing = -5 911 077)
"""
if not raw_x or not raw_y:
return None
# Нормализуем запятую как десятичный разделитель
x_str = raw_x.replace(",", ".")
y_str = raw_y.replace(",", ".")
# Берём первое число из (возможно) многоточечной строки
mx = _COORD_FIRST_RE.search(x_str)
my = _COORD_FIRST_RE.search(y_str)
if not mx or not my:
return None
try:
msk_x = float(mx.group())
msk_y = float(my.group())
except ValueError:
return None
# Грубая проверка диапазона МСК-66 X (false_easting=1500000, типичные отклонения ±50 км)
if not (1_450_000 <= msk_x <= 1_600_000):
logger.debug("msk66_to_wgs84: X=%s вне ожидаемого диапазона МСК-66", msk_x)
return None
try:
lon, lat = _MSK66_TO_WGS84.transform(msk_x, msk_y)
except Exception as exc: # pyproj может поднять ProjError
logger.debug("msk66_to_wgs84: transform failed X=%s Y=%s: %s", msk_x, msk_y, exc)
return None
if not (_EKB_LON_MIN <= lon <= _EKB_LON_MAX and _EKB_LAT_MIN <= lat <= _EKB_LAT_MAX):
logger.debug("msk66_to_wgs84: результат (%.4f, %.4f) вне ЕКБ-региона", lon, lat)
return None
return lon, lat
EKBURG_PERMITS_URLS: dict[int, str] = {
2026: "https://xn--80acgfbsl1azdqr.xn--p1ai/file/70bf01bf31538ee9dd82dadfc47192a0",
2025: "https://xn--80acgfbsl1azdqr.xn--p1ai/file/6a0a18c9ee327e6e4f76c32a5385a6bd",
@ -87,6 +170,8 @@ class PermitRow:
rve_date: date | None # только РВЭ
raw_coord_x: str | None
raw_coord_y: str | None
geocoded_lat: float | None # WGS84 широта (из МСК-66, None если не сконвертировано)
geocoded_lon: float | None # WGS84 долгота (из МСК-66, None если не сконвертировано)
source_year: int
source_url: str
raw_row: dict[str, Any] = field(default_factory=dict)
@ -214,6 +299,12 @@ def _parse_row(
raw: dict[str, Any] = {str(i): str(v) for i, v in enumerate(row) if v is not None}
raw_x = _to_str(row[6]) if len(row) > 6 else None
raw_y = _to_str(row[7]) if len(row) > 7 else None
coords = msk66_to_wgs84(raw_x, raw_y)
geocoded_lon: float | None = coords[0] if coords else None
geocoded_lat: float | None = coords[1] if coords else None
return PermitRow(
permit_type=permit_type,
permit_number=permit_number,
@ -227,11 +318,15 @@ def _parse_row(
cadastral_number=_to_str(row[5]) if len(row) > 5 else None,
total_area_sqm=_to_float(row[12]) if len(row) > 12 else None,
living_area_sqm=_to_float(row[13]) if len(row) > 13 else None,
living_area_fact_sqm=_to_float(row[14]) if permit_type == "RVE" and len(row) > 14 else None,
living_area_fact_sqm=(
_to_float(row[14]) if permit_type == "RVE" and len(row) > 14 else None
),
rve_number=_to_str(row[15]) if permit_type == "RVE" and len(row) > 15 else None,
rve_date=_to_date(row[16]) if permit_type == "RVE" and len(row) > 16 else None,
raw_coord_x=_to_str(row[6]) if len(row) > 6 else None,
raw_coord_y=_to_str(row[7]) if len(row) > 7 else None,
raw_coord_x=raw_x,
raw_coord_y=raw_y,
geocoded_lat=geocoded_lat,
geocoded_lon=geocoded_lon,
source_year=year,
source_url=source_url,
raw_row=raw,

View file

@ -17,16 +17,31 @@ from app.services.scrapers.ekburg_permits import (
EKBURG_PERMITS_URLS,
EkburgPermitsClient,
PermitRow,
msk66_to_wgs84,
)
from app.workers.celery_app import celery_app
logger = logging.getLogger(__name__)
# Размер батча для backfill_geom (SAVEPOINT per batch)
_BACKFILL_BATCH_SIZE = 200
def _upsert_permit(db: Session, row: PermitRow) -> None:
"""UPSERT одной строки разрешения в ekburg_construction_permits."""
"""UPSERT одной строки разрешения в ekburg_construction_permits.
Если geocoded_lat/lon заполнены записывает geom как ST_SetSRID(ST_MakePoint(lon,lat),4326).
"""
if row.geocoded_lon is not None and row.geocoded_lat is not None:
geom_expr = (
"ST_SetSRID(ST_MakePoint(CAST(:geocoded_lon AS double precision),"
" CAST(:geocoded_lat AS double precision)), 4326)"
)
else:
geom_expr = "NULL"
db.execute(
text("""
text(f"""
INSERT INTO ekburg_construction_permits (
permit_type, permit_number,
issue_date, expiry_date,
@ -36,6 +51,7 @@ def _upsert_permit(db: Session, row: PermitRow) -> None:
total_area_sqm, living_area_sqm, living_area_fact_sqm,
rve_number, rve_date,
raw_coord_x, raw_coord_y,
geocoded_lat, geocoded_lon, geom,
source_year, source_url, raw_row
)
VALUES (
@ -47,6 +63,7 @@ def _upsert_permit(db: Session, row: PermitRow) -> None:
:total_area_sqm, :living_area_sqm, :living_area_fact_sqm,
:rve_number, :rve_date,
:raw_coord_x, :raw_coord_y,
:geocoded_lat, :geocoded_lon, {geom_expr},
:source_year, :source_url, CAST(:raw_row AS jsonb)
)
ON CONFLICT (permit_type, permit_number) DO UPDATE SET
@ -65,6 +82,9 @@ def _upsert_permit(db: Session, row: PermitRow) -> None:
rve_date = EXCLUDED.rve_date,
raw_coord_x = EXCLUDED.raw_coord_x,
raw_coord_y = EXCLUDED.raw_coord_y,
geocoded_lat = EXCLUDED.geocoded_lat,
geocoded_lon = EXCLUDED.geocoded_lon,
geom = EXCLUDED.geom,
raw_row = EXCLUDED.raw_row,
fetched_at = NOW()
"""),
@ -86,6 +106,8 @@ def _upsert_permit(db: Session, row: PermitRow) -> None:
"rve_date": row.rve_date,
"raw_coord_x": row.raw_coord_x,
"raw_coord_y": row.raw_coord_y,
"geocoded_lat": row.geocoded_lat,
"geocoded_lon": row.geocoded_lon,
"source_year": row.source_year,
"source_url": row.source_url,
"raw_row": json.dumps(row.raw_row, ensure_ascii=False),
@ -93,6 +115,107 @@ def _upsert_permit(db: Session, row: PermitRow) -> None:
)
def backfill_geom(db: Session | None = None) -> dict[str, int]:
"""Бэкфилл geom/geocoded_lat/lon для существующих строк с raw_coord_x/y но без geom.
Выполняется разово (admin-trigger или параметр refresh). Использует батч-UPDATE
с SAVEPOINT per batch для устойчивости.
Возвращает {"updated": N, "skipped": N, "errors": N}.
"""
own_session = db is None
if own_session:
db = SessionLocal()
updated = 0
skipped = 0
errors = 0
try:
# Выбираем строки с координатами, у которых geom пуст
rows = db.execute( # type: ignore[union-attr]
text("""
SELECT id, raw_coord_x, raw_coord_y
FROM ekburg_construction_permits
WHERE raw_coord_x IS NOT NULL
AND raw_coord_y IS NOT NULL
AND geom IS NULL
ORDER BY id
""")
).fetchall()
logger.info("ekburg_permits backfill_geom: %d строк для обработки", len(rows))
batch: list[dict] = []
def _flush_batch(b: list[dict]) -> tuple[int, int]:
u = e = 0
try:
with db.begin_nested(): # type: ignore[union-attr]
for params in b:
db.execute( # type: ignore[union-attr]
text("""
UPDATE ekburg_construction_permits
SET
geocoded_lat = CAST(:lat AS double precision),
geocoded_lon = CAST(:lon AS double precision),
geom = ST_SetSRID(
ST_MakePoint(
CAST(:lon AS double precision),
CAST(:lat AS double precision)
),
4326
),
fetched_at = NOW()
WHERE id = CAST(:id AS bigint)
"""),
params,
)
u += 1
except Exception as exc:
logger.warning("backfill_geom: batch failed: %s", exc)
e += len(b)
u = 0
return u, e
for row in rows:
coords = msk66_to_wgs84(row.raw_coord_x, row.raw_coord_y)
if coords is None:
skipped += 1
continue
lon, lat = coords
batch.append({"id": row.id, "lon": lon, "lat": lat})
if len(batch) >= _BACKFILL_BATCH_SIZE:
u, e = _flush_batch(batch)
updated += u
errors += e
batch = []
if batch:
u, e = _flush_batch(batch)
updated += u
errors += e
db.commit() # type: ignore[union-attr]
except Exception as exc:
logger.error("backfill_geom: unexpected error: %s", exc)
raise
finally:
if own_session:
db.close() # type: ignore[union-attr]
logger.info(
"ekburg_permits backfill_geom done: updated=%d skipped=%d errors=%d",
updated,
skipped,
errors,
)
return {"updated": updated, "skipped": skipped, "errors": errors}
@celery_app.task(name="tasks.ekburg_permits_sync.refresh_year", queue="celery")
def refresh_year(year: int) -> dict[str, int]:
"""Скачать + распарсить + upsert РНС/РВЭ за один год.
@ -149,3 +272,13 @@ def refresh_all() -> dict[str, dict[str, int]]:
results[str(year)] = refresh_year(year)
logger.info("ekburg_permits_sync: refresh_all done: %s", results)
return results
@celery_app.task(name="tasks.ekburg_permits_sync.backfill_geom_task", queue="celery")
def backfill_geom_task() -> dict[str, int]:
"""Разовый бэкфилл geom для уже загруженных строк без geom.
Запускать вручную через admin или celery CLI:
celery -A app.workers.celery_app call tasks.ekburg_permits_sync.backfill_geom_task
"""
return backfill_geom()

View file

@ -23,6 +23,7 @@ from app.services.scrapers.ekburg_permits import (
_to_date,
_to_float,
_to_str,
msk66_to_wgs84,
)
# ── helpers ───────────────────────────────────────────────────────────────────
@ -263,6 +264,11 @@ class TestParseXlsx:
assert row.total_area_sqm == pytest.approx(12345.67)
assert row.living_area_sqm == pytest.approx(8000.0)
assert row.source_year == 2024
# Тестовые координаты X=1530000, Y=380000 — должны попасть в ЕКБ-регион
assert row.geocoded_lat is not None
assert row.geocoded_lon is not None
assert 56.3 < row.geocoded_lat < 57.5
assert 59.5 < row.geocoded_lon < 61.5
def test_parse_rve_yields_permit_row_with_rve_fields(self) -> None:
content = _make_xlsx_bytes({"реестр разрешений на ввод": _rve_sheet_rows()})
@ -353,3 +359,163 @@ class TestSslConfiguration:
assert client._client._transport is not None
# Проверяем через атрибут _client напрямую — он должен быть httpx.Client
assert isinstance(client._client, httpx.Client)
# ── unit tests: msk66_to_wgs84 ───────────────────────────────────────────────
# Параметры допуска для cad_parcels cross-check (≈50м по Хаверсину)
# На широте ЕКБ (~56.9°): 1° широты ≈ 111320м, 1° долготы ≈ 60730м
# 50м ÷ 111320 ≈ 0.00045° по широте; 50м ÷ 60730 ≈ 0.00082° по долготе
_TOL_LAT_DEG = 0.00045 # ~50м по широте
_TOL_LON_DEG = 0.00082 # ~50м по долготе
class TestMsk66ToWgs84:
"""Тесты трансформации МСК-66 → WGS84.
Реальные данные из ekburg_construction_permits (сверены с cad_parcels ground-truth).
Калиброванный proj4: x_0=1496944, y_0=-5911077
avg_err14м, median7м, p9558м на 738 точках с кадастровым номером.
Точки для cross-check получены из cad_parcels.geom (ST_Centroid) и проверены SQL:
X=1534814.7997, Y=394813.2001 lon60.619637, lat56.871948
(cad 66:41:0110016:1046, погрешность <50м от центроида участка)
X=1544026.7997, Y=381585.4401 lon60.768255, lat56.752330
(cad 66:41:0503018:248, аэропорт Кольцово)
X=1526550.1397, Y=387129.3801 lon60.483224, lat56.803547
(cad 66:41:0306057:89, западная часть ЕКБ)
X=1531655.2797, Y=398430.2401 lon60.568320, lat56.904674
(cad 66:41:0106051:10, ул. Республиканская 1а; err_m=12)
"""
def test_kosmonavtov_cad_crosscheck(self) -> None:
"""X=1534814.7997/Y=394813.2001 — сверка с centroid cad_parcels 66:41:0110016:1046.
Ожидаемые координаты получены из PostGIS: lon60.619637, lat56.871948.
Допуск ±0.00082° по долготе (~50м), ±0.00045° по широте (~50м).
"""
result = msk66_to_wgs84("1534814.7997", "394813.2001")
assert result is not None
lon, lat = result
assert (
abs(lon - 60.619637) < _TOL_LON_DEG
), f"lon={lon:.6f} далеко от эталонного 60.619637 (Δ={abs(lon - 60.619637):.6f}°)"
assert (
abs(lat - 56.871948) < _TOL_LAT_DEG
), f"lat={lat:.6f} далеко от эталонного 56.871948 (Δ={abs(lat - 56.871948):.6f}°)"
def test_koltsovo_cad_crosscheck(self) -> None:
"""X=1544026.7997/Y=381585.4401 — сверка с centroid cad_parcels 66:41:0503018:248.
Ожидаемые координаты: lon60.768255, lat56.752330 (аэропорт Кольцово).
"""
result = msk66_to_wgs84("1544026.7997", "381585.4401")
assert result is not None
lon, lat = result
assert (
abs(lon - 60.768255) < _TOL_LON_DEG
), f"lon={lon:.6f} далеко от эталонного 60.768255 (Δ={abs(lon - 60.768255):.6f}°)"
assert (
abs(lat - 56.752330) < _TOL_LAT_DEG
), f"lat={lat:.6f} далеко от эталонного 56.752330 (Δ={abs(lat - 56.752330):.6f}°)"
def test_smorodinovaya_cad_crosscheck(self) -> None:
"""X=1526550.1397/Y=387129.3801 — сверка с centroid cad_parcels 66:41:0306057:89.
Ожидаемые координаты: lon60.483224, lat56.803547 (западная часть ЕКБ).
"""
result = msk66_to_wgs84("1526550.1397", "387129.3801")
assert result is not None
lon, lat = result
assert (
abs(lon - 60.483224) < _TOL_LON_DEG
), f"lon={lon:.6f} далеко от эталонного 60.483224 (Δ={abs(lon - 60.483224):.6f}°)"
assert (
abs(lat - 56.803547) < _TOL_LAT_DEG
), f"lat={lat:.6f} далеко от эталонного 56.803547 (Δ={abs(lat - 56.803547):.6f}°)"
def test_respublikanskaya_cad_crosscheck(self) -> None:
"""X=1531655.2797/Y=398430.2401 — ул. Республиканская 1а, cad 66:41:0106051:10.
Ожидаемые координаты: lon60.568320, lat56.904674 (err_m=12 по cad_parcels).
Этот тест ловит регрессию easting-сдвига: при x_0=1500000 computed был бы ~60.519°.
"""
result = msk66_to_wgs84("1531655.2797", "398430.2401")
assert result is not None
lon, lat = result
assert (
abs(lon - 60.568320) < _TOL_LON_DEG
), f"lon={lon:.6f} далеко от эталонного 60.568320 (Δ={abs(lon - 60.568320):.6f}°)"
assert (
abs(lat - 56.904674) < _TOL_LAT_DEG
), f"lat={lat:.6f} далеко от эталонного 56.904674 (Δ={abs(lat - 56.904674):.6f}°)"
def test_kosmonavtov_in_ekb(self) -> None:
"""Пр. Космонавтов — попадает в ЕКБ-регион."""
result = msk66_to_wgs84("1534814.7997", "394813.2001")
assert result is not None
lon, lat = result
assert 56.5 < lat < 57.2, f"lat={lat} вне ЕКБ"
assert 60.0 < lon < 61.0, f"lon={lon} вне ЕКБ"
def test_koltsovo_in_ekb(self) -> None:
"""Район аэропорта Кольцово — попадает в ЕКБ-регион."""
result = msk66_to_wgs84("1544026.7997", "381585.4401")
assert result is not None
lon, lat = result
assert 56.5 < lat < 57.2, f"lat={lat} вне ЕКБ"
assert 60.3 < lon < 61.2, f"lon={lon} вне ЕКБ"
def test_smorodinovaya_in_ekb(self) -> None:
"""Ул. Смородиновая, западная часть ЕКБ."""
result = msk66_to_wgs84("1526550.1397", "387129.3801")
assert result is not None
lon, lat = result
assert 56.5 < lat < 57.2, f"lat={lat} вне ЕКБ"
assert 59.5 < lon < 61.0, f"lon={lon} вне ЕКБ"
def test_multi_point_cell_takes_first(self) -> None:
"""Если ячейка содержит несколько точек через пробел — берём первую."""
# Реальный пример: "1528834.9797 1528819.0197"
result = msk66_to_wgs84("1528834.9797 1528819.0197", "389190.2401 389136.0601")
assert result is not None
lon, lat = result
assert 56.3 < lat < 57.5
assert 59.5 < lon < 61.5
def test_comma_decimal_separator(self) -> None:
"""Запятая как десятичный разделитель (встречается в xlsx)."""
# "1534655,7197" — реальный случай
result = msk66_to_wgs84("1534655,7197", "387724.5201")
assert result is not None
lon, lat = result
assert 56.3 < lat < 57.5
assert 59.5 < lon < 61.5
def test_none_input_returns_none(self) -> None:
assert msk66_to_wgs84(None, None) is None
def test_empty_string_returns_none(self) -> None:
assert msk66_to_wgs84("", "") is None
def test_non_numeric_returns_none(self) -> None:
assert msk66_to_wgs84("нет данных", "нет данных") is None
def test_out_of_range_x_returns_none(self) -> None:
"""X вне диапазона МСК-66 (например, случайное число) → None."""
# X=500000 — это не МСК-66 (false_easting ~1500000, диапазон 1450000-1600000)
result = msk66_to_wgs84("500000.0", "380000.0")
assert result is None
def test_result_is_lon_lat_order(self) -> None:
"""Результат — (lon, lat), а не (lat, lon).
ЕКБ: lon ~60.5E, lat ~56.9N lon > lat числово.
"""
result = msk66_to_wgs84("1534814.7997", "394813.2001")
assert result is not None
lon, lat = result
# Проверяем что первый элемент — долгота (60.x), второй — широта (56.x)
assert 56.0 < lat < 58.0, f"ожидали широту, получили {lat}"
assert 59.0 < lon < 62.0, f"ожидали долготу, получили {lon}"
assert lon > lat # lon ~60.6 > lat ~56.9