Identifizierung aufeinanderfolgender Punkte innerhalb eines bestimmten Puffers

Ich habe eine Punktedatei, stündliche Umzüge eines Tieres, und ich möchte in der Lage sein, einen Puffer um jeden Punkt zu legen und die Anzahl der nachfolgenden Punkte zu berechnen, die sich innerhalb der Pufferzone befinden. Ich suche nach einer Methode, die entlang der Punktdatei funktioniert, wie ein sich bewegendes Fenster, das nur die nachfolgenden Punkte zählt, die sich innerhalb der Pufferzone befinden.

Zum Beispiel platziere ich an Punkt 10 einen Puffer von 1500 m und möchte wissen, ob Punkt 11 innerhalb des Puffers liegt und wenn ja, ob Punkt 12 innerhalb des Puffers liegt und so weiter. Ich möchte nicht wissen, ob sich Punkt 52 innerhalb der Pufferzone befindet, es sei denn, alle vorherigen Punkte befanden sich innerhalb des Puffers. Ich möchte auch nicht wissen, ob sich die Punkte 9 oder 8 usw. im Puffer befinden.

Ich habe eine zusätzliche Toolbox namens "Moving Window Point Analysis" gefunden und ausprobiert, die als Moving Window für die Punktdatei fungiert. Dies funktioniert gut, aber sehr langsam und schließt alle Punkte ein, die sich innerhalb der Pufferzone befinden, auch wenn es sich nicht um aufeinanderfolgende Punkte handelt. Ich kann keinen Weg finden, um nur aufeinanderfolgende Punkte zu betrachten.

Ich möchte eine Methode, die eine Ausgabetabelle bereitstellt, da ich viele Datenpunkte auf diese Weise betrachten kann.

Ich verwende ArcGIS 10. Jede Hilfe, die jemand leisten kann, wäre sehr dankbar.

Bei einer Liste von Punktpositionen (vorzugsweise in projizierten Koordinaten, damit Entfernungen leicht berechnet werden können) kann dieses Problem mit fünf einfacheren Operationen gelöst werden :

1. Berechnen Sie Punkt-Punkt-Abstände.

2. Identifizieren Sie für jeden Punkt i, i = 1, 2, ... die Indizes dieser Punkte in Abständen, die kleiner als der Pufferradius sind (z. B. 1500).

3. Beschränken Sie diese Indizes auf i oder höher.

4. Behalten Sie nur die erste aufeinanderfolgende Gruppe von Indizes bei, die keine Unterbrechung aufweisen.

5. Geben Sie die Anzahl dieser Gruppe aus.

In entspricht Rjede dieser Operationen einer Operation. Um gelten für jeden Punkt dieser Sequenz ist es bequem, die meisten der Arbeit innerhalb einer Funktion kapseln wir definieren , also:

#
# forward(j, xy, r) counts how many contiguous rows in array xy, starting at index j,
#                   are within (Euclidean) distance r of the jth row of xy.
#
forward <- function(j, xy, r) {
  # Steps 1 and 2: compute an array of indexes of points within distance r of point j.
  i <- which(apply(xy, 1, function(x){sum((x-xy[j,])^2) <= r^2}))
  # Step 3: select only the indexes at or after j.
  i <- i[i >= j]
  # Steps 4 and 5: retain only the first consecutive group and count it.
  length(which(i <= (1:length(i) + j)))
}

(Eine effizientere Version dieser Funktion finden Sie weiter unten.)

Ich habe diese Funktion flexibel genug gemacht, um verschiedene Punktlisten ( xy) und Pufferabstände ( r) als Parameter zu akzeptieren .

Normalerweise lesen Sie eine Datei mit Punktpositionen (und sortieren Sie sie gegebenenfalls nach Zeit). Um dies in Aktion zu zeigen, werden wir nur einige Beispieldaten zufällig generieren :

# Create sample data
n<-16                                     # Number of points
set.seed(17)                              # For reproducibility
xy <- matrix(rnorm(2*n) + 1:n, n, 2) * 300
#
# Display the track.
plot(xy, xlab="x", ylab="y")
lines(xy, col="Gray")

Ihr typischer Abstand beträgt 300 * Sqrt (2) = ungefähr 500. Wir führen die Berechnung durch, indem wir diese Funktion auf die Punkte im Array anwendenxy (und dann die Ergebnisse wieder anheften xy, da dies ein praktisches Format für den Export in ein GIS wäre ):

radius <- 1500
z <- sapply(1:n, function(u){forward(u,xy,radius)})
result <- cbind(xy, z)                              # List of points, counts

Anschließend analysieren Sie das resultArray weiter, Rindem Sie es entweder in eine Datei schreiben oder in eine andere Software importieren. Hier ist das Ergebnis für die Beispieldaten :

                        z
  [1,]   -4.502615  551.5413 4
  [2,]  576.108979  647.8110 3
  [3,]  830.103893 1087.7863 4
  [4,]  954.819620 1390.0754 3
...
 [15,] 4977.361529 4146.7291 2
 [16,] 4783.446283 4511.9500 1

(Denken Sie daran , dass die Zählungen sind die Punkte , an denen sie basieren, so dass jede Zählung 1 oder größer sein müssen.)

Wenn Sie viele tausend Punkte haben, ist diese Methode zu ineffizient : Sie berechnet viel zu viele unnötige Punkt-zu-Punkt-Abstände. Da wir jedoch die Arbeit in die forwardFunktion eingekapselt haben , lässt sich die Ineffizienz leicht beheben. Hier ist eine Version, die besser funktioniert, wenn mehr als ein paar hundert Punkte beteiligt sind:

forward <- function(j, xy, r) {
   n <- dim(xy)[1]     # Limit the search to the number of points in xy
   r2 <- r^2           # Pre-compute the squared distance threshold
   z <- xy[j,]         # Pre-fetch the base point coordinates
   i <- j+1            # Initialize an index into xy (just past point j)

   # Advance i while point i remains within distance r of point j.
   while(i <= n && sum((xy[i,]-z)^2) <= r2) i <- i+1

   # Return the count (including point j).
   i-j
}

Um dies zu testen, habe ich wie zuvor zufällige Punkte erstellt, aber zwei Parameter variiert: n(die Anzahl der Punkte) und ihre Standardabweichung (oben als 300 fest codiert). Die Standardabweichung bestimmt die durchschnittliche Anzahl von Punkten in jedem Puffer ("Durchschnitt" in der folgenden Tabelle): Je mehr Punkte vorhanden sind, desto länger dauert die Ausführung dieses Algorithmus. (Bei komplexeren Algorithmen hängt die Laufzeit nicht so stark davon ab, wie viele Punkte sich in jedem Puffer befinden.) Hier einige Zeitpunkte:

 Time (sec)    n    SD  Average  Distances checked per minute
 1.30       10^3     3  291      13.4 million
 1.72       10^4    30   35.7    12.5
 2.50       10^5   300    3.79    9.1
16.4        10^6  3000    1.04    3.8

Ich denke, Ihre beste Wette ist es, eine kleine Routine mit ArcPy zu schreiben. Ich würde so etwas wie diesen Pseudocode erstellen:

Select all points sorted by location-id    
Iterate for each point:
    Select points by location using a distance (no need to create buffer geometry)
    Sort points by location-id
    Set a variable to the value of your reference id
    consecutive-counter = 0
        Iterate your selection:
            Is the location-id of the first (or next) point equal to variable + 1?
            if 'yes' increment consecutive-counter by 1
            Repeat until 'no'

Ich bin nicht sicher, was Sie mit den Informationen tun möchten, aber ich nehme an, Sie könnten ein Feld in Ihrer Tabelle erstellen und es mit der Anzahl aufeinanderfolgender Punkte aktualisieren (wenn ja, fügen Sie das Feld zuerst hinzu).

Ich würde empfehlen, Feature-Layer zu erstellen (wie eine Datenbanktabellenansicht, jedoch für Features in Arc). Machen Sie zwei aus den Originaldaten und öffnen Sie einen Aktualisierungscursor auf dem ersten Satz, der Ihre Gesamtsortierung angibt (da ESRI keine vollständigen SQL-Abfragen berücksichtigt). Verwenden Sie den zweiten, um nach Position auszuwählen, und öffnen Sie einen Suchcursor auf dem resultierenden Auswahlsatz.

[BEARBEITEN NACH James ANFRAGE] Mit Model Builder ausarbeiten. Wenn Sie Model Builder noch nicht verwendet haben, müssen Sie nur mit der rechten Maustaste in arcToolbox klicken. Wählen Sie 'Toolbox hinzufügen'. Klicken Sie mit der rechten Maustaste auf die neue Toolbox und klicken Sie auf "Neu-> Modell". Wenn Sie ein neues Modellfenster haben, ziehen Sie die benötigten Werkzeuge und Daten per Drag & Drop in das Fenster und verknüpfen Sie sie visuell (mit dem kleinen Pfeilwerkzeug). Wenn Sie so weit wie möglich gekommen sind (Sie können Ihre Cursor hier nicht hinzufügen), verwenden Sie die Option im Dateimenü des Model Builder, um nach Python zu exportieren. Das bringt Sie den größten Teil des Weges dorthin. Es ist automatisch generierter Code, wird also ein bisschen böse, aber funktional sein. Verwenden Sie dann die Links in meiner Antwort oben, um die Cursor zu verstehen und hinzuzufügen.

Wenn Sie Python noch nicht kennen, haben Sie keine Angst davor, Code zu schreiben! Python ist eine sehr einfache Skriptsprache, mit der Sie schnell Ergebnisse erzielen können. Esri hat auch eine Anleitung dazu.

Wenn Sie mit Ihrem Code nicht weiterkommen, veröffentlichen Sie ihn in diesem Forum und bitten Sie um Hilfe. Es gibt viele Leute hier, die helfen können. Eine Warnung: Stellen Sie sicher, dass Sie die richtige Hilfe von ESRI verwenden. Sie haben die Python-API zwischen den Versionen 9.x und 10 (arcgisscripting bzw. arcpy) massiv geändert. Wenn Sie also ArcGIS 9.x verwenden, finden Sie die entsprechenden Links zu meinen!

Sie können den Modellbauer in ArcGIS verwenden, um aufeinanderfolgende ID-Werte zu finden. Ich habe mein Modell als Python-Skript exportiert. Der Code generiert einen neuen shp mit aufeinanderfolgenden ID-Werten. !ICH WÜRDE! ist das Basis-ID-Feld. Sie müssen den Pfad, den Namen und den ID-Feldnamen von point2.shp entsprechend Ihrem Fall aktualisieren.

# Import arcpy module
import arcpy


# Local variables:
point2_shp = "C:\\temp\\point2.shp"
point2_shp__2_ = "C:\\temp\\point2.shp"
point2_shp__4_ = "C:\\temp\\point2.shp"
Freq_dbf__2_ = "C:\\temp\\Freq.dbf"
point2_shp__5_ = "C:\\temp\\point2.shp"
point2__2_ = "C:\\temp\\point2.shp"
point2__4_ = "C:\\temp\\point2.shp"
Freq_dbf = "C:\\temp\\Freq.dbf"
PointConsecutive_shp = "C:\\temp\\PointConsecutive.shp"

# Process: Add Field
arcpy.AddField_management(point2_shp, "AUTOID", "LONG", "", "", "", "", "NON_NULLABLE", "NON_REQUIRED", "")

# Process: Calculate Field
arcpy.CalculateField_management(point2__2_, "AUTOID", "autoIncrement()", "PYTHON_9.3", "rec=0\\ndef autoIncrement():\\n global rec\\n pStart = 1 #adjust start value, if req'd \\n pInterval = 1 #adjust interval value, if req'd\\n if (rec == 0): \\n  rec = pStart \\n else: \\n  rec = rec + pInterval \\n return rec\\n")

# Process: Add Field (2)
arcpy.AddField_management(point2_shp, "DIFF", "LONG", "", "", "", "", "NON_NULLABLE", "NON_REQUIRED", "")

# Process: Calculate Field (2)
arcpy.CalculateField_management(point2__4_, "DIFF", "!ID! - !AUTOID!", "PYTHON_9.3", "")

# Process: Frequency
arcpy.Frequency_analysis(point2_shp__2_, Freq_dbf, "DIFF", "")

# Process: Join Field
arcpy.JoinField_management(point2_shp__4_, "DIFF", Freq_dbf__2_, "DIFF", "")

# Process: Select
arcpy.Select_analysis(point2_shp__5_, PointConsecutive_shp, "\"FREQUENCY\" >1")