Bagaimana cara saya mengelompokkan titik dekat dengan posisi GPS?

10

Saya orang IT, jadi saya tidak tahu terlalu banyak tentang proyeksi dan sebagainya, saya harap Anda dapat membantu saya.

Saya telah membuat aplikasi untuk Android yang mengumpulkan posisi GPS, jadi saya memiliki lintang dan bujur pada waktu tertentu. Saya ingin menyimpan bersama elemen-elemen yang berdekatan satu sama lain, dalam kelompok area medan dengan ukuran fisik yang sama; masalahnya adalah saya tidak tahu poinnya sebelumnya dan mereka bisa datang dari posisi mana pun di dunia .

Ide pertama saya (untuk menjelaskan sedikit lebih banyak masalah) adalah menggunakan desimal lintang dan bujur untuk pengelompokan; misalnya, satu grup akan menjadi posisi dengan garis lintang antara 35.123 dan 35.124, dan garis bujur antara 60.101 dan 60.102. Jadi jika saya mendapatkan posisi seperti lat = 35.1235647 dan lon = 60.1012254598, poin ini akan masuk ke grup itu.

Solusi ini akan OK untuk representasi 2D kartesian, karena saya akan memiliki luas 0,001 unit lebar dan tinggi; Namun, karena ukuran 1degree dari bujur berbeda pada lintang yang berbeda, saya tidak dapat menggunakan pendekatan ini.

Ada ide?

Kuu
sumber
Mengapa Anda tidak dapat menyimpan posisi seperti pada, dan kemudian melakukan pemrosesan nanti? Juga di khatulistiwa, 1 derajat garis bujur adalah sekitar 111 km, jadi 0,001 derajat akan sedikit lebih dari 1 Km. Apakah Anda benar-benar ingin tempat sampah Anda begitu besar?
Devdatta Tengshe
Tingkat 0,001 hanyalah contoh dari ide saya. Tentu saja saya harus menyesuaikannya dengan persyaratan. Saya tidak dapat melakukan pemrosesan pos karena ini adalah aplikasi waktu nyata, dan saya tidak dapat mengatakan kepada pengguna saya "tunggu sampai besok, karena saya harus mengelompokkan poin" Bagaimanapun juga, terima kasih atas ide-idenya;)
Kuu

Jawaban:

6

Salah satu cara cepat dan kotor menggunakan pembagian bola rekursif . Dimulai dengan triangulasi permukaan bumi, secara rekursi membelah setiap segitiga dari sebuah sudut melintang ke tengah sisi terpanjangnya. (Idealnya Anda akan membagi segitiga menjadi dua bagian dengan diameter yang sama atau bagian dengan luas yang sama, tetapi karena itu melibatkan beberapa perhitungan fiddly, saya hanya membagi sisi menjadi dua: ini menyebabkan berbagai segitiga pada akhirnya sedikit berbeda ukurannya, tetapi yang sepertinya tidak penting untuk aplikasi ini.)

Tentu saja Anda akan mempertahankan pembagian ini dalam struktur data yang memungkinkan dengan cepat mengidentifikasi segitiga di mana titik arbitrer berada. Pohon biner (berdasarkan pada panggilan rekursif) berfungsi dengan baik: setiap kali sebuah segitiga terbelah, pohon itu terbelah pada simpul segitiga itu. Data mengenai bidang yang terbelah dipertahankan, sehingga Anda dapat dengan cepat menentukan di sisi mana saja titik arbitrer berada: yang akan menentukan apakah akan melakukan perjalanan ke kiri atau ke kanan melewati pohon.

(Apakah saya mengatakan pemisahan "bidang"? Ya - jika memodelkan permukaan bumi sebagai bola dan menggunakan koordinat geosentris (x, y, z), maka sebagian besar perhitungan kami dilakukan dalam tiga dimensi, di mana sisi segitiga adalah potongan persimpangan bola dengan pesawat melalui asalnya. Ini membuat perhitungan cepat dan mudah.)


Saya akan menggambarkan dengan menunjukkan prosedur pada satu oktan bola; tujuh oktan lainnya diproses dengan cara yang sama. Oktan seperti itu adalah segitiga 90-90-90. Dalam grafik saya, saya akan menggambar segitiga Euclidean yang mencakup sudut yang sama: mereka tidak terlihat sangat bagus sampai mereka menjadi kecil, tetapi mereka dapat dengan mudah dan cepat digambar. Berikut adalah segitiga Euclidean yang sesuai dengan oktan: ini adalah awal dari prosedur.

Gambar 1

Karena semua pihak memiliki panjang yang sama, satu dipilih secara acak sebagai "terpanjang" dan dibagi:

Gambar 2

Ulangi ini untuk masing-masing segitiga baru:

Gambar 3

Setelah n langkah kita akan memiliki 2 ^ n segitiga. Berikut adalah situasi setelah 10 langkah, menunjukkan 1024 segitiga dalam oktan (dan 8192 pada keseluruhan bola):

Gambar 4

Sebagai ilustrasi lebih lanjut, saya menghasilkan titik acak dalam oktan ini dan melakukan perjalanan pohon pembagian sampai sisi terpanjang segitiga mencapai kurang dari 0,05 radian. Segitiga (cartesian) ditunjukkan dengan titik probe berwarna merah.

Gambar 5

Secara kebetulan, untuk mempersempit lokasi titik ke satu derajat garis lintang (kurang-lebih), Anda akan mencatat bahwa ini sekitar 1/60 radian dan mencakup sekitar (1/60) ^ 2 / (Pi / 2) = 1/6000 dari permukaan total. Karena setiap subdivisi kira-kira membagi ukuran segitiga, sekitar 13 hingga 14 subdivisi dari oktan akan melakukan trik. Itu tidak banyak perhitungan - seperti yang akan kita lihat di bawah ini - membuatnya efisien untuk tidak menyimpan pohon sama sekali, tetapi untuk melakukan pembagian dengan cepat. Pada awalnya Anda akan mencatat di mana oktan titik terletak - yang ditentukan oleh tanda-tanda tiga koordinatnya, yang dapat dicatat sebagai angka biner tiga digit - dan pada setiap langkah Anda ingin mengingat apakah titik terletak di sebelah kiri (0) atau kanan (1) dari segitiga. Itu memberi angka biner 14 digit lagi. Anda dapat menggunakan kode ini untuk mengelompokkan poin arbitrer.

(Secara umum, ketika dua kode ditutup sebagai angka biner aktual, titik-titik yang sesuai dekat; namun, poin masih dapat ditutup dan memiliki kode yang sangat berbeda. Pertimbangkan dua titik terpisah satu meter terpisah dengan Khatulistiwa, misalnya: kode mereka harus berbeda bahkan sebelum titik biner, karena mereka berada di oktan yang berbeda. Hal semacam ini tidak dapat dihindari dengan pembagian ruang yang tetap.)


Saya menggunakan Mathematica 8 untuk mengimplementasikan ini: Anda dapat mengambil apa adanya atau sebagai pseudocode untuk implementasi di lingkungan pemrograman favorit Anda.

Tentukan pada sisi mana bidang 0-ab titik p terletak pada:

side[p_, {a_, b_}] := If[Det[{p, a, b}] >=  0, left, right];

Perbaiki segitiga abc berdasarkan titik p.

refine[p_, {a_, b_, c_}] := Block[{sides, x, y, z, m},
  sides = Norm /@ {b - c, c - a, a - b} // N;
  {x, y, z} = RotateLeft[{a, b, c}, First[Position[sides, Max[sides]]] - 1];
  m = Normalize[Mean[{y, z}]];
  If[side[p, {x, m}] === right, {y, m, x}, {x, m, z}] 
  ]

Gambar terakhir diambil dengan menampilkan oktan dan, di atas itu, dengan menjadikan daftar berikut ini sebagai satu set poligon:

p = Normalize@RandomReal[NormalDistribution[0, 1], 3]        (* Random point *)
{a, b, c} = IdentityMatrix[3] . DiagonalMatrix[Sign[p]] // N (* First octant *)
NestWhileList[refine[p, #] &, {a, b, c}, Norm[#[[1]] - #[[2]]] >= 0.05 &, 1, 16]

NestWhileListberulang kali menerapkan operasi ( refine) ketika suatu kondisi berkaitan (segitiga besar) atau sampai jumlah operasi maksimum telah tercapai (16).

Untuk menampilkan triangulasi penuh oktan, saya mulai dengan oktan pertama dan mengulangi penyempurnaan sepuluh kali. Ini dimulai dengan sedikit modifikasi refine:

split[{a_, b_, c_}] := Module[{sides, x, y, z, m},
  sides = Norm /@ {b - c, c - a, a - b} // N;
  {x, y, z} = RotateLeft[{a, b, c}, First[Position[sides, Max[sides]]] - 1];
  m = Normalize[Mean[{y, z}]];
  {{y, m, x}, {x, m, z}}
  ]

Perbedaannya adalah bahwa splitmengembalikan kedua bagian dari segitiga inputnya daripada segitiga di mana titik tertentu berada. Triangulasi penuh diperoleh dengan iterasi ini:

triangles = NestList[Flatten[split /@ #, 1] &, {IdentityMatrix[3] // N}, 10];

Untuk memeriksanya, saya menghitung ukuran setiap segitiga dan melihat kisarannya. ("Ukuran" ini proporsional dengan sosok berbentuk piramida yang digantikan oleh setiap segitiga dan pusat bola; untuk segitiga kecil seperti ini, ukuran ini pada dasarnya sebanding dengan area bulatnya.)

Through[{Min, Max}[Map[Round[Det[#], 0.00001] &, triangles[[10]] // N, {1}]]]

{0,00523, 0,00739}

Dengan demikian ukuran bervariasi naik atau turun sekitar 25% dari rata-rata mereka: yang tampaknya masuk akal untuk mencapai cara yang hampir seragam untuk mengelompokkan poin.

Dalam pemindaian kode ini, Anda tidak akan melihat trigonometri : satu-satunya tempat yang diperlukan, jika sama sekali, adalah konversi bolak-balik antara koordinat bola dan koordinat kartesius. Kode juga tidak memproyeksikan permukaan bumi ke peta apa pun, sehingga menghindari distorsi yang menyertainya. Kalau tidak, itu hanya menggunakan rata-rata ( Mean), Teorema Pythagoras ( Norm) dan penentu 3 oleh 3 ( Det) untuk melakukan semua pekerjaan. (Ada beberapa perintah manipulasi daftar sederhana seperti RotateLeftdan Flatten, juga, bersama dengan pencarian sisi terpanjang dari setiap segitiga.)

whuber
sumber
1

Ini adalah yang sulit, karena proyeksi memperkenalkan distorsi ketika beralih dari sistem koordinat geografis WGS84 3D ke proyeksi 2D datar. Pada skala global, Anda pasti akan mengalami distorsi di suatu tempat dalam sistem.

Saya pikir taruhan terbaik Anda adalah memproyeksikan ke proyeksi Universal Transverse Mercator . Sejauh yang saya tahu ini adalah yang terdekat dengan proyeksi global dengan paling sedikit distorsi.

Jika Anda dapat melonggarkan persyaratan bahwa kelompok harus didefinisikan di area dengan ukuran yang sama persis serta persyaratan pemrosesan waktu nyata, ada algoritma pengelompokan seperti DBSCAN , dan keluarga turunannya, yang dapat membantu menghasilkan pengelompokan.

Sean Barbeau
sumber
1
UTM bukan "proyeksi global": demonstrasi adalah bahwa hampir setiap pasangan koordinat yang valid, seperti (500000, 5000000), sesuai dengan setidaknya 120 titik terpisah yang berbeda dalam sistem UTM. Dan, sayangnya, algoritma pengelompokan tidak memenuhi kebutuhan OP untuk dapat mengelompokkan poin secara real time hanya berdasarkan lokasi mereka (dan bukan pada kedekatannya dengan titik lain).
whuber
@whuber re: "proyeksi global" - benar. Itu sebabnya saya mengatakan "yang terdekat dengan proyeksi global". Jika Anda mengetahui sistem proyeksi yang lebih baik yang lebih sesuai silakan tinggalkan di komentar dan saya akan mengedit jawaban saya. Selain itu, OP tidak memiliki persyaratan waktu nyata di pos awal. Saya akan mengedit jawaban saya untuk memperhitungkan ini.
Sean Barbeau
Sean, (1) Solusi saya untuk masalah proyeksi global adalah dengan tidak menggunakannya. Ada ada ada proyeksi global tanpa singularitas. (2) Benar, klarifikasi waktu-nyata muncul dalam komentar. Teks yang mengikuti "ide pertama saya" memang melakukan pekerjaan dengan baik, meskipun, dalam menekankan bahwa masalah ini adalah memartisi permukaan bumi daripada mengelompokkan satu set lokasi. Itulah poin yang saya (tidak terlalu efektif) coba sampaikan dalam komentar saya sebelumnya kepada Anda.
whuber