Can anyone help me what's the algorythm to calculate points of intersection whet is given one point p1 [lat,long] , bearing x (let say 137 degrees), arc center [lat,long] and radius? I woulde like to get p2[lat,long] and p3[lat,long].
For better imagination what I mean something like on this picture but on sphere with only one point and bearing:
In the Specification for the Wide Area Augmentation System (WAAS) on page 99 there are two formulae for calculating the "pierce point" of a GPS signal through the ionosphere. This may not be your application area, but may be able to adapt these to your application and they happen to work with longitude, latitude, azimuth, and elevation angle. Here's a picture from the document:
On the site http://paulbourke.net/geometry/circlesphere/ there is a section called "Intersection of a Line and a Sphere (or circle)" with more general equations for intersections. At the end of that section you will find some implementations in C, VB, Python, and LISP.
As a bonus (this is a limited time offer, you are one of a select few to receive this offer ;-) ) here is a short Python script for visualizing line-sphere intersections in Blender (adapted from this blog post. Note that is takes Earth-centered, Earth-fixed coordinates, not longitude/latitude/azimuth/elevation angles.
import bpy
EARTH_EQUATOR = 6378137.0 # meters
IONO_HEIGHT = 350000.0 # meters
# scale down the radius by 1E6 so that it fits in blender's world
radius = (EARTH_EQUATOR + IONO_HEIGHT) / 1E6
bpy.ops.mesh.primitive_uv_sphere_add(size=radius)
from mathutils import Vector
# weight
w = 10
def MakePolyLine(objname, curvename, cList):
curvedata = bpy.data.curves.new(name=curvename, type='CURVE')
curvedata.dimensions = '3D'
objectdata = bpy.data.objects.new(objname, curvedata)
objectdata.location = (0,0,0) #object origin
bpy.context.scene.objects.link(objectdata)
polyline = curvedata.splines.new('POLY')
polyline.points.add(len(cList)-1)
for num in range(len(cList)):
polyline.points[num].co = (cList[num])+(w,)
polyline.order_u = len(polyline.points)-1
polyline.use_endpoint_u = True
# here's a sample vector
listOfVectors = [(0.81,3.443,-5.723),(3.723,0.111,-5.603)]
MakePolyLine("NameOfMyCurveObject", "NameOfMyCurve", listOfVectors)
With more vectors, this would output something like the image below:
