Where communities thrive


  • Join over 1.5M+ people
  • Join over 100K+ communities
  • Free without limits
  • Create your own community
People
Activity
    Walter Lorenzetti
    @wlorenzetti
    I just try with pyproj 2.2.2 with Transformer object: 
    # https://pyproj4.github.io/pyproj/stable/gotchas.html#upgrading-to-pyproj-2-from-pyproj-1

from pyproj import Transformer
import time

start = time.time()
x = 61759.831
y = -41125.585
transformer = Transformer.from_crs("+proj=cass +lat_0=43.318293 +lon_0=11.332212 +x_0=0.000000000000000 +y_0=0.000000000000000 +datum=WGS84 +units=m +no_defs", "epsg:3003")
ret = transformer.transform(y, x)
stop = time.time()
print "TIME: {}".format(stop - start)
print ret
    but timing in not good again.. 
    TIME: 0.0580408573151
(1646302.7175560538, 4859396.61741957)
    I' tryee inside a 20 time cycle to test memory effect: 
    TIME: 0.0537190437317
(1646302.7175560538, 4859396.61741957)
TIME: 0.0591440200806
(1646302.7175560538, 4859396.61741957)
TIME: 0.0580699443817
(1646302.7175560538, 4859396.61741957)
TIME: 0.0595941543579
(1646302.7175560538, 4859396.61741957)
TIME: 0.0600090026855
(1646302.7175560538, 4859396.61741957)
TIME: 0.059278011322
(1646302.7175560538, 4859396.61741957)
TIME: 0.0588028430939
(1646302.7175560538, 4859396.61741957)
TIME: 0.0584590435028
(1646302.7175560538, 4859396.61741957)
TIME: 0.0585939884186
(1646302.7175560538, 4859396.61741957)
TIME: 0.0572869777679
(1646302.7175560538, 4859396.61741957)
TIME: 0.0574839115143
(1646302.7175560538, 4859396.61741957)
TIME: 0.0569059848785
(1646302.7175560538, 4859396.61741957)
TIME: 0.0589108467102
(1646302.7175560538, 4859396.61741957)
TIME: 0.0600910186768
(1646302.7175560538, 4859396.61741957)
TIME: 0.0592269897461
(1646302.7175560538, 4859396.61741957)
TIME: 0.058445930481
(1646302.7175560538, 4859396.61741957)
TIME: 0.0585761070251
(1646302.7175560538, 4859396.61741957)
TIME: 0.0573010444641
(1646302.7175560538, 4859396.61741957)
TIME: 0.067663192749
(1646302.7175560538, 4859396.61741957)
TIME: 0.0576071739197
(1646302.7175560538, 4859396.61741957)
    But I can see a better timing...
    Joris Van den Bossche
    @jorisvandenbossche
    So what I said above: it will only be faster if you use the transform object repeatedly. For a single transformation, it will still be slower compared to pyproj 1.9.6
    Walter Lorenzetti
    @wlorenzetti
    Oh sorry I have putted 'transformer' instance creation inside the cycle, now the timing is good! Thanks: 
    TIME: 6.41345977783e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.50203704834e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.21593475342e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.31130218506e-05
(1646302.7175560538, 4859396.61741957)
TIME: 9.77516174316e-06
(1646302.7175560538, 4859396.61741957)
TIME: 1.09672546387e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.38282775879e-05
(1646302.7175560538, 4859396.61741957)
TIME: 9.05990600586e-06
(1646302.7175560538, 4859396.61741957)
TIME: 1.00135803223e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.09672546387e-05
(1646302.7175560538, 4859396.61741957)
TIME: 8.10623168945e-06
(1646302.7175560538, 4859396.61741957)
TIME: 7.86781311035e-06
(1646302.7175560538, 4859396.61741957)
TIME: 1.69277191162e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.09672546387e-05
(1646302.7175560538, 4859396.61741957)
TIME: 8.10623168945e-06
(1646302.7175560538, 4859396.61741957)
TIME: 1.4066696167e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.4066696167e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.09672546387e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.71661376953e-05
(1646302.7175560538, 4859396.61741957)
TIME: 1.69277191162e-05
(1646302.7175560538, 4859396.61741957)
    donglx2018
    @donglx2018
    hi
    pyproj fails to compille when I pip install basemap
    donglx2018
    @donglx2018
    the error message is as follows:
    src/_proj.c:6754:24: error: 'PyThreadState' {aka 'struct _ts'} has no member named 'exc_type'; did you mean 'curexc_type'?
    tmp_type = tstate->exc_type;
    ^~~~
    curexc_type
    any idea about this error?
    Alan D. Snow
    @snowman2
    What is your operating system version? Python version?
    David Hoese
    @djhoese
    @snowman2 I'm bugging you here since it seems like the quickest way to find an answer, sorry. I'm hacking together a custom MapServer where I want to serve data in the +proj=geos projection. I was thinking of defining my own EPSG code for this, but with the changes to PROJ where it uses a sqlite DB now, is it possible to do this any more?
    Alan D. Snow
    @snowman2
    You could potentially just modify the DB directly. I have seen some connect to it here: https://github.com/ioos/compliance-checker/blob/a5b457a44b609b9da104897a7ea46c4f62cc42d4/compliance_checker/cf/cf.py#L4244-L4253
    Joshua Arnott
    @snorfalorpagus
    Is pyproj thread safe? There is mention in the changelog for 1.9.0, but it doesn't explicitly state in the documentation that it's OK to call in multiple threads.
    Alan D. Snow
    @snowman2
    It should be. I would recommend version 2.4 for that. If you see any issues with threading, please raise and issue.
    Joris Van den Bossche
    @jorisvandenbossche
    It's threadsafe as long as you are using different objects in different threads, right? (I seem to remember some recent discussion about that)
    Alan D. Snow
    @snowman2

    It's threadsafe as long as you are using different objects in different threads, right?

    Sounds right. Each object has a context attached and only one context per thread is considered safe. IIRC

    Joshua Arnott
    @snorfalorpagus

    Each object has a context attached and only one context per thread is considered safe.

    OK. That's probably my problem as I'm creating transforms as globals. I'll move them into a thread-local context.

    David Hoese
    @djhoese
    @snowman2 Does pyproj and/or rasterio ever check the other projection definitions defined in the PROJ library? I'm wondering if they ever read the other.extra file (https://github.com/OSGeo/PROJ/blob/master/data/other.extra). From what I can tell they don't and I don't think they are included in the PROJ database. So are the files just there for backwards compatibility?
    Alan D. Snow
    @snowman2
    I am not sure about rasterio/GDAL, but pyproj would only use it if PROJ uses it. pyproj itself does not directly access it. I would recommend asking the PROJ devs about it.
    Ilya Ashikhmin
    @ilalex

    I'm a little bit confused about Proj initialization. I have a following code for geometry transformation:

    from functools import partial

import pyproj
from pyproj import Proj
from shapely.geometry import Point
from shapely.ops import transform

if __name__ == '__main__':
    p = Point(400_000, 6_000_000)
    old_proj = partial(pyproj.transform, Proj(init="epsg:32649"), Proj(init="epsg:4326"))
    new_proj = partial(pyproj.transform, Proj("epsg:32649"),      Proj("epsg:4326"))
    print(transform(old_proj, p))
    print(transform(new_proj, p))

    It results in: 

    /home/ilya/.cache/pypoetry/virtualenvs/pyprojbug-Ra-0DuD_-py3.6/lib/python3.6/site-packages/pyproj/crs.py:77: FutureWarning: '+init=<authority>:<code>' syntax is deprecated. '<authority>:<code>' is the preferred initialization method.
  return _prepare_from_string(" ".join(pjargs))
/home/ilya/.cache/pypoetry/virtualenvs/pyprojbug-Ra-0DuD_-py3.6/lib/python3.6/site-packages/pyproj/crs.py:77: FutureWarning: '+init=<authority>:<code>' syntax is deprecated. '<authority>:<code>' is the preferred initialization method.
  return _prepare_from_string(" ".join(pjargs))
POINT (109.4692987436448 54.13837331779994)
POINT (54.13837331779994 109.4692987436448)

    Why is Proj(init="epsg:32649") and Proj("epsg:32649") giving different results?

    Alan D. Snow
    @snowman2
    I believe you are running into axis order changes http://pyproj4.github.io/pyproj/stable/gotchas.html#axis-order-changes-in-proj-6. I would recommend reading the link. I would also recommend using the Transformer class with always_xy for a simple fix.
    Topi Filppula
    @tfilppula
    How should I use pipelines? For example a really simple one:
    string = "+proj=pipeline +ellps=GRS80 +step +proj=cart"
pipe = Transformer.from_pipeline(string)
pipe.transform(50, 25, errcheck=True)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/usr/local/lib/python3.7/site-packages/pyproj/transformer.py", line 446, in transform
    errcheck=errcheck,
  File "pyproj/_transformer.pyx", line 463, in pyproj._transformer._Transformer._transform
pyproj.exceptions.ProjError: transform error: latitude or longitude exceeded limits
    Whereas cct seems to work:
    echo 50 25 0 | cct +proj=pipeline +ellps=GRS80 +step +proj=cart
 3717892.6072   4430811.8715  2679074.4629           inf
    Topi Filppula
    @tfilppula
    Python version is 3.7.6, pyproj 2.5.0
    Alan D. Snow
    @snowman2
    Mind raising an issue for this on GitHub?
    Topi Filppula
    @tfilppula
    Issue raised, hope it fulfils (at least most) requirements.
    David Hoese
    @djhoese

    @snowman2 Coming with another shot in the dark question: Have you had any issues pop up recently with unicode issues with pyproj or PROJ? We have a Satpy user on Windows that is having random (I think) errors with a particular PROJ dict being passed to pyproj. It errors and says either:

    pyproj.exceptions.ProjError: Invalid projection b'\x02\x01roj=geos +sweep=x +lon_0=-75 +h=35786023 +x_0=0 +y_0=0 +ellps=GRS80 +units=m +no_defs'.: (Internal Proj Error: proj_create: unrecognized format / unknown name)

    Or:

      File "pyproj/_crs.pyx", line 24, in pyproj._crs.cstrdecode
  File "C:\ProgramData\Miniconda3\lib\site-packages\pyproj\compat.py", line 21, in pystrdecode
    return cstr.decode("utf-8")
UnicodeDecodeError: 'utf-8' codec can't decode bytes in position 8-9: invalid continuation byte

    I'm currently working with them to narrow this down, but it seems related to a particular type of processing with various dask workers, etc. Note the PROJ string they are using stays the same and will sometimes work. Any ideas?

    Alan D. Snow
    @snowman2
    Unfortunately, no. This is all new to me.
    David Hoese
    @djhoese
    @snowman2 I noticed in pyproj 2.5+ that to_cf() produces a reference_ellipsoid_name with a space in it. At least for WGS 84. Is this intended? Any idea if this is PROJ doing it or pyproj?
    of course, I could file a bug but wanted to know if it was intended first
    Alan D. Snow
    @snowman2
    @djhoese, this is intended. Since it is WKT based, it uses the form of the ellipsoid name in the EPSG database. The old version used the PROJ string version of the ellipsoid name.
    David Hoese
    @djhoese
    ah great. That makes sense. Thanks.
    David Hoese
    @djhoese
    @snowman2 I'm coming back to the unicode issue I mentioned above. I've now hit it in my own processing. I've started documenting things here: pytroll/satpy#1114 and I'm working on reproducing it. If you look at the tracebacks in that issue that I posted, my guess is that pyproj's CRS object isn't being hashed/serialized properly. Or at least it doesn't seem to be thread safe. If you have any ideas what might be going on it would be great if you could let me know.
    Alan D. Snow
    @snowman2
    Huh, not sure. Mind opening up an issue on github?
    David Hoese
    @djhoese
    I can, I thought I'd wait until I could reproduce it, but so far I can't. I'll make an issue in a bit
    georg1me
    @georg1me
    Hello all. I am trying to create a data visualization on a map, and then publish it as a django app on my portfolio website. I was able to successfully use cartopy and proj in spyder. However, once I created a virtual environment in vscode to finish the 'production' version, I can't get conda in vscode to recognize that I have already installed the dependencies for cartopy, including proj. I am pretty comfortable using virtual environments, but I have never encountered a package that I can't install directly through pip. Any advice would be greatly appreciated. Thanks!
    Alan D. Snow
    @snowman2
    I think that you are looking for: https://gitter.im/SciTools/cartopy
    georg1me
    @georg1me
    Thank you @snowman2
    David Hoese
    @djhoese
    @snowman2 Just because you hang out with a different crowd than me, I'm wondering if you know the best options for finding intersections between polygons on two separate projections? And of the solutions you know about, do they handle things if points on the polygons don't exist in the other projection? Is there something that should work for all Earth-based CRSes
    Alan D. Snow
    @snowman2
    Just seeing this now. There are several options, but in general the idea is to get the polygons in the same projection. Shapely is a good resource for intersects: https://shapely.readthedocs.io/en/latest/manual.html#object.intersects
    Reproject example: https://shapely.readthedocs.io/en/latest/manual.html#shapely.ops.transform
    If you have many polygons, you could look into https://geopandas.org/ and https://github.com/pygeos/pygeos/
    David Hoese
    @djhoese
    Thanks @snowman2. I think the other satpy developer working on this was using shapely and transforming already. I think the problem was what to do if one of the vertices of the polygon can't be transformed to the other projection.
    Alan D. Snow
    @snowman2
    Ah, I see what you mean. Does the problem occur in both transform directions? Or is it possible to try transforming to the other CRS to see if all points can be transformed?