Geodesic distance between geometry shapes Python
Having a dfA with a column called geometry with the following geometrical shapes:
d = {'id': [1, 2], 'geometry': ['POINT (70.66000 33.45000)', 'POINT (74.08000 4.60000)']}
dfA = pd.DataFrame(data=d)
dfA
  id  geometry 

 0  1  POINT (70.66 33.45) 
 1  2  POINT (74.08 4.6) 
I would like to calculate the minimum geodesic distance with each of the geometric shapes of the dfB's geometry column:
d = {'id': [1, 2, 3], 'geometry': ['LINESTRING (58.66000 34.58000, 59.66000 35.58000)', 'LINESTRING (47.91000 15.78000, 48.91000 16.78000)', 'POINT (66.86000 10.48000)']}
dfB = pd.DataFrame(data=d)
dfB
  id  geometry 

 0  1  LINESTRING (58.66 34.58, 59.66 35.58) 
 1  2  LINESTRING (47.91 15.78, 48.91 16.78) 
 2  3  POINT (66.86 10.48) 
I have tried to do this calculation using the Python shapely and geopandas libraries by following the steps below:
from shapely import wkt
import geopandas as gpd
dfA['geometry'] = dfA['geometry'].apply(wkt.loads)
dfA = gpd.GeoDataFrame(dfA, geometry='geometry')
dfB['geometry']= dfB['geometry'].apply(wkt.loads)
for i, value in dfB.iterrows():
e = dfB.iloc[i]['id']
dfA[str(e)] = dfA['geometry'].distance(dfB.iloc[i]['geometry'])
dfA
  id  geometry  1  2  3 

 0  1  POINT (70.66 33.45)  11,20432506  27,40349248  44,09404608 
 1  2  POINT (74.08 4.6)  42,10521108  33,0247377  9,311433832 
Unfortunately, shapely distance function calculates the Euclidean Distance and not the geodesic distance.
Another strategy to follow would be to use a function that calculates the geodesic distance from point A to all points on line B [B1, B2, B3,...] and keep the minimum distance. That is to say: dist_AB = min(geodist(A, B1), geodist(A, B2), geodist(A, B3), ....)
This solution works but computationally it is very expensive since we are talking about a calculation from thousands of points against thousands of lines. Any other more optimal way to perform this calculation will be of a lot of help.
See also questions close to this topic

How to apply different text formats to different text in python tkinter Text widget
I am writing a python code to create a RICH TEXT EDITOR in pythonmy python text editor. It has a bold button and all formating buttons like italic,fonts etc.
These buttons would format the selected text.
But the problem with all these buttons is that First when I bold the selected text, it is working fine but at the second time it applies the same format to the previously selected text
This is my code for the bold and italic function :
def boldtxt(): bldbtn['activebackground']="#ffad33" contboxfont['weight']="bold" contbox.tag_add("bold",SEL_FIRST,SEL_LAST) contbox.tag_config("bold",font=contboxfont) def italictxt(): itlbtn['activebackground']="#ffad33" contboxfont['slant']="italic" contbox.tag_add("italic",SEL_FIRST,SEL_LAST) contbox.tag_config("italic",font=contboxfont) contbox=Text() contbox.pack() bldbtn=Button(comand=boldtxt) bldbtn.pack() itlbtn=Button(command=italictxt) itlbtn.pack()
I want that every word (or selected item should be formatted according to the user. For ex. if the bold button is pressed the selected text should be bold and if the italic button is pressed then the selected text would be italic. Thanks in advance :)

A Sympy/IPython method to display function label with expression
I've been searching for a way to display both the function label as well as the function expression when working in Jupyter or any other IPython terminal, but I never managed to find a solution online.
In short I wanted to say something like
display(f)
and it should return:
f(x) = a_{2}x^{2} + a_{1}x + a_{0}
rather than just the terse:
a_{2}x^{2} + a_{1}x + a_{0}So I've decided to share what I came up with in case others have similar issues.
I find it particularly useful especially when you want to hide the code in your Jupyter notebook as described in How to hide code in cells in IPython Notebook
Here is the function:
from sympy import latex from IPython.display import Math def parade(expr, func, y=None, label=None): """Returns a latex expression to display using IPython.display.display() Parameters: repr: the representation of the function e.g. 'f(x)' func: the sympy functions symbol e.g. f y (optional): The value or expression this function is equal to. label (optional): Append a latex formula label Returns: the latex expression to display. """ expression = ('{expr} = {func}'.format(expr=expr, func=latex(func))) if y is not None: expression += ' = {y}'.format(y=latex(y)) if label is not None: expression += '\\tag{label}'.format(label=latex(label)) return Math(expression)
To use the function type:
from sympy import symbols from IPython.display import display x, a_0, a_1, a_2 = symbols('x a_0 a_1 a_2') f = a_2 * x**2 + a_1 * x + a_0 display(parade('f(x)', f))
and it outputs:
f(x) = a_{2}x^{2} + a_{1}x + a_{0}If anyone has a way to improve on this, I'd appreciate it if you could share it.

try/except issue in Python
I ran the following code as a practice purpose from one of the Python book and I want the output as shown in the URL provided below. So, when I am running the first activity in the question (from book  shown in the photo) to verify if he has done it correctly, it's working but when he ran the second activity, it's not working properly. Ideally, I should get except statement printed if I am entering any string in the first user prompted input but I am being asked for the second user prompted input i.e., rate. Why is that so? And also please tell if it is advisable to put complete code in try portion or ideally we should only put the portion in a try statement for which we are sure that it wouldn't work out and it would jump to except block?
When I enter forty, I am not given the error message which is in the book which means that there can be some changes made in the code. However, it is working fine when I am including all lines of codes under try instead of just 2 lines which are currently there i.e., (fh = float(sh) and fr = float(sr)). How can I correct my existing written code by just writing two statements in a try portion?
Your help would be appreciated.
Issue: Enter Hours: fortyError, please enter numeric input
Image of the required output is below:
Following is the code:
sh = input("Enter hours:") sr = input("Enter Rate:") try: fh = float(sh) fr = float(sr) except: print('Error, please enter numeric input') quit() print(fh, fr) if fh > 40: reg = fr * fh otp = (fh  40.0) * (fr * 0.5) xp = reg + otp else: xp = fh * fr print("Pay:",xp)

Confusion with Arrays
I'm new to JavaScript and although I have looked at many forums/tutorials I'm struggling with exactly what I need to change to achieve what I want.
I currently have a Google spreadsheet that is populated with UK post codes, these postcodes are automatically added at various points using a combination of google forms and Microsoft forms/flow.
What I want to do is to automatically calculate the distance in miles between 2 post codes that are in each row.
For example, in column C there would be a pickup post code, and in column D there would be a Drop off post code  in column E I would like the mileage to automatically display when a new row with new post codes is added.
Currently, I am using the below customer function to calculate the distance in miles, but I need to be able to wrap it in the
ARRAYFORMULA
function in google sheets but when I do only the first cell displays the correct mileage (the others remain blank)The code below seems to work very well I just can't figure out how to be able to allow the
GOOGLEDISTANCE
function to accept an array of Pick up post codes and an array of drop off post codes and then return an array of the calculated mileage in miles.Code:
function GOOGLEDISTANCE(pick_up,drop_off,return_type) { Utilities.sleep(1000); var travelMode = Maps.DirectionFinder.Mode.BICYCLING; var directions = Maps.newDirectionFinder() .setRegion('UK') .setLanguage('enGB') .setOrigin(pick_up) .setDestination(drop_off) // .addWaypoint(waypoint1) .setMode(travelMode) .getDirections(); if (directions.status !== "OK") return "Error: " + directions.status; var route = directions.routes[0].legs[0]; var time = route.duration.value; var distance = route.distance.value; var steps = route.steps.map(function(step) { return step.html_instructions.replace(/<[^>]+>/g, ""); }).join("\n"); switch(return_type) { case "MILES": case "miles": return distance * 0.000621371; break; case "KILOMETERS": case "kilometers": return distance / 1000; break; case "MINUTES": case "minutes": return time / 60; break; case "HOURS": case "hours": return time / 60 / 60; break; case "STEPS": case "steps": return steps; break; default: // Default to miles return distance * 0.000621371; } }

Difference between distm function or the distVincentyEllipsoid in R
Could you fully explain the big difference in using the
distm
function or thedistVincentyEllipsoid
function to calculate the distance of geodesic coordinates in R?I noticed that using distm for this calculation, it takes much longer. Could you please explain to me beyond the difference, why does this happen?
Thank you!

Float distance for 80 bit long double
For
float
,double
, and__float128
, we have the following code which computes the number of representables between them:int32_t float_distance(float x, float y) { static_assert(sizeof(float) == sizeof(int32_t), "float is incorrect size."); int32_t xi = *reinterpret_cast<int32_t*>(&x); int32_t yi = *reinterpret_cast<int32_t*>(&y); return yi  xi; } int64_t float_distance(double x, double y) { static_assert(sizeof(double) == sizeof(int64_t), "double is incorrect size."); int64_t xi = *reinterpret_cast<int64_t*>(&x); int64_t yi = *reinterpret_cast<int64_t*>(&y); return yi  xi; } __int128_t float_distance(__float128 x, __float128 y) { static_assert(sizeof(__float128) == sizeof(__int128_t), "quad is incorrect size."); __int128_t xi = *reinterpret_cast<__int128_t*>(&x); __int128_t yi = *reinterpret_cast<__int128_t*>(&y); return yi  xi; }
(This code works for x,y > 0, for legibility we aren't dealing with the general case.)
There is no
int80_t
, so what is the analogous code forlong double
? Out of desperation, I tried both__int128_t
and__int64_t
return types, but no joy.Edit: It seems like this is not a wellknown property of IEEE754 floating point numbers. Here's a good read for those who are wondering why this works.

How to find intersections between datasets in Geopandas
I have two datasets in Projected CRS: EPSG:32618, I am trying to find the intersection of the 1000m buffer on one point with the points in other datasets, I have plotted that one buffer polygon over the points I am trying to intersect the plot shows that the intersection exists while the overlay and within functions of geopandas show no intersection.
Han anybody ever faced this issue and is there a workaround? Thanks.

Can you store a constant global object in Django?
I am working on a small webpage that uses geospacial data. I did my initial analysis in Python using GeoPandas and Shapely and I am attempting to build a webpage from this. The problem is, when using Django, I can't seem to find a way to keep the shape file stored as a constant object. Each time a request is made to do operations on the shapefile, I need to load the data from source. This takes something like 6 seconds while a standard dataframe deep copy
df.copy()
would take fractions of a second. Is there a way I can store a dataframe in Django that can be accessed and deep copied by the views without rereading the shapefile? 
Exception has occurred: SSLError when adding basemap with Contextily to GeoPandas map
I'm trying to use Contextily to add a basemap to a map created with GeoPandas, but am getting the following error:
HTTPSConnectionPool(host='a.tile.openstreetmap.org', port=443): Max retries exceeded with url: /14/8452/7869.png (Caused by SSLError(SSLError("bad handshake: Error([('SSL routines', 'tls_process_server_certificate', 'certificate verify failed')])"))) File "D:\generate_udf_maps.py", line 42, in generate_map ctx.add_basemap(ax, source=ctx.providers.OpenStreetMap.Mapnik) File "D:\generate_udf_maps.py", line 60, in <module> generate_map(gdb_path, None, operations_file)
Here is my code:
import geopandas as gpd import pandas as pd import matplotlib.pyplot as plt from shapely.geometry import Point, Polygon import contextily as ctx fig, ax = plt.subplots() p1 = Point(640864.516, 784724.566) p2 = Point(638594.804, 788492.544) p3 = Point(644076.197, 788141.637) p4 = Point(641089.953, 789539.343) df = gpd.GeoDataFrame(geometry=[p1, p2, p3, p4]) df.plot(ax=ax) ctx.add_basemap(ax, source=ctx.providers.OpenStreetMap.Mapnik) # Error happens here. If this line is not present I get a map with just the points. plt.show()
I have tried different map providers to no avail and saw several links giving suggestions related to a reference package, which I'm not sure are applicable here. Thanks in advance!

using qgis and shaply error: GEOSGeom_createLinearRing_r returned a NULL pointer
I tried to create a polygon shapefile in QGIS and read it in python by shapely. An example code looks like this:
import fiona from shapely.geometry import shape multipolys = fiona.open(somepath) multi = multipolys[0] coord = shape(multi['geometry'])
The EOSGeom_createLinearRing_r returned a NULL pointer I checked if the polygon is valid in QGIS and no error was reported. Actually, it does not work for even a simple triangle generated in QGIS. Anyone knows how to solve it?
Thank you

Why does zip change input floats?
I have list of shapely Points in GeoSeries.
coords.head():
0 POINT (26.17690 80.81700) 1 POINT (15.54390 80.61700) 2 POINT (20.67690 80.36700) 3 POINT (6.10610 80.83300) 4 POINT (17.63910 79.88300) Name: geometry, dtype: geometry
When trying to get list of zipped coordinates using command
pd.Series(zip(coords.geometry.x, coords.geometry.y)).head()
I get next sample:0 (26.1769, 80.817) 1 (15.5439, 80.617) 2 (20.6769, 80.367) 3 (6.1061000000000005, 80.833) 4 (17.63909999999999, 79.883) dtype: object
By the way example of
coords.geometry.x.head()
:0 26.1769 1 15.5439 2 20.6769 3 6.1061 4 17.6391 dtype: float64
Also strange thing is that when I try to reproduce results:
new_coords = [(26.17690, 80.81700), (15.54390, 80.61700), (20.67690, 80.36700), (6.10610, 80.83300), (17.63910, 79.88300)] new_coords = gpd.GeoSeries([Point(p) for p in new_coords]) pd.Series(zip(new_coords.geometry.x, new_coords.geometry.y)) new_coords
Zip doesn't behave oddly:
0 POINT (26.17690 80.81700) 1 POINT (15.54390 80.61700) 2 POINT (20.67690 80.36700) 3 POINT (6.10610 80.83300) 4 POINT (17.63910 79.88300) dtype: geometry
Main goal here is to get accurate value of coordinates in order to merge frames of data, thus it is not acceptable that zip returns similar value but not the same.

Transform GPS coordinates data of a point to new point with different coordinates by only adding distance
If I have a point X1(lat1,lon1) and a distance d. I we want a new point by adding distance d horizontally or vertically to X1. How can I do this in MATLAB?

Converting latitidue, longitude to image x/y over a maptilegl generated map
I forked tileservergl and modified the static endpoint to return the center lat/lng and zoom in the response headers with the image. My app always calls this endpoint as styles/klokantechbasic/static/auto/...?path=...
I then use those parameters to draw new markers on the image at some selected lat/lng coordinates, according to my application requirements.
In order to project the lat/lng coordinates to pixel (x,y) offsets from upper left corner, I use the following transformation, given in Pythonic Pseudocode:
Let (h, w) be the height and width of the image returned by tileservergl. px_radius = (256 * (2 ** zoom)) / (2 * PI) # pixels dy_rad = radians(lat  center_lat) # angle from center on yaxis x = w/2 + px_radius * radians(lng  center_lng) y = h/2  px_radius * 1/2 * ln[(1 + sin(dy_rad)) / (1  sin(dy_rad))]
To test this, I called the static endpoint using the outline of Utah as a path. Then my code uses the above transformation to draw a marker at the upper left and lower right corners of the state boundary. As you can in the image, the x coordinate of each point aligns with the east and west boundary lines, but the y coordinates fall short on either end.
Since the scale is smaller than it should be only in the y direction, I suspect that the difference might result from the fact that tileservergl and mercator are using a projection based on the earth being an elipsoid (fatter around the equator). Please don't hold me to this  I am pulling at straws now.
I spent several hours studying the code in tileservergl and the code in mapbox/sphericalmercator. It appears that sphericalmercator is transforming the (lat, lng) coordinates before converting them to pixels, using this code:
var to3857 = proj4('EPSG:3857'); var toDataProj = proj4(info.proj4); dataProjWGStoInternalWGS = function(xy) { return to3857.inverse(toDataProj.forward(xy)); };
To test this, I tried using the Python library pyproj (for "Python Projection") as so:
proj = Proj(proj='merc', ellps='WSG84', lat_ts=lat_ctr) (x_center, y_center) = proj(lng_center, lat_center) # meters (x, y) = proj(lng, lat) # meters x = w/2 + (x  x_center) * px_radius / 6378137.0 y = h/2  (y  y_center) * px_radius / 6378137.0
This time, both the x and y resulting pixel coordinates fall short of the path coordiantes.
It's clear that I don't know how to use the pyproj library to mimic what tileservergl is doing. What other values for 'ellps' might I use?