SQL使用PostgreSQL

使用場景為生成馬賽克風(fēng)格的圖片，即把某個區(qū)域分割成若干個小的區(qū)域，計算每個小的區(qū)域中包含的某些體量，用不同顏色顯示。

1. 矩形（rectangle）

原理：每一行/列小矩形的中心點連線都是平行的，根據(jù)橫縱切割的數(shù)量可以算出所有的中心點，再根據(jù)每個中心點計算出四個頂點的坐標(biāo)。（如果想展現(xiàn)正方形可以1.手動控制橫縱切割的比;2.簡單修改代碼邏輯，只傳入橫向切割數(shù)量，用邊長計算縱向數(shù)量。）

def cal_rasterized_points_from_boundary_v2(obj, boundary_li, lng_count, lat_count):
    """
    計算邊界內(nèi)柵格化后的中心點，返回若干個矩形的坐標(biāo)數(shù)組 -- 馬賽克凹凸有致版
    :param obj: 多邊形的object
    :param boundary_li: multipolygon.boundary.coords轉(zhuǎn)化后的list
    :param lng_count: 橫向切割數(shù)量
    :param lat_count: 縱向切割數(shù)量
    :return: 
    """
    # 目前只考慮MultiPolygon包含一個多邊形的情況。
    boundary_list = boundary_li[0]
    lngs = []
    lats = []
    for bl in boundary_list:
        lngs.append(bl[0])
        lats.append(bl[1])
    max_lng = float(max(lngs))
    min_lng = float(min(lngs))
    max_lat = float(max(lats))
    min_lat = float(min(lats))
    max_lng_amap, max_lat_amap = wgs84togcj02(max_lng, max_lat)
    min_lng_amap, min_lat_amap = wgs84togcj02(min_lng, min_lat)

    lng_gap = (max_lng_amap - min_lng_amap) / lng_count
    lat_gap = (max_lat_amap - min_lat_amap) / lat_count
    data_lngs = []
    data_lats = []
    cps = []
    result = []
    for i in range(lng_count):
        lng = min_lng_amap + lng_gap * 0.5 + lng_gap * i
        if lng < max_lng_amap:
            data_lngs.append(lng)
    for i in range(lat_count):
        lat = min_lat_amap + lat_gap * 0.5 + lat_gap * i
        if lat < max_lat_amap:
            data_lats.append(lat)
    for lat in data_lats:
        for lng in data_lngs:
            cps.append([lng, lat])
    for cp in cps:
        # TODO: To judge cp whether in boundary or not
        lng, lat = gcj02towgs84(cp[0], cp[1])
        point = Point((lng, lat))
        li = obj.filter(boundary__contains=point)
        if not li:
            continue
        result.append([[cp[0] - lng_gap * 0.5, cp[1] - lat_gap * 0.5],
                       [cp[0] + lng_gap * 0.5, cp[1] - lat_gap * 0.5],
                       [cp[0] + lng_gap * 0.5, cp[1] + lat_gap * 0.5],
                       [cp[0] - lng_gap * 0.5, cp[1] + lat_gap * 0.5]])

    return result

2. 六邊形（hexagonal）

原理：

20170124_01_pic_004.png

六邊形可以拆分成六個等邊三角形。

20170124_01_pic_010.png

每個中心點的位置每兩行一個循環(huán)。同行兩個中心點間距為三倍的邊長，行間距為√3/2倍邊長（經(jīng)緯度間距離轉(zhuǎn)換存在差異，現(xiàn)采用緯度間*0.88）

def cal_rasterized_points_from_boundary(obj, boundary_li, lng_count):
    """
    計算邊界內(nèi)柵格化后的中心點，返回若干個六邊形的坐標(biāo)數(shù)組 -- 馬賽克凹凸有致版  
    蜂巢六邊形切割 -- 粗略用lng的差來表示邊長，因為分割后數(shù)值很小，可以忽略坐標(biāo)和距離轉(zhuǎn)換的誤差
    :param obj: 多邊形的object
    :param boundary_li: multipolygon.boundary.coords轉(zhuǎn)化后的list
    :param lng_count: 橫向切割數(shù)量
    :return: 
    """
    # 目前只考慮MultiPolygon包含一個多邊形的情況。
    boundary_list = boundary_li[0]
    lngs = []
    lats = []
    for bl in boundary_list:
        lngs.append(bl[0])
        lats.append(bl[1])
    max_lng = float(max(lngs))
    min_lng = float(min(lngs))
    max_lat = float(max(lats))
    min_lat = float(min(lats))
    max_lng_amap, max_lat_amap = wgs84togcj02(max_lng, max_lat)
    min_lng_amap, min_lat_amap = wgs84togcj02(min_lng, min_lat)

    lng_gap = (max_lng_amap - min_lng_amap) / lng_count
    side_length = lng_gap / 3
    lat_gap = side_length * math.sqrt(3) / 2 * 0.88
    lat_count = int((max_lat_amap - min_lat_amap) / lat_gap) + 1
    data_lngs_first = []
    data_lngs_second = []
    data_lats = []
    cps = []
    result = []
    for i in range(lng_count):
        lng = min_lng_amap + lng_gap * 0.5 + lng_gap * i
        if lng < max_lng_amap:
            data_lngs_first.append(lng)
    for lng in data_lngs_first:
        sec_lng = lng - side_length * 3 / 2
        if sec_lng < max_lng_amap:
            data_lngs_second.append(sec_lng)

    for i in range(lat_count):
        lat = min_lat_amap + lat_gap * 0.5 + lat_gap * i
        if lat < max_lat_amap:
            data_lats.append(lat)

    for i, lat in enumerate(data_lats):
        # 六邊形單雙行中心點不對齊
        if i % 2 == 0:
            for lng in data_lngs_first:
                cps.append([lng, lat])
        else:
            for lng in data_lngs_second:
                cps.append([lng, lat])
    for cp in cps:
        lng, lat = gcj02towgs84(cp[0], cp[1])
        point = Point((lng, lat))
        li = obj.filter(boundary__contains=point)
        if not li:
            continue
        # 左上角為第一個點
        result.append([[cp[0] - side_length * 0.5, cp[1] + lat_gap],
                       [cp[0] + side_length * 0.5, cp[1] + lat_gap],
                       [cp[0] + side_length, cp[1]],
                       [cp[0] + side_length * 0.5, cp[1] - lat_gap],
                       [cp[0] - side_length * 0.5, cp[1] - lat_gap],
                       [cp[0] - side_length, cp[1]]])

    return result


def generate_comb_mosaic_item(bc_id, lng_count):
    area = BusinessCircle.objects.filter(business_circle_id=bc_id)
    if not area:
        return ""

    mp = area[0].boundary
    if len(mp) > 1:
        return "返回的是多個多邊形，這樣不ok"

    boundary_li = [[list(list(location)) for location in list(polygon)] for polygon in list(mp.boundary.coords)]

    result = cal_rasterized_points_from_boundary(area, boundary_li, lng_count)
    # data = [[j, i] for i in range(lat_count) for j in range(lng_count)]
    data = []
    cursor = connection.cursor()
    SQL_string_list = []
    dp_result = copy.deepcopy(result)
    for index, points in enumerate(dp_result):
        MULTIPOLYGON_txt = ""
        points.append(points[0])
        for point in points:
            lng, lat = gcj02towgs84(point[0], point[1])
            MULTIPOLYGON_txt += ", {} {}".format(lng, lat)
        MULTIPOLYGON_txt = MULTIPOLYGON_txt[2:]
        SQL_RAW = """SELECT avg(forecast_price) FROM economic_dqchina_loupaninfo
        WHERE st_contains(st_geomfromtext('MULTIPOLYGON((({})))', 4326) :: geometry, center_point :: geometry) UNION ALL """.format(MULTIPOLYGON_txt)

        SQL_string_list.append(SQL_RAW)
    SQL_RAWS = "".join(SQL_string_list)
    cursor.execute(SQL_RAWS[:-10])
    rows = cursor.fetchall()
    for index, row in enumerate(rows):
        data.append(row[0])
    item = {
        "data": [float(d) if d else 0 for d in data],
        "result": result,
    }
    # print(item["data"])
    # print(item["result"])
    return item