def normal(*vecs):
    # 3~4個のVectorのNormalを求める
    if len(vecs) == 3:
        return geom.normal(*vecs)
    elif len(vecs) == 4:
        n1 = geom.normal(vecs[0], vecs[1], vecs[3])
        n2 = geom.normal(vecs[1], vecs[2], vecs[3])
        if n1.dot(n2) < 0:
            n2.negate()
        return (n1 + n2).normalized()

def get_color_from_normal(dvk, pol, num_verts, vectorlight, colo):
    if num_verts <= 4:
        normal_no = normal(dvk[pol[0]], dvk[pol[1]], dvk[pol[2]])
    else:
        normal_no = normal(dvk[pol[0]], dvk[pol[1]], dvk[pol[2]], dvk[pol[3]])

    normal_no = (normal_no.angle(vectorlight, 0)) / pi

    r = (normal_no * colo[0]) - 0.1
    g = (normal_no * colo[1]) - 0.1
    b = (normal_no * colo[2]) - 0.1
    return (r+0.2, g+0.2, b+0.2)

 def __init__(self, v1, v2, v3):
     self.verts = [v1, v2, v3]
     self.normal = geom.normal(v1.co, v2.co, v3.co)
     self.edge_keys = [tuple(sorted((self.verts[i - 1], self.verts[i]),
                                    key=lambda v: v.index))
                       for i in range(3)]
     self.outer_verts = []

def generate_3PT(pts, obj, nv, mode=1):
    mw = obj.matrix_world
    V = Vector
    nv = max(3, nv)

    # construction
    v1, v2, v3, v4 = V(pts[0]), V(pts[1]), V(pts[1]), V(pts[2])
    edge1_mid = v1.lerp(v2, 0.5)
    edge2_mid = v3.lerp(v4, 0.5)
    axis = geometry.normal(v1, v2, v4)
    mat_rot = mathutils.Matrix.Rotation(math.radians(90.0), 4, axis)

    # triangle edges
    v1_ = ((v1 - edge1_mid) * mat_rot) + edge1_mid
    v2_ = ((v2 - edge1_mid) * mat_rot) + edge1_mid
    v3_ = ((v3 - edge2_mid) * mat_rot) + edge2_mid
    v4_ = ((v4 - edge2_mid) * mat_rot) + edge2_mid

    r = geometry.intersect_line_line(v1_, v2_, v3_, v4_)
    if r:
        p1, _ = r
        cp = mw * p1
        bpy.context.scene.cursor_location = cp

        if mode == 0:
            pass

        elif mode == 1:
            generate_bmesh_repr(p1, v1, axis, nv)

    else:
        print('not on a circle')

def _binary_write(filepath, faces):
    with open(filepath, 'wb') as data:
        fw = data.write
        # header
        # we write padding at header beginning to avoid to
        # call len(list(faces)) which may be expensive
        fw(struct.calcsize('<80sI') * b'\0')

        # 3 vertex == 9f
        pack = struct.Struct('<9f').pack

        # number of vertices written
        nb = 0

        for face in faces:
            # calculate face normal
            # write normal + vertexes + pad as attributes
            fw(struct.pack('<3f', *normal(*face)) + pack(*itertools.chain.from_iterable(face)))
            # attribute byte count (unused)
            fw(b'\0\0') 
            nb += 1

        # header, with correct value now
        data.seek(0)
        fw(struct.pack('<80sI', _header_version().encode('ascii'), nb))

 def calc_vert_normal(self, vert, looptris, fallback=Vector((0, 0, 0))):
     normal = Vector()
     num = 0
     for tri in looptris:
         normal += geom.normal(*[loop.vert.co for loop in tri])
         num += 1
     return normal / num

def generate_3PT_mode_1_(pts, obj):
    origin = obj.location
    transform_matrix = obj.matrix_local
    V = Vector

    # construction
    v1, v2, v3, v4 = V(pts[0]), V(pts[1]), V(pts[1]), V(pts[2])
    edge1_mid = v1.lerp(v2, 0.5)
    edge2_mid = v3.lerp(v4, 0.5)
    axis = geometry.normal(v1, v2, v4)
    mat_rot = mathutils.Matrix.Rotation(math.radians(90.0), 4, axis)

    # triangle edges
    v1_ = ((v1 - edge1_mid) * mat_rot) + edge1_mid
    v2_ = ((v2 - edge1_mid) * mat_rot) + edge1_mid
    v3_ = ((v3 - edge2_mid) * mat_rot) + edge2_mid
    v4_ = ((v4 - edge2_mid) * mat_rot) + edge2_mid

    r = geometry.intersect_line_line(v1_, v2_, v3_, v4_)
    if r:
        p1, _ = r
        # cp = transform_matrix * (p1 + origin)
        cp = transform_matrix * p1
        bpy.context.scene.cursor_location = cp
        # generate_gp3d_stroke(cp, axis, obj, radius=(p1-v1).length)
    else:
        print('not on a circle')

def createMirror(angle):
    mirror[0].xyz = [0,-1,-1]
    mirror[1].xyz = [-1,1,1]
    mirror[2].xyz = [1,1,1]
    mirror[0].rotate( Matrix.Rotation( angle, 4, norm_mirror) )
    mirror[1].rotate( Matrix.Rotation( angle, 4, norm_mirror)  )
    mirror[2].rotate( Matrix.Rotation( angle, 4, norm_mirror)  )
    n = geometry.normal(mirror[2],mirror[1],mirror[0])
    return [mirror[2],mirror[1],mirror[0],n]

def getPolyNormal(poly):
	""" Returns the normal of poligon based on the position of their vertex. It calculates the normal, it doesn't return manually modified normals.

	:param poly: The poligon.
	:type poly: |KX_PolyProxy|
	"""

	mesh = poly.getMesh()
	s = poly.getNumVertex()
	v1 = mesh.getVertex(0, poly.v1)
	v2 = mesh.getVertex(0, poly.v2)
	v3 = mesh.getVertex(0, poly.v3)
	if s == 4: v4v = mesh.getVertex(0, poly.v4).XYZ
	else: v4v = None

	if v4v: normal = geometry.normal(v1.XYZ, v2.XYZ, v3.XYZ, v4v)
	else: normal = geometry.normal(v1.XYZ, v2.XYZ, v3.XYZ)
	return normal

def flatFace(verts):
    #print("    Flatface: original vertices %s %s %s"%(verts[0],verts[1],verts[2]))
    quatr = mg.normal(verts).rotation_difference(Vector((0, 0, 1)))
    eul = quatr.to_euler()
    eul.z = 0.0
    #print("    Flatface called Euler is x=%2.2g deg y=%2.2g deg"%(degrees(eul.x),degrees(eul.y)))
    for v in verts:
        v.rotate(eul)
        v.z = 0

    def execute(self, context):
        obj = context.edit_object
        me = obj.data
        bm = bmesh.from_edit_mesh(me)
        bm.select_history.validate()
        vertCount = 0
        selectedVerts = []
        for vert in (vert for vert in bm.verts if vert.select):  #count loop verts
            vertCount+=1
            if len(selectedVerts)<3:
                selectedVerts.append(vert)

        if vertCount < 1:
            self.report({'INFO'}, 'Select at least one vertex')
            return {'CANCELLED'}


        vertsMedianLoc = mathutils.Vector((0.0 ,0.0 ,0.0))
        vertsMedianNorm = mathutils.Vector((0.0 ,0.0 ,0.0))
        vertsCoordList = []
        vectOffset = mathutils.Vector((0.0, 0.0, self.Depth))
        offsetMatrix = mathutils.Matrix.Translation(vectOffset)
        yawMatrix = mathutils.Matrix.Rotation(radians(self.TrimRoll), 4, 'X')
        pitchMatrix = mathutils.Matrix.Rotation(radians(self.TripPitch), 4, 'Y')

        for vert in selectedVerts:
            vertsMedianLoc+=vert.co
            vertsCoordList.append(vert.co)
            vertsMedianNorm+=vert.normal
        vertsMedianLoc /= len(selectedVerts)
        if vertCount>2:
            trisNorm=normal(vertsCoordList)
            if degrees(vertsMedianNorm.angle(trisNorm))>90:
                vertsMedianNorm= trisNorm*(-1)
            else:
                vertsMedianNorm = trisNorm
        if vertCount==2:  #fixes normal for 2 verts
            vectorV1_V2=selectedVerts[0].co-selectedVerts[1].co
            perpendicular= vectorV1_V2.cross(vertsMedianNorm)
            vertsMedianNorm = perpendicular.cross(vectorV1_V2)


        normToRotation = vertsMedianNorm.to_track_quat('Z', 'X').to_euler()
        selectionLocalMatrix = mathutils.Matrix.Translation(vertsMedianLoc)
        selectionLocalMatrix *= selectionLocalMatrix.Rotation(normToRotation[2], 4, 'Z')
        selectionLocalMatrix *= selectionLocalMatrix.Rotation(normToRotation[1], 4, 'Y')
        selectionLocalMatrix *= selectionLocalMatrix.Rotation(normToRotation[0], 4, 'X')

        spaceMatrixTransform = context.object.matrix_world*selectionLocalMatrix*offsetMatrix*yawMatrix*pitchMatrix

        loc = spaceMatrixTransform.to_translation()
        norm = spaceMatrixTransform.to_3x3()*mathutils.Vector((0.0, 0.0, 1.0))
        bpy.ops.mesh.select_all(action='SELECT')
        bpy.ops.mesh.bisect(plane_co=loc, plane_no=norm, use_fill=True, clear_inner= not self.clearOut, clear_outer=self.clearOut)
        bpy.ops.mesh.select_all(action='DESELECT')
        bmesh.update_edit_mesh(me, True)
        return {"FINISHED"}

def _ascii_write(filepath, faces):
    with open(filepath, 'w') as data:
        fw = data.write
        header = _header_version()
        fw('solid %s\n' % header)

        for face in faces:
            # calculate face normal
            fw('facet normal %f %f %f\nouter loop\n' % normal(*face)[:])
            for vert in face:
                fw('vertex %f %f %f\n' % vert[:])
            fw('endloop\nendfacet\n')

        fw('endsolid %s\n' % header)

def build_mirror_pivot(f,verts,ob):
  v0 = mathutils.Vector(verts[f[0]])
  v1 = mathutils.Vector(verts[f[1]])
  v2 = mathutils.Vector(verts[f[2]])
  if len(f) == 4:
    v3 = mathutils.Vector(verts[f[3]])
    norm = normal(v0,v1,v2,v3)
  else:
    v3 = None
    norm = normal(v0,v1,v2)
  rot = norm.to_track_quat('X','Z')
  r = rot.to_euler('XYZ')
  p0 = v0.lerp(v1,0.5)
  if v3 == None: p1 = v2
  else: p1 = v2.lerp(v3,0.5)
  center = p0.lerp(p1,0.5)

  dummy = bpy.data.objects.new(ob.name + "_mirror_pivot",None)
  dummy.location = center.to_tuple()
  print("ROTATING TO ",r.x," ",r.y," ",r.z)
  dummy.rotation_euler = (r.x,r.y,r.z)
  bpy.context.scene.objects.link(dummy)
  return dummy

def generate_3PT_mode_1(pts=None, num_verts=20, make_edges=False):
    '''
    Arc from start - throught - Eend
    - call this function only if you have 3 pts,
    - do your error checking before passing to it.
    '''
    num_verts -= 1
    verts, edges = [], []
    V = Vector

    # construction
    v1, v2, v3, v4 = V(pts[0]), V(pts[1]), V(pts[1]), V(pts[2])
    edge1_mid = v1.lerp(v2, 0.5)
    edge2_mid = v3.lerp(v4, 0.5)
    axis = geometry.normal(v1, v2, v4)
    mat_rot = mathutils.Matrix.Rotation(math.radians(90.0), 4, axis)

    # triangle edges
    v1_ = ((v1 - edge1_mid) * mat_rot) + edge1_mid
    v2_ = ((v2 - edge1_mid) * mat_rot) + edge1_mid
    v3_ = ((v3 - edge2_mid) * mat_rot) + edge2_mid
    v4_ = ((v4 - edge2_mid) * mat_rot) + edge2_mid

    r = geometry.intersect_line_line(v1_, v2_, v3_, v4_)
    if r:
        # do arc
        p1, _ = r

        # find arc angle.
        a = (v1 - p1).angle((v4 - p1), 0)
        s = (2 * math.pi) - a

        interior_angle = (v1 - v2).angle(v4 - v3, 0)
        if interior_angle > 0.5 * math.pi:
            s = math.pi + 2 * (0.5 * math.pi - interior_angle)

        for i in range(num_verts + 1):
            mat_rot = mathutils.Matrix.Rotation(((s / num_verts) * i), 4, axis)
            vec = ((v4 - p1) * mat_rot) + p1
            verts.append(vec[:])
    else:
        # do straight line
        step_size = 1 / num_verts
        verts = [v1_.lerp(v4_, i * step_size)[:] for i in range(num_verts + 1)]

    if make_edges:
        edges = [(n, n + 1) for n in range(len(verts) - 1)]

    return verts, edges

def compute_distances_mu(plane, pts, result, gates, tolerance):
    plane_origin = V(plane[0])
    plane_a, plane_b = V(plane[1]), V(plane[2])
    norm = normal([plane_origin, plane_a, plane_b])
    if norm.length == 0:
        print("Error: the three points of the plane are aligned. Not valid plane")
    local_result = [[] for res in result]
    for p in pts:
        data = compute_point_tri_dist(V(p), plane_origin, plane_a, plane_b, norm, tolerance)
        for i, r in enumerate(local_result):
            r.append(data[i])

    for i, res in enumerate(result):
        if gates[i]:
            res.append(local_result[i])

def generate_3PT_mode_1(pts, obj, nv):
    origin = obj.location
    mw = obj.matrix_world
    V = Vector

    nv = max(3, nv)

    # construction
    v1, v2, v3, v4 = V(pts[0]), V(pts[1]), V(pts[1]), V(pts[2])
    edge1_mid = v1.lerp(v2, 0.5)
    edge2_mid = v3.lerp(v4, 0.5)
    axis = geometry.normal(v1, v2, v4)
    mat_rot = mathutils.Matrix.Rotation(math.radians(90.0), 4, axis)

    # triangle edges
    v1_ = ((v1 - edge1_mid) * mat_rot) + edge1_mid
    v2_ = ((v2 - edge1_mid) * mat_rot) + edge1_mid
    v3_ = ((v3 - edge2_mid) * mat_rot) + edge2_mid
    v4_ = ((v4 - edge2_mid) * mat_rot) + edge2_mid

    r = geometry.intersect_line_line(v1_, v2_, v3_, v4_)
    if r:
        p1, _ = r
        cp = mw * p1
        bpy.context.scene.cursor_location = cp
        layer = get_layer()
        generate_gp3d_stroke(layer, p1, v1, axis, mw, origin, nv)

        ''' 
        # f = [i for i in dir(bpy.context) if 'gpencil' in i]
        active_gpencil_frame
        active_gpencil_layer
        editable_gpencil_layers
        editable_gpencil_strokes
        gpencil_data
        gpencil_data_owner
        visible_gpencil_layers
        '''
        
        #bpy.context.active_gpencil_layer = layer
        #print(bpy.context.gpencil_data)
        scn = bpy.context.scene
        scn.grease_pencil = bpy.data.grease_pencil['tc_circle_000']

    else:
        print('not on a circle')

def generate_3PT(pts, obj, nv, mode=0):
    origin = obj.location
    mw = obj.matrix_world
    V = Vector

    nv = max(3, nv)

    # construction
    v1, v2, v3, v4 = V(pts[0]), V(pts[1]), V(pts[1]), V(pts[2])
    edge1_mid = v1.lerp(v2, 0.5)
    edge2_mid = v3.lerp(v4, 0.5)
    axis = geometry.normal(v1, v2, v4)
    mat_rot = mathutils.Matrix.Rotation(math.radians(90.0), 4, axis)

    # triangle edges
    v1_ = ((v1 - edge1_mid) * mat_rot) + edge1_mid
    v2_ = ((v2 - edge1_mid) * mat_rot) + edge1_mid
    v3_ = ((v3 - edge2_mid) * mat_rot) + edge2_mid
    v4_ = ((v4 - edge2_mid) * mat_rot) + edge2_mid

    r = geometry.intersect_line_line(v1_, v2_, v3_, v4_)
    if r:
        p1, _ = r
        cp = mw * p1
        bpy.context.scene.cursor_location = cp

        if mode == 0:
            layer = get_layer()
            generate_gp3d_stroke(layer, p1, v1, axis, mw, origin, nv)

            scn = bpy.context.scene
            scn.grease_pencil = bpy.data.grease_pencil['tc_circle_000']

        elif mode == 1:
            generate_bmesh_repr(p1, v1, axis, nv)

    else:
        print('not on a circle')

    def process(self):
        if self.outputs['Centers'].is_linked or self.outputs['Normals'].is_linked or \
                self.outputs['Origins'].is_linked or self.outputs['Norm_abs'].is_linked:
            if 'Polygons' in self.inputs and 'Vertices' in self.inputs \
                and self.inputs['Polygons'].is_linked and self.inputs['Vertices'].is_linked:

                pols_ = SvGetSocketAnyType(self, self.inputs['Polygons'])
                vers_tupls = SvGetSocketAnyType(self, self.inputs['Vertices'])
                vers_vects = Vector_generate(vers_tupls)
                
                # make mesh temp утилитарно - удалить в конце
                mat_collect = []
                normals_out = []
                origins = []
                norm_abs_out = []
                for verst, versv, pols in zip(vers_tupls, vers_vects, pols_):
                    normals = []
                    centrs = []
                    norm_abs = []
                    p0_xdirs = []
                    for p in pols:
                        v0 = versv[p[0]]
                        v1 = versv[p[1]]
                        v2 = versv[p[2]]
                        # save direction of 1st point in polygon
                        p0_xdirs.append(v0)
                        # normals
                        norm = geometry.normal(v0, v1, v2)
                        normals.append(norm)                       
                        # centrs
                        x,y,z = zip(*[verst[poi] for poi in p])
                        x,y,z = sum(x)/len(x), sum(y)/len(y), sum(z)/len(z)
                        current_center = Vector((x,y,z))
                        centrs.append(current_center)
                        # normal absolute !!!
                        # это совершенно нормально!!! ;-)
                        norm_abs.append(current_center+norm)
                        
                        if self.Separate:
                            norm_abs_out.append(norm_abs)    
                            origins.append(centrs)
                            normals_out.append(normals)
                        else:
                            norm_abs_out.extend(norm_abs)    
                            origins.extend(centrs)
                            normals_out.extend(normals)
                    mat_collect_ = []

                    for cen, nor, p0 in zip(centrs, normals, p0_xdirs):
                        zdir = nor
                        xdir = (Vector(p0) - cen).normalized()
                        ydir = zdir.cross(xdir)
                        lM = [(xdir[0], ydir[0], zdir[0], cen[0]),
                              (xdir[1], ydir[1], zdir[1], cen[1]),
                              (xdir[2], ydir[2], zdir[2], cen[2]),
                              (0.0, 0.0, 0.0, 1.0)]
                        mat_collect_.append(lM)
                    mat_collect.extend(mat_collect_)

                if not self.Separate:
                    SvSetSocketAnyType(self, 'Centers', mat_collect)
                    SvSetSocketAnyType(self, 'Norm_abs', Vector_degenerate([norm_abs_out]))
                    SvSetSocketAnyType(self, 'Origins', Vector_degenerate([origins]))
                    SvSetSocketAnyType(self, 'Normals', Vector_degenerate([normals_out]))
                else:
                    SvSetSocketAnyType(self, 'Centers', mat_collect)
                    SvSetSocketAnyType(self, 'Norm_abs', Vector_degenerate(norm_abs_out))
                    SvSetSocketAnyType(self, 'Origins', Vector_degenerate(origins))
                    SvSetSocketAnyType(self, 'Normals', Vector_degenerate(normals_out))

 def calc_normal(self):
     """:rtype: Vector"""
     return geom.normal(*self.coords)

    def process(self):
        verts_socket, poly_socket = self.inputs
        norm_socket, norm_abs_socket, origins_socket, centers_socket = self.outputs

        if not any([s.is_linked for s in self.outputs]):
            return

        if not (verts_socket.is_linked and poly_socket.is_linked):
            return

        pols_ = poly_socket.sv_get()
        vers_tupls = verts_socket.sv_get()
        vers_vects = Vector_generate(vers_tupls)

        # make mesh temp утилитарно - удалить в конце
        mat_collect = []
        normals_out = []
        origins = []
        norm_abs_out = []
        for verst, versv, pols in zip(vers_tupls, vers_vects, pols_):
            normals = []
            centrs = []
            norm_abs = []
            p0_xdirs = []
            for p in pols:
                v0 = versv[p[0]]
                v1 = versv[p[1]]
                v2 = versv[p[2]]
                # save direction of 1st point in polygon
                p0_xdirs.append(v0)
                # normals
                norm = geometry.normal(v0, v1, v2)
                normals.append(norm)
                # centrs
                x,y,z = zip(*[verst[poi] for poi in p])
                x,y,z = sum(x)/len(x), sum(y)/len(y), sum(z)/len(z)
                current_center = Vector((x,y,z))
                centrs.append(current_center)
                # normal absolute !!!
                # это совершенно нормально!!! ;-)
                norm_abs.append(current_center+norm)

            if self.Separate:
                norm_abs_out.append(norm_abs)
                origins.append(centrs)
                normals_out.append(normals)
            else:
                norm_abs_out.extend(norm_abs)
                origins.extend(centrs)
                normals_out.extend(normals)
                
            mat_collect_ = []

            for cen, nor, p0 in zip(centrs, normals, p0_xdirs):
                zdir = nor
                xdir = (Vector(p0) - cen).normalized()
                ydir = zdir.cross(xdir)
                lM = [(xdir[0], ydir[0], zdir[0], cen[0]),
                      (xdir[1], ydir[1], zdir[1], cen[1]),
                      (xdir[2], ydir[2], zdir[2], cen[2]),
                      (0.0, 0.0, 0.0, 1.0)]
                mat_collect_.append(Matrix(lM))
            mat_collect.extend(mat_collect_)

        if not self.Separate:
            norm_abs_out = [norm_abs_out]
            origins = [origins]
            normals_out = [normals_out]

        centers_socket.sv_set(mat_collect)
        norm_abs_socket.sv_set(Vector_degenerate(norm_abs_out))
        origins_socket.sv_set(Vector_degenerate(origins))
        norm_socket.sv_set(Vector_degenerate(normals_out))