geolib.py 文件源码

def cT_helper(x, y, z, in_srs, out_srs):
    """Helper function that wraps osr CoordinatTransformation
    """
    x, y, z = np.atleast_1d(x), np.atleast_1d(y), np.atleast_1d(z)
    #Handle cases where z is 0 - probably a better way to use broadcasting for this
    if x.shape[0] != z.shape[0]:
        #Watch out for masked array input here
        orig_z = z[0]
        z = np.zeros_like(x)
        z[:] = orig_z
    orig_shape = x.shape
    cT = osr.CoordinateTransformation(in_srs, out_srs)
    #x2, y2, z2 = zip(*[cT.TransformPoint(*xyz) for xyz in zip(x, y, z)])
    x2, y2, z2 = list(zip(*[cT.TransformPoint(*xyz) for xyz in zip(np.ravel(x),np.ravel(y),np.ravel(z))]))
    if len(x2) == 1:
        x2, y2, z2 = x2[0], y2[0], z2[0] 
    else:
        x2 = np.array(x2).reshape(orig_shape)
        y2 = np.array(y2).reshape(orig_shape)
        z2 = np.array(z2).reshape(orig_shape)
    return x2, y2, z2