convert python code to C++ code:   import numpy as np
import torch
import torch
from torch_scatter import scatter
import time
from matplotlib import pyplot as plt
import os
from dv import AedatFile
import scipy
import cv2

N_CHANNELS = 12 # 不同 channels 对应不同 representations 

def create_window(event, stacking_type="SBN"):
    # x = event[:, 0]
    # y = event[:, 1]
    # t = event[:, 2]
    # p = event[:, 3]
    windows = []
    num_of_events = event.shape[0] // 3   
    # 应该在这个文件里改
    equispaced_factor = 1 / 3
    cur_num_of_events = len(event[:, 2])# len(event[0]) # 这块可能有问题
    windows.append(event) # 0
    # windows.append((x, y, p, t))
# 稍等
    if stacking_type == "SBN":
        for i in range(3):
            event_=event[i * num_of_events: (i + 1) * num_of_events,:]  # [0:3, :] [3:6, :] [6:9, :]
            windows.append(event_.copy()) # 1 2 3 
            # x_ = x[i * num_of_events : (i + 1) * num_of_events]
            # y_ = y[i * num_of_events : (i + 1) * num_of_events]
            # p_ = p[i * num_of_events : (i + 1) * num_of_events]
            # t_ = t[i * num_of_events : (i + 1) * num_of_events]
            # windows.append((x_, y_, p_, t_))

        # [1/2,1],[3/4,1],[7/8,1]
        for _ in range(3):
            cur_num_of_events = cur_num_of_events // 2
            event=event[cur_num_of_events:,:]
            # cur_num_of_events = cur_num_of_events // 2
            # x = x[cur_num_of_events:]
            # y = y[cur_num_of_events:]
            # p = p[cur_num_of_events:]
            # t = t[cur_num_of_events:]
            # windows.append((x, y, p, t))
            windows.append(event.copy()) # 4 5 6

    elif stacking_type == "SBT":
        for i in range(3):
            event_=event[(event[:,2]<=(i + 1) * equispaced_factor)&(event[:,2]>= i * equispaced_factor)]
            windows.append(event_.copy())

        factor = 1
        for _ in range(4):
            factor = factor / 2
            event=event[(event[:,2]<=factor)]
            windows.append(event.copy())

    return windows


def run(src, index,aggregation,height,width):
    if aggregation == "variance":
        event_surface = scatter(
            src=src, index=index, dim_size=height * width, reduce="mean"
        )
        event_surface_squared = scatter(
            src=src ** 2,
            index=index,
            dim_size= height * width,
            reduce="mean",
        )
        result = event_surface_squared - event_surface ** 2
        event_surface = result
    else:
        event_surface = scatter(
            src=src,
            index=index,
            dim_size= height * width,
            reduce= aggregation,
        )
    event_surface = event_surface.reshape(height, width)

    return event_surface

def get_processed_stack_event(window,func,aggregation,height,width):
    index=torch.tensor(window[:,0]+window[:,1]*width,dtype=torch.int64)
    # print("event:{}".format(window[936]))
    # print("index:{}".format(np.where(index>684800)))

    if  func== "timestamp":
        event_surface=run(torch.tensor(window[:,2]), index,aggregation=aggregation,height=height,width=width)
    elif func == "polarity":
        event_surface = run(torch.tensor(window[:,3]), index,aggregation=aggregation,height=height,width=width)
    elif func == "count":
        event_surface = run(
            torch.tensor(window).new_ones([window.shape[0]]), index, aggregation=aggregation,height=height,width=width 
        )
    elif func == "timestamp_pos":
        positive_events = window[window[:, 3] == 1]
        positive_index = torch.tensor(
            positive_events[:, 0] + positive_events[:, 1] * width
        ).to(torch.int64)
        event_surface = run(
            torch.tensor(positive_events[:, 2]), positive_index, aggregation=aggregation,height=height,width=width
        )
    elif func == "timestamp_neg":
        negative_events = window[window[:, 3] == -1]
        if len(list(negative_events)) == 0:
            negative_events = window[window[:, 3] == 0]
        negative_index = torch.tensor(
            negative_events[:, 0] + negative_events[:,1] * width
        ).to(torch.int64)
        event_surface = run(
            torch.tensor(negative_events[:, 2]), negative_index, aggregation=aggregation,height=height,width=width
        )
    elif func == "count_pos":
        positive_events = window[window[:, 3] == 1]
        positive_index = torch.tensor(
            positive_events[:, 0] + positive_events[:, 1] * width
        ).to(torch.int64)
        event_surface = run(
            torch.tensor(positive_events).new_ones([positive_events.shape[0]]),
            positive_index, aggregation=aggregation,height=height,width=width
        )
    elif func == "count_neg":
        negative_events = window[window[:, 3] == -1]
        if len(list(negative_events)) == 0:
            negative_events = window[window[:, 3] == 0]
        negative_index = torch.tensor(
            negative_events[:, 0] + negative_events[:, 1] * width
        ).to(torch.int64)
        event_surface = run(
            torch.tensor(negative_events).new_ones([negative_events.shape[0]]),
            negative_index, aggregation=aggregation,height=height,width=width
        )

    return event_surface

def get_optimized_representation(event,height,width,num_bins):
    
    
    torch.cuda.synchronize()
    start_time = time.time()
    # 和timestamp相关的可视化有问题 
  
    window_indexes = [0, 3, 2, 6, 5, 6, 2, 5, 1, 0, 4, 1]
    functions = [
        "polarity",         # 0
        "timestamp_neg",    # 1 
        "count_neg",        # 2
        "polarity",         # 3
        "count_pos",        # 4
        "count",            # 5
        "timestamp_pos",    # 6
        "count_neg",        # 7  
        "timestamp_neg",    # 8
        "timestamp_pos",    # 9
        "timestamp",        # 10
        "count",            # 11
    ]
    aggregations = [
        "variance",
        "variance",
        "mean",
        "sum",
        "mean",
        "sum",
        "mean",
        "mean",
        "max",
        "max",
        "max",
        "mean",
    ]

    event[:,2] = event[:,2]-event[0,2]

    delta_t=event[-1,2]-event[0,2]
    if delta_t!=0:
        event[:,2]=event[:,2]/delta_t # 时间归一化 0 - 1.0
    else:
        event[:,2] = event[:,2]

    windows=create_window(event)

    #你看这个，就是要表征的图，是不分极性的，只能背景颜色
    # 我给你看一下论文 中的图
    representation=torch.zeros([num_bins,height,width]) 
    

    for i in range(len(aggregations)):
        stacked_event=get_processed_stack_event(window=windows[window_indexes[i]],func=functions[i],
                                                aggregation=aggregations[i],height=height,width=width)
        representation[i]=stacked_event
    
    
    torch.cuda.synchronize()
    end_time = time.time()

    process_time = end_time - start_time
    # print("representation:{}".format(representation.shape))

    return process_time, representation


# main functions 

if __name__ == '__main__':
    save_voxel = 5000
    use_mode = 'frame_interval_time'
    device = torch.device("cuda:0")
    data_path = r"/data2/local_userdata/heyuting/Projects1/COESOT-main/CEUTrack/data/COESOT/train/"

    save_path = r"/data2/local_userdata/heyuting/Projects1/COESOT-main/CEUTrack/data/COESOT/train/"
    video_files = os.listdir(data_path)
    dvs_img_interval = 1
    for videoID in range(len(video_files)):
        foldName = video_files[videoID]
        if foldName not in ['dvSave-2022_02_07_22_19_02']:
            continue

        print("==>> enter in: ", foldName)
        fileLIST = os.listdir(os.path.join(data_path, foldName))
        mat_save = os.path.join(save_path, foldName, foldName+'_voxel/')

        if not os.path.exists(os.path.join(save_path,foldName,foldName+'_OPR')):
            os.mkdir(os.path.join(save_path,foldName,foldName+'_OPR'))
        new_voxel_save=os.path.join(save_path,foldName,foldName+'_OPR/')
        read_aedat4_path=os.path.join(data_path,foldName,foldName +'.aedat4')
        frame_exposure_time = []
        frame_interval_time = []
        with AedatFile(read_aedat4_path) as f:
            for frame in f['frames']:
                if use_mode=='frame_exposure_time':
                    frame_exposure_time.append([frame.timestamp_start_of_exposure,frame.timestamp_end_of_exposure])
                elif use_mode=='frame_interval_time':
                    frame_interval_time.append([frame.timestamp_start_of_frame,
                                            frame.timestamp_end_of_frame]) 
                use_mode=='frame_interval_time'
        frame_timestamp = frame_interval_time
        frame_num = len(frame_timestamp)
        print(frame_num)
        events = np.hstack([packet for packet in f['events'].numpy()])
        cost_time = 0.0
        for frame_no in range(0, int(frame_num / dvs_img_interval) - 1):
            # start_time = time.time() 
            start_idx = np.where(events['timestamp'] >= frame_timestamp[frame_no][0])[0][0]
            end_idx = np.where(events['timestamp'] >= frame_timestamp[frame_no][1])[0][0]
            sub_event = events[start_idx:end_idx]
            idx_length = end_idx - start_idx
            t = events[start_idx:end_idx]['timestamp']
            x_all = events[start_idx:end_idx]['x']
            y_all = events[start_idx:end_idx]['y']
            p_all = events[start_idx:end_idx]['polarity']
            if start_idx == end_idx:
                time_length = 0
                # t = ((t-t) / time_length) * 1000
                scipy.io.savemat(mat_save+'frame{:04d}.mat'.format(frame_no),mdict={'coor':np.zeros([100,3]),
                                'features':np.zeros([100,16])}) 
                print('empty event frame:', frame_no)
                continue
            else:
                time_length = t[-1] - t[0]
                # t = ((t-t[0]) / time_length) * 1000
            all_idx = np.where(sub_event['polarity'] != -1)      # all event
            neg_idx = np.where(sub_event['polarity'] == 0)      # neg event
            t = t[all_idx]
            x = x_all[all_idx]
            y = y_all[all_idx]
            p = p_all[all_idx]
            # t = t
            # x = x_all
            # y = y_all
            # p = p_all
            p[neg_idx] = -1    
            num_events=end_idx-start_idx
            num_bins = 12 
            height=346
            width=260
            events_tensor = np.stack([x, y, t, p], axis=-1)
            grid_voxel, re = get_optimized_representation(events_tensor, 260, 346, num_bins) 
            cost_time = cost_time + grid_voxel 
            img = re.numpy()
            
            fig, axes = plt.subplots(nrows=3, ncols=4, figsize=(12, 9))

            # 在每个子图中画一些示例内容
            img_idx = 0
            for i in range(3):
                for j in range(4):
                    # 在每个子图中绘制一些内容，这里以文本标注子图位置
                    tmp_img = img[img_idx, :, :]
                    im = axes[i, j].imshow(tmp_img, cmap = 'inferno')
                    img_idx = img_idx + 1
            plt.tight_layout()
            fig.suptitle('Event Frame')
            # plt.imshow(img, cmap='inferno')
            plt.xlabel("x [pixels]")
            plt.ylabel("y [pixels]")
            plt.colorbar(im)
            # print(mat_save)
            plt.savefig(new_voxel_save+'new_voxel{:04d}.png'.format(frame_no))

        cost_time = cost_time / frame_num
        print('mean cost time',str(cost_time))