clc,clear
[A,B] = xlsread('CityTags.xlsx');
[Geo,GeoName] = xlsread('Region.xlsx');
Geo = Geo(:,1:4);
num = length(A);
pre = B(1);
x = [];
y = [];
idx = [];
C = [];
D = [];
Region = [];
Tags = [];
Range = [];
% 单词方向 ：水平 = 1，竖直 = 0
Dir = [];
% 第 cnt 个聚类
cnt = 1;
% 为输入信息 打上聚类标签
cluster_points = [];
for i = 1:num
    cur = B(i);
    if ~isequal(pre,cur)
        X = [x' y'];
        Region = unique(Region);
        Rsz = size(Region);
        Region = [];
        K = Rsz(2);
        [idx,D] = kmeans(X,K);
        for cluster = 1:K
            CurCluster = [X(idx==cluster,1) X(idx==cluster,2)];
            CurCluster_Size = size(CurCluster);
            tmp = zeros(CurCluster_Size(1),1);
            tmp = tmp + cnt;
            cluster_points = [cluster_points;X(idx==cluster,1) X(idx==cluster,2) tmp];
            direction = abs(PCA_Rotation(CurCluster'));
            if direction > 45
                Dir = [Dir;0];
            else
                Dir = [Dir;1];
            end
            sz = sum(idx==cluster);
            tmpX = sort(X(idx==cluster,:));
            min_x = tmpX(ceil(0.05*sz),1);
            max_x = tmpX(ceil(0.95*sz),1);
            min_y = tmpX(ceil(0.05*sz),2);
            max_y = tmpX(ceil(0.95*sz),2);
            Range = [Range;min_x,max_x,min_y,max_y];
            C = [C;D(cluster,1) D(cluster,2) sz cnt];
            cnt = cnt + 1;
            x = [];
            y = [];
            Tags = [Tags;string(pre)];
        end
    end
    x = [x A(i,1)];
    y = [y A(i,2)];
    Region = [Region A(i,3)];
    pre=cur;
end
X = [x' y'];
Region = unique(Region);
Rsz = size(Region);
K = Rsz(2);
[idx,D] = kmeans(X,K);
for cluster = 1:K
    CurCluster = [X(idx==cluster,1) X(idx==cluster,2)];
    CurCluster_Size = size(CurCluster);
    tmp = zeros(CurCluster_Size(1),1);
    tmp = tmp + cnt;
    cluster_points = [cluster_points;X(idx==cluster,1) X(idx==cluster,2) tmp];
    direction = abs(PCA_Rotation(CurCluster'));
    if direction > 45
        Dir = [Dir;0];
    else
        Dir = [Dir;1];
    end
    sz = sum(idx==cluster);
    Range = [Range;min(X(idx==cluster,1)),max(X(idx==cluster,1)),min(X(idx==cluster,2)),max(X(idx==cluster,2))];
    C = [C;D(cluster,1) D(cluster,2) sz cnt];
    cnt = cnt + 1;
    x = [];
    y = [];
    Tags = [Tags;string(cur)];
end
Map = [];
GeoSize = size(Geo);
InfoSize = size(C);
for i = 1:InfoSize(1)
    for j = 1:GeoSize(1)
        cx = C(i,1);
        cy = C(i,2);
        if cx >= Geo(j,1) && cx <= Geo(j,2) && cy >= Geo(j,3) && cy <= Geo(j,4)
            Map = [Map;string(GeoName(j))];
            break;
        end
    end
end
xlswrite('CityTagsCluster1.xlsx', [Tags,C,Dir Map,Range], 'A1')
xlswrite('CityTagsWithIndex1.xlsx',cluster_points,'A1')

function [Dir,res] = PCA(Data)
% 以[-pi/2,pi/2]之间的角度pi表示方向
% 1.去除均值
[buf K] = size(Data);
miu = mean(Data')';
for k=1:K
    Data(:,k) = Data(:,k)-miu; 
end;
sigma = zeros(2,2);
% 2.计算协方差
for k=1:K
    x = Data(:,k);
    sigma = sigma+x*x';
end;
sigma = sigma/K;  
% 3.特征分解
[V,D] = eig(sigma);
if (D(1,1)<D(2,2)) % 把较大的特征值对应的向量挪到第一行
     buf = V(1,:);
     V(1,:) = V(2,:);
     V(2,:) = buf;
end;
% 4.求解主方向向量
res = inv(V)*[1;0];
Dir = atan(res(2)/res(1))/pi*180;

# Load Geo Map
dataf = pd.read_excel('GeoRegion.xlsx')
wordmap = dataf.values
sz = wordmap.shape
geomap = dict()
for i in range(0,sz[0]):
    geomap[wordmap[i][0]] =[wordmap[i][1],wordmap[i][2],wordmap[i][3],wordmap[i][4]]

# Load Data Set
# Load Cluster Info
df = pd.read_excel('CityTagsClusterVer4.xlsx')
words = df.values
# Load Tags labeled by cluster index
df2 = pd.read_excel('CityTagsWithIndex2.xlsx')
points = df2.values
cluster_size = words.shape[0] + 10
cluster_points = [[] for i in range(cluster_size)]

for i in range(0, points.shape[0]):
    cluster_points[points[i][2]].append((points[i][0],points[i][1]))

class WordCluster(object):
    def __init__(self,word,x,y,num,direction,region,min_x,max_x,min_y,max_y,idx):
        self.word = word
        self.x = x
        self.y = y
        self.num = num
        self.direction = direction
        self.region = region
        self.min_x = min_x
        self.max_x = max_x
        self.min_y = min_y
        self.max_y = max_y
        self.idx = idx

    def __lt__(self, other):
        return self

    def __lt__(self, other):
        return self.num > other.num

    def __str__(self):
        return '(' + self.word + ', ' + str(self.x) + ', ' + str(self.y) + ', ' + str(self.num) + ', '+str(self.direction)+')'

# compute the area of M in pixel : Am
img = cv2.imread("SD.png")
gray = cv2.cvtColor(img, cv2.COLOR_BGRA2GRAY)
h, w = gray.shape[:2]
m = np.reshape(gray, [1, w * h])
mean = m.sum() / (w * h)
ret, binary = cv2.threshold(gray, mean, 255, cv2.THRESH_BINARY)
Am = len(binary[binary == 0])

font_size = int((cmath.sqrt((freq * Am)/n).real))

def sample_position_v2(self, size_x, size_y, b_x, b_y, e_x, e_y, center_x, center_y):
    dis_to_center = 99999999
    ans_x = 0
    ans_y = 0
    for y in range(max(1, b_y), min(e_y - size_y, self.height - size_y)):
        for x in range(max(1,b_x), min(e_x - size_x,self.width - size_x)):
            area = self.integral[y - 1, x - 1] + self.integral[y + size_y - 1, x + size_x - 1]
            area -= self.integral[y - 1, x + size_x - 1] + self.integral[y + size_y - 1, x - 1]
            if not area:
                tmpdis = Euclidean_distance(x + size_x/2 , center_x, y + size_y/2, center_y).real
                if tmpdis < dis_to_center:
                    dis_to_center = tmpdis
                    ans_x = x
                    ans_y = y
    if dis_to_center == 99999999:
        return None
    else:
        return ans_y, ans_x

weight_aver = 0.5
penalty = 9999999999.9
scaling_factor = 1.0
# find the best placement
b_x = center_x - (box_size[0] + self.margin) // 2
b_y = center_y - (box_size[1] + self.margin) // 2
for s_f in np.arange(0.05,1.00,0.05):
    tmp_penalty = weight_aver * (1 - s_f)
    new_x = (int)(center_x - box_size[0] * s_f / 2)
    new_y = (int)(center_y - box_size[1] * s_f / 2)
    distance_error = cmath.sqrt((b_x - new_x ) * (b_x - new_x ) + (b_y - new_y) * (b_y - new_y))
    tmp_penalty = tmp_penalty + (1 - weight_aver) * distance_error / cmath.sqrt(A)
    if tmp_penalty <= penalty:
        scaling_factor = s_f
        penalty = tmp_penalty
freq = (int)(freq * scaling_factor)
que.put(WordCluster(word, center_x, center_y , freq, curword.direction, CurRegion, min_x, max_x, min_y, max_y, curword.idx))

def Color_Distribution(self, i):
    Hue_Spectrum = ['rgb(29, 122, 202)', 'rgb(0, 139, 182)', 'rgb(1, 146, 149)', 'rgb(0, 146, 107)',
                    'rgb(1, 142, 46)', 'rgb(78, 134, 0)', 'rgb(131, 123, 0)', 'rgb(161, 111, 0)',
                    'rgb(182, 94, 46)', 'rgb(197, 78, 110)', 'rgb(197, 69, 152)', 'rgb(180, 77, 186)',
                    'rgb(135, 98, 204)']
    l = len(Hue_Spectrum)
    c = 3
    return Hue_Spectrum[(i + c) % l]

def GeoGenerate(self, words, A , Am ,geomap, cluster_points, PointsNum, max_font_size=None):
     '''
     :param words: Preprocessed data sets
     :param A: the area of the rectangle bounding region M
     :param Am: the area of M in pixel
     :param geomap : the region coordinate
     :param cluster_points : the points that labeled with the cluster index
     :param max_font_size: Max Font Size
     :return: Geo Word Cloud

     '''
     # Coverage error
     Measure1 = 0.0
     # Words not represent
     Measure2 = 0.0
     # difference between M and M'
     Measure3 = 0.0
     MM = np.zeros((800,500))
     # the total frequency
     n = 0
     que = Q.PriorityQueue()
     min_freq = 999999999
     for i in range(0, words.shape[0], 1):
         que.put(WordCluster(words[i][0], words[i][1], words[i][2], words[i][3], words[i][4], words[i][5],words[i][6],words[i][7],words[i][8],words[i][9], i+1))
         n = n + words[i][3]
         if min_freq > words[i][3]:
             min_freq = words[i][3]

     # the number of the cluster
     CountWord  = que.qsize()

     '''
     #users can diy the fill word to complete the shape of the map
     
     if que.qsize() < self.max_words:
         for i in range(0, self.max_words - que.qsize(), 1):
             # users can diy fill words
             que.put(WordCluster('fill', self.width//2, self.height//2, (int)(min_freq * 0.5) , 1))
     '''

     frequencies = list()

     if self.random_state is not None:
         random_state = self.random_state
     else:
         random_state = Random()

     # set the mask
     if self.mask is not None:
         boolean_mask = self._get_bolean_mask(self.mask)
         width = self.mask.shape[1]
         height = self.mask.shape[0]
     else:
         boolean_mask = None
         height, width = self.height, self.width

     occupancy = IntegralOccupancyMap(height, width, boolean_mask)

     # create image
     img_grey = Image.new("L", (width, height))
     draw = ImageDraw.Draw(img_grey)
     img_array = np.asarray(img_grey)
     font_sizes, positions, orientations, colors = [], [], [], []

     last_freq = 1.

     if max_font_size is None:
         # if not provided use default font_size
         max_font_size = self.max_font_size

     curidx = 0

     # start drawing grey image
     while not que.empty():
         curword = que.get()
         freq = curword.num
         word = curword.word
         CurRegion = curword.region
         # the cluster center
         center_x = (int)(curword.x)
         center_y = (int)(curword.y)
         # the cluster range
         min_x = (int)(curword.min_x)
         max_x = (int)(curword.max_x)
         min_y = (int)(curword.min_y)
         max_y = (int)(curword.max_y)

         # select the font size
         rs = self.relative_scaling

         if rs != 0:
             font_size = int((cmath.sqrt((freq * Am)/n).real))

         # select the rotation
         if curword.direction == 1:
             orientation = None
         else:
             orientation = Image.ROTATE_90

         tried_other_orientation = False
         #Origin_orientation = orientation

         while True:
             # try to find a position
             font = ImageFont.truetype(self.font_path, font_size)
             # transpose font optionally
             transposed_font = ImageFont.TransposedFont(
                 font, orientation=orientation)
             # get size of resulting text
             box_size = draw.textsize(word, font=transposed_font)
             # find possible places using integral image:
             '''
             if freq < min_freq:
                 result = occupancy.sample_position_v2(box_size[0] + self.margin, box_size[1] + self.margin,
                                                       1, 1,self.width ,self.height )
             else:
             '''
             '''
             b_x = geomap[CurRegion][0]
             b_y = geomap[CurRegion][2]
             e_x = geomap[CurRegion][1]
             e_y = geomap[CurRegion][3]
             '''

             result = occupancy.sample_position_v2(box_size[0] + self.margin, box_size[1] + self.margin,
                                                  min_x, min_y, max_x, max_y,center_x,center_y)

             if result is not None or font_size < self.min_font_size:
                 # either we found a place or font-size went too small
                 break
             # if we didn't find a place, make font smaller
             # but first try to rotate
             if not tried_other_orientation :
                 orientation = (Image.ROTATE_90 if orientation is None else Image.ROTATE_270)
                 tried_other_orientation = True
                 continue
             # scale the word down
             else:
                 weight_aver = 0.5
                 penalty = 9999999999.9
                 scaling_factor = 1.0
                 # find the best placement
                 b_x = center_x - (box_size[0] + self.margin) // 2
                 b_y = center_y - (box_size[1] + self.margin) // 2
                 for s_f in np.arange(0.05,1.00,0.05):
                     tmp_penalty = weight_aver * (1 - s_f)
                     new_x = (int)(center_x - box_size[0] * s_f / 2)
                     new_y = (int)(center_y - box_size[1] * s_f / 2)
                     distance_error = cmath.sqrt((b_x - new_x ) * (b_x - new_x ) + (b_y - new_y) * (b_y - new_y))
                     tmp_penalty = tmp_penalty + (1 - weight_aver) * distance_error / cmath.sqrt(A)
                     if tmp_penalty <= penalty:
                         scaling_factor = s_f
                         penalty = tmp_penalty
                 freq = (int)(freq * scaling_factor)
                 que.put(WordCluster(word, center_x, center_y , freq, curword.direction, CurRegion, min_x, max_x, min_y, max_y, curword.idx))
                 break

         # can place
         if result is not None:
             x, y = np.array(result) + self.margin // 2
             # actually draw the text
             draw.text((y, x), word, fill="white", font=transposed_font)
             tmpmeasure = Hausdorff_Distance(y,x,box_size[0] + self.margin,box_size[1] + self.margin,curword.idx,cluster_points)
             tmp1, tmp2 = np.array(tmpmeasure)
             Measure1 = Measure1 + tmp1
             Measure2 = Measure2 + tmp2
             Fill_Map(MM, y, x,box_size[0] + self.margin, box_size[1] + self.margin)
             # append attributes
             frequencies.append((word, freq))
             positions.append((x, y))
             orientations.append(orientation)
             font_sizes.append(font_size)
             '''
             colors.append(self.color_func(word, font_size=font_size,
                                           position=(x, y),
                                           orientation=orientation,
                                           random_state=random_state,
                                           font_path=self.font_path))
             '''
             colors.append(self.Color_Distribution(curidx))
             curidx = curidx + 1
             # print(colors)
             # recompute integral image
             if self.mask is None:
                 img_array = np.asarray(img_grey)
             else:
                 img_array = np.asarray(img_grey) + boolean_mask
             # recompute bottom right
             # the order of the cumsum's is important for speed ?!
             occupancy.update(img_array, x, y)
             last_freq = freq

     # layout
     self.layout_ = list(zip(frequencies, font_sizes, positions,
                             orientations, colors))


     # Three Measures
     Measure1 = (Measure1/PointsNum).real
     Measure2 = (PointsNum - Measure2)/PointsNum
     Measure3 = (162384 - len(MM[MM == 1]))/162384
     print ("Measure1:")
     print (Measure1)
     print ("Measure2:")
     print (Measure2)
     print("Measure3:")
     print(Measure3)

     return self

def Euclidean_distance(x1, x2, y1, y2):
    return cmath.sqrt((x2 - x1) * (x2 - x1)  + (y2 - y1) * (y2 - y1))

def Hausdorff_Distance(x, y, w, h, idx, cluster_points):
    sz = len(cluster_points[idx])
    min_x = x
    max_x = x + w
    min_y = y
    max_y = y + h
    dis = 0
    num = 0
    for i in range(0, sz):
        cur_x = cluster_points[idx][i][0]
        cur_y = cluster_points[idx][i][1]
        if cur_x >= min_x and cur_x <= max_x and cur_y >= min_y and cur_y <= max_y:
            num = num + 1
            continue
        mindis = 99999999
        for x2 in range(min_x,max_x + 1):
            tmpdis1 = Euclidean_distance(cur_x, x2, cur_y, min_y).real
            tmpdis2 = Euclidean_distance(cur_x, x2, cur_y, max_y).real
            tmpdis = min(tmpdis1, tmpdis2)
            if mindis > tmpdis:
                mindis = tmpdis
        for y2 in range(min_y, max_y + 1):
            tmpdis1 = Euclidean_distance(cur_x, min_x, cur_y, y2).real
            tmpdis2 = Euclidean_distance(cur_x, max_x, cur_y, y2).real
            tmpdis = min(tmpdis1, tmpdis2)
            if mindis > tmpdis:
                mindis = tmpdis
        dis = dis + mindis
    print((dis,num))
    return dis,num

def Fill_Map(MM, x, y, size_x, size_y):
    for i in range(x - 1, x+size_x):
        for j in range(y - 1, y+size_y):
            MM[i][j] = 1
           
.
.
.
# After placing one word
tmpmeasure = Hausdorff_Distance(y,x,box_size[0] + self.margin,box_size[1] + self.margin,curword.idx,cluster_points)
tmp1, tmp2 = np.array(tmpmeasure)
Measure1 = Measure1 + tmp1
Measure2 = Measure2 + tmp2
Fill_Map(MM, y, x,box_size[0] + self.margin, box_size[1] + self.margin)
.
.
.
# At the end of the algorithm 
Measure1 = (Measure1/PointsNum).real
Measure2 = (PointsNum - Measure2)/PointsNum
Measure3 = (162384 - len(MM[MM == 1]))/162384