Descubra o ponto que fica entre duas linhas paralelas

Estou enfrentando um problema no ArcGIS. Eu trabalho em um banco de dados de navegação. Em nosso banco de dados, as ruas de pista única são representadas por uma única linha, enquanto uma rua de pista múltipla (rua com divisor no centro) é representada por duas linhas paralelas (linhas vermelhas na imagem).

Eu tenho um shapefile de ponto com alguns pontos caindo dentro da rua com várias pistas e outros fora.

Eu quero criar um script do ArcPy que encontre os pontos que se enquadram nas ruas de várias faixas. ou seja, entre essas linhas paralelas (marcadas na figura).

Eu não sei como conseguir isso, alguém pode me ajudar?

Eu fiz alguns exercícios e descobri que a criação de buffer em um lado da linha pode criar dentro do polígono de várias faixas (mostrado na figura).

mas agora o problema é que o polígono está realmente cruzando a linha (ou seja, sobrepondo-se ao limite de várias faixas). então ele vai pegar pontos desnecessários. existe alguma maneira de alinhar esse polígono com a linha da rua?

Nota: integrar não funcionará aqui, porque também move a linha da rua. Eu preciso apenas alinhar o polígono ao longo da linha da rua.

Eu tentaria abaixo do algoritmo arcpy (mesmo manual!)

1. Encontre a largura adequada das ruas de duas faixas - aqui você pode precisar agrupar ruas com a mesma largura e seguir o procedimento abaixo para cada cluster.
2. Crie um buffer de ambas as linhas para as duas direções (direita e esquerda) com essa largura (ou um pouco menos que isso - para garantir a área da estrada).
3. Execute a ferramenta Intersecção para obter a região Sobreposta.
4. Execute Selecionar por local para selecionar pontos que se enquadram dentro deste polígono.

Eu diria que este é um exercício geométrico.

PSEUDO-CÓDIGO:

• Para cada ponto (ponto preto), encontre a estrada mais próxima e localize a projeção do ponto nessa estrada (ponto vermelho).
• Desenhe uma linha curta (tracejada) na direção oposta, começando no ponto preto
• Descubra se há interseção entre a linha curta e o mesmo nome da estrada, estrela azul. Se houver, o ponto preto é o que procuramos.

Como se pode ver, existem casos especiais - pontos pretos circulados:

1. Estrada muito sinuosa de 1 linha. Isso pode ser eliminado por a) trabalhando apenas com estradas de 2 linhas ou b) certificando-se de que os FIDs das estradas que cruzam ponto vermelho e estrela sejam diferentes. No entanto, se a estrada dobrada tiver uma junção com outra estrada de 1 linha, isso pode não funcionar.
2. O ponto preto está localizado na extensão da estrada de 1 linha exatamente perpendicular. Nesse caso, existe a possibilidade de uma estrada de 1 faixa ser escolhida como o vizinho mais próximo.
3. Ponto preto fica na linha.

Todos os casos acima são muito improváveis, no entanto, parece que a opção mais segura é trabalhar apenas com duas estradas de linha, ou seja, exportá-las para uma classe de recurso separada. O caso 3 é engraçado, vamos deixar ao acaso, porque a menor distância da linha nunca é zero, portanto, a direção 'oposta' do raio que conecta 2 pontos pode ser encontrada.

Implementação do Python:

import arcpy, traceback, os, sys
from arcpy import env
env.overwriteoutput=True

# things to change ---------
maxD=30
mxd = arcpy.mapping.MapDocument("CURRENT")
pointLR = arcpy.mapping.ListLayers(mxd,"NODES")[0]
lineLR = arcpy.mapping.ListLayers(mxd,"LINKS")[0]
sjOneToMany=r'D:\scratch\sj2.shp'
RDNAME='street'
# -------------------------
dDest=arcpy.Describe(lineLR)
SR=dDest.spatialReference

try:
    def showPyMessage():
        arcpy.AddMessage(str(time.ctime()) + " - " + message)
    g = arcpy.Geometry()
    geometryList=arcpy.CopyFeatures_management(pointLR,g)
    n=len(geometryList)
    endPoint=arcpy.Point()

    arcpy.SpatialJoin_analysis(pointLR, lineLR,sjOneToMany,"JOIN_ONE_TO_MANY","KEEP_COMMON","","WITHIN_A_DISTANCE",maxD)
    initFidList=(-1,)
    for fid in range(n):
        query='"TARGET_FID" = %s' %str(fid)
        nearTable=arcpy.da.TableToNumPyArray(sjOneToMany,("TARGET_FID","JOIN_FID"),query)
        if len(nearTable)<2:continue
        fidLines=[int(row[1]) for row in nearTable]
        query='"FID" in %s' %str(tuple(fidLines))
        listOfLines={}
        blackPoint=geometryList[fid]
        with arcpy.da.SearchCursor(lineLR,("FID", "Shape@","STREET"),query) as rows:
            dMin=100000
            for row in rows:
                shp=row[1];dCur=blackPoint.distanceTo(shp)
                listOfLines[row[0]]=row[-2:]
                if dCur<dMin:
                    fidNear,lineNear, roadNear=row
                    dMin=dCur
            chainage=lineNear.measureOnLine(blackPoint)
            redPoint=lineNear.positionAlongLine (chainage).firstPoint
            smallD=blackPoint.distanceTo(redPoint)
            fp=blackPoint.firstPoint
            dX=(redPoint.X-fp.X)*(maxD-smallD)/smallD
            dY=(redPoint.Y-fp.Y)*(maxD-smallD)/smallD
            endPoint.X=fp.X-dX;endPoint.Y=fp.Y-dY
            dashLine=arcpy.Polyline(arcpy.Array([fp,endPoint]),SR)

            for n in listOfLines:
                if n==fidNear:continue
                line, road=listOfLines[n]
                if road!=roadNear:continue
                blueStars=dashLine.intersect(line,1)
                if blueStars.partCount==0:continue
                initFidList+=(fid,); break
    query='"FID" in %s' %str(initFidList)
    arcpy.SelectLayerByAttribute_management(pointLR, "NEW_SELECTION", query)
    arcpy.AddMessage ('\n %i point(s) found' %(len(initFidList)-1))
except:
    message = "\n*** PYTHON ERRORS *** "; showPyMessage()
    message = "Python Traceback Info: " + traceback.format_tb(sys.exc_info()[2])[0]; showPyMessage()
    message = "Python Error Info: " +  str(sys.exc_type)+ ": " + str(sys.exc_value) + "\n"; showPyMessage()            

Existe outra solução possível, talvez mais elegante. Envolve triangulação. Informe-me se for de seu interesse e atualizarei minha resposta

Como as ruas são paralelas, presumi que elas foram criadas com a Copy Parallelferramenta na barra de ferramentas Editar, fazendo com que o par de linhas tenha a mesma direção. Podemos então iterar sobre as coordenadas da primeira linha e adicioná-las a um polígono e, em seguida, iterar no verso da segunda linha. Definitivamente, há uma maneira melhor de abordar os pares de linhas de agarrar; a abordagem OID funciona, mas não é muito bonita.

import collections
import arcpy

FC = "fc"
points = "points"
pgons = "pgons"
arcpy.env.overwriteOutput = True

def buildpoly(oid_coords):
    #create ddict of the form OID:<x1y1, x2y2, ..., xn-1yn-1, xnyn>
    ddict = collections.defaultdict(list)    
    for k,v in oid_coords:
        ddict[k].append(v)

    line1,line2 = ddict.keys()    

    #Assume that the parallel lines have same direction, so reverse the second
    arr = arcpy.Array()
    arr.extend(arcpy.Point(*pt) for pt in ddict[line1])    
    arr.extend(arcpy.Point(*pt) for pt in ddict[line2][::-1])

    return arcpy.Polygon(arr)

#id is an integer field that pairs parallel lines together
unique = list(set(t[0] for t in arcpy.da.SearchCursor(FC, "id")))
polygons = []
for uni in unique:
    polygons.append(buildpoly([r for r in row] for row in arcpy.da.SearchCursor(FC,
                                                                                ["OID@", "SHAPE@XY"],
                                                                                "id={}".format(uni),
                                                                                explode_to_points=True)))


arcpy.CopyFeatures_management(polygons, pgons)

A partir daí, é uma ligação para Intersect / Select Layer por local / o que você tem. Observe que o Spolígono em forma não é perfeito, pois eu o desenhei à mão livre e há alguns arcos que explode_to_pointsnão são adequados. Basta executar Densifyou equivalente.