text: # link: https://archive.ics.uci.edu/ml/datasets/Taxi+Service+Trajectory+-+Prediction+Challenge%2C+ECML+PKDD+2015
  import json
  import pandas as pd
  import os
  import numpy as np
  import time
  import math
  
  old_time_format = '%Y-%m-%d %H:%M:%S'
  new_time_format = '%Y-%m-%dT%H:%M:%SZ'
  
  # MIN_TIME = '2013-07-01 00:00:00'
  # MAX_TIME = '2014-07-01 00:00:00'
  # MIN_TIMESTAMP = float(
  #     datetime.timestamp(
  #         pd.to_datetime(MIN_TIME, utc=True, format=old_time_format)))
  # MAX_TIMESTAMP = float(
  #     datetime.timestamp(
  #         pd.to_datetime(MAX_TIME, utc=True, format=old_time_format)))
  
  
  def judge_id(value, dividing_points, equally=True):
      if equally:
          min_v = dividing_points[0]
          interval = dividing_points[1] - dividing_points[0]
          idx = int((value - min_v) / interval)
          max_id = len(dividing_points) - 2
          return min(max_id, idx)
      else:
          for i, num in enumerate(dividing_points):
              if value <= num:
                  return i - 1
          return len(dividing_points)
  
  
  def gen_config_geo():
      geo = {"including_types": [
          "Polygon"
      ],
          "Polygon": {
              "row_id": "num",
              "column_id": "num"
          }
      }
      return geo
  
  
  def gen_config_grid(row_num, column_num):
      grid = {
          "including_types": [
              "state"
          ],
          "state": {
              "row_id": row_num,
              "column_id": column_num,
              "inflow": "num",
              "outflow": "num"
          }
      }
      return grid
  
  
  def gen_config_info(file_name, interval):
      info = \
          {
              "data_col": [
                  "inflow",
                  "outflow"
              ],
              "data_files": [
                  file_name
              ],
              "geo_file": file_name,
              "output_dim": 2,
              "init_weight_inf_or_zero": "inf",
              "set_weight_link_or_dist": "dist",
              "calculate_weight_adj": False,
              "weight_adj_epsilon": 0.1,
              "time_intervals": interval,
          }
      return info
  
  
  def gen_config(output_dir_flow, file_name, row_num, column_num, interval):
      print('begin config')
      config = {}
      data = json.loads(json.dumps(config))
      data["geo"] = gen_config_geo()
      data["grid"] = gen_config_grid(row_num, column_num)
      data["info"] = gen_config_info(file_name, interval)
      # config = json.dumps(data)
      with open(output_dir_flow + "/config.json", "w") as f:
          json.dump(data, f, ensure_ascii=False, indent=1)
      # print(config)
      print('finish config')
  
  
  if __name__ == '__main__':
      start_time = time.time()
      # 参数
      # 时间间隔
      interval = 3600
      # 行数
      row_num = 20
      # 列数
      column_num = 10
      # 开始年月日
      (start_year, start_month, start_day) = (2013, 7, 1)
      # 结束年月日
      (end_year, end_month, end_day) = (2013, 9, 30)
      # 输入文件名称
      # input_file_name = 'Porto.csv'
      input_file_name = 'train.csv'
      # 输出文件名称
      file_name = 'PORTO%d%02d-%d%02d' % (start_year, start_month, end_year, end_month)
      # 输出文件夹名称
      output_dir_flow = 'output/PORTO%d%02d-%d%02d' % (start_year, start_month, end_year, end_month)
      # 输入文件夹名称
      input_dir_flow = 'input/Porto/'
  
      data_url = input_dir_flow + input_file_name
  
      # 创建输出文件夹
      if not os.path.exists(output_dir_flow):
          os.makedirs(output_dir_flow)
  
      # 读入文件
      data_set_Porto = pd.read_csv(
          data_url,
          usecols=[4, 5, 8],
          dtype={'TAXI_ID': np.int64, 'TIMESTAMP': np.int64}
      )  #
      data_set_Porto.reset_index(drop=True, inplace=True)
      data_set_Porto.columns = ['taxi_id', 'timestamp', 'polyline']
      print('finish reading csv')
      print(data_set_Porto.dtypes)
      print('shape:', data_set_Porto.shape)
  
      # 筛选时间范围（由6个参数决定）
      dt1 = '%d-%02d-%02dT00:00:00Z' % (start_year, start_month, start_day)
      dt2 = '%d-%02d-%02dT23:59:59Z' % (end_year, end_month, end_day)
      # print(dt1, dt2)
      stimes = time.mktime(time.strptime(dt1, '%Y-%m-%dT%H:%M:%SZ'))
      etimes = time.mktime(time.strptime(dt2, '%Y-%m-%dT%H:%M:%SZ'))
      # print(stimes, etimes)
      data_set_Porto = data_set_Porto[(data_set_Porto['timestamp'] <= etimes)
                                      & (data_set_Porto['timestamp'] >= stimes)]
      print('select time, shape:', data_set_Porto.shape)
      # # 去掉有空值列
      # data_set_Porto.dropna(subset=['taxi_id', 'timestamp', 'polyline'], inplace=True)
      # print('drop na, shape:', data_set_Porto.shape)
  
      # 生成时间序列
      print('begin to generate dividing points')
      # 时间戳的最小最大值，按interval对齐
      min_timestamp = stimes
      max_timestamp = etimes
      min_timestamp = float(math.floor(min_timestamp / interval) * interval)
      max_timestamp = float(math.ceil(max_timestamp / interval) * interval)
      # 以interval为颗粒度生成时间分割序列
      time_dividing_point = \
          list(np.arange(min_timestamp, max_timestamp, interval))
      # print(time_dividing_point)
      convert = []
      for t in time_dividing_point:
          convert.append(time.strftime('%Y-%m-%dT%H:%M:%SZ', time.localtime(t)))
  
      # 生成经纬度范围
      # 遍历polyline，维护经纬度的最大最小值
      LAT_MAX, LAT_MIN = float("-inf"), float("inf")
      LON_MAX, LON_MIN = float("-inf"), float("inf")
      for index, row in data_set_Porto.iterrows():
          poly_lst = eval(row['polyline'])
          for lon, lat in poly_lst:
              LON_MAX = max(LON_MAX, lon)
              LON_MIN = min(LON_MIN, lon)
              LAT_MAX = max(LAT_MAX, lat)
              LAT_MIN = min(LAT_MIN, lat)
  
      # 生成经纬度列表
      lat_diff = LAT_MAX - LAT_MIN
      lat_dividing_points = \
          [round(LAT_MIN + lat_diff / row_num * i, 3)
           for i in range(row_num + 1)]
      lon_diff = LON_MAX - LON_MIN
      lon_dividing_points = \
          [round(LON_MIN + lon_diff / column_num * i, 3)
           for i in range(column_num + 1)]
  
      # 生成 grid
      # generate gird data (.geo)
      print('begin to generate geo')
      geo_data = pd.DataFrame(
          columns=['geo_id', 'type', 'coordinates', 'row_id', 'column_id'])
      for i in range(row_num):
          for j in range(column_num):
              index = i * column_num + j
              coordinates = [[
                  [lon_dividing_points[j], lat_dividing_points[i]],
                  [lon_dividing_points[j + 1], lat_dividing_points[i]],
                  [lon_dividing_points[j + 1], lat_dividing_points[i + 1]],
                  [lon_dividing_points[j], lat_dividing_points[i + 1]],
                  [lon_dividing_points[j], lat_dividing_points[i]]
              ]]  # list of list of [lon, lat]
              geo_data.loc[index] = [index, 'Polygon', coordinates, i, j]
      geo_data.to_csv(output_dir_flow + file_name + '.geo', index=False)
      print('finish geo')
  
      print('size of grid data:', len(time_dividing_point), len(lon_dividing_points) - 1, len(lat_dividing_points) - 1, 2)
  
      # 存储结果的数据结构，最后一维代表inflow 和 outflow
      grid_data = np.zeros((len(time_dividing_point), len(lon_dividing_points) - 1, len(lat_dividing_points) - 1, 2))
      print("begin to calculate grid")
      # 逐行统计in/out_flow
      for index, row in data_set_Porto.iterrows():
          # 初始化为第一个轨迹点的数据
          trajectory_list = eval(row['polyline'])
          if len(trajectory_list) == 0:
              continue
          timestamp = row['timestamp']
          time_index = judge_id(timestamp, time_dividing_point)
          trajectory0 = trajectory_list[0]
          lon_index = judge_id(trajectory0[0], lon_dividing_points)
          lat_index = judge_id(trajectory0[1], lat_dividing_points)
          # 对于一辆taxi的轨迹数据
          for lon, lat in trajectory_list[1:]:
              timestamp += 15 * 1000  # 15s
              if timestamp > etimes:
                  break
              time_index_new = judge_id(timestamp, time_dividing_point)
              lon_index_new = judge_id(lon, lon_dividing_points)
              lat_index_new = judge_id(lat, lat_dividing_points)
              if lon_index == lon_index_new and lat_index == lat_index_new:
                  time_index = time_index_new
              else:
                  # out += 1
                  grid_data[time_index][lon_index][lat_index][1] = grid_data[time_index][lon_index][lat_index][1] + 1
                  # in += 1
                  grid_data[time_index_new][lon_index_new][lat_index_new][0] =\
                      grid_data[time_index_new][lon_index_new][lat_index_new][0] + 1
                  time_index = time_index_new
                  lon_index = lon_index_new
                  lat_index = lat_index_new
  
      # grid_data.to_csv('flow_data.csv')
      print('finish calculating grid')
  
      print('begin to write .grid')
      grid_file = open(output_dir_flow + file_name + '.grid', 'w')
      # 'dyna_id', 'type', 'time', 'row_id', 'column_id', 'inflow', 'outflow'
      grid_file.write('dyna_id,' + 'type,' + 'time,' + 'row_id,' + 'column_id,' + 'inflow,' + 'outflow' + '\n')
      dyna_id = 0
      for lat_index in range(grid_data.shape[2]):
          for lon_index in range(grid_data.shape[1]):
              for time_index in range(grid_data.shape[0]):
                  times = convert[time_index]
                  grid_file.write(str(dyna_id) + ',' + 'state' + ',' + str(times)
                                  + ',' + str(lat_index) + ',' + str(lon_index)
                                  + ',' + str(grid_data[time_index][lon_index][lat_index][0]) +
                                  ',' + str(grid_data[time_index][lon_index][lat_index][1]) + '\n')
                  dyna_id = dyna_id + 1
                  if dyna_id % 10000 == 0:
                      print(str(dyna_id) + '//' + str(grid_data.shape[0] * grid_data.shape[1] * grid_data.shape[2]))
      grid_file.close()
      print('finish grid file')
  
      # 生成config.json文件
      gen_config(output_dir_flow, file_name, row_num, column_num, interval)
      print('finish')
      end_time = time.time()
      print(end_time - start_time)
  