用python进行快速坐标匹配

经常碰到的一个问题,很多个目标多次观测组成的列表,需要对这个列表进行分类。由于工作中主要用python,如果使用Python通过两次循环来对列表进行坐标匹配,一旦列表的长度稍长(如万的量级),时间就长的让人难以忍受。Python在处理循环时的速度实在让人不敢恭维。

所幸有一些其他方法可以直接做这件事,这里介绍一个可以快速坐标匹配的包,esutil中的htm模块。这个模块使用Hierarchical Triangular Mesh方法来进行坐标匹配(参见Hierarchical Triangular Mesh, Grid (spatial index)等),可以极快的完成大量数据的坐标匹配工作。

下面写一个通过htm进行坐标匹配的例子

import sys
import numpy as np
import pandas as pd
from esutil import htm

def main():
    name = sys.argv[1]
    df = pd.read_csv(name)
    ralst = np.array(df.ra)
    declst = np.array(df.dec)
    mat = htm.Matcher(15, ralst, declst)
    k0, k1, err = mat.match(ralst, declst, radius=0.00028)
    print (k0, k1, err)
    unique_objID = np.unique(k1)
    print unique_objID

if __name__ == '__main__':
    main()

这里需要先指定被匹配的坐标列表,我们是列表自己和自己匹配,两者设置成一样就行。默认匹配到的个数是1,一旦搜索到匹配的对象就会停止搜索,这样用来被匹配的列表中目标即使有重复也没有关系,因为总是只会搜索到这些重复目标中的第一个就停止。返回的元素k0是匹配的列表下表,k1是其在mat中的对应目标的下标列表, err是两者之间的偏差(肯定小于radius)。

某次执行的输出结果如下

(array([ 0, 1, 2, ..., 41598, 41599, 41600]), array([ 37, 37, 37, ..., 41589, 41589, 41589]), array([0.00000000e+00, 0.00000000e+00, 0.00000000e+00, ...,
8.53773646e-07, 8.53773646e-07, 8.53773646e-07]))
[ 37 86 135 184 233 282 331 380 429 478 527 576
625 674 723 772 821 870 919 968 1017 1066 1115 1164
1213 1262 1311 1360 1409 1458 1507 1556 1605 1654 1703 1752
1801 1850 1899 1948 1997 2046 2095 2144 2193 2242 2291 2340
2389 2438 2487 2536 2585 2634 2683 2732 2781 2830 2879 2928
2977 3026 3075 3124 3173 3222 3271 3320 3369 3418 3467 3516
3565 3614 3663 3712 3761 3810 3859 ... ]

对41000个元素的列表匹配费时0.76s,可以说很快速了。

Visits: 336

发表回复

您的电子邮箱地址不会被公开。 必填项已用 * 标注

*