mesh_tissue/weight_tools.py

   1 # SPDX-License-Identifier: GPL-2.0-or-later
   2
   3 #-------------------------- COLORS / GROUPS EXCHANGER -------------------------#
   4 #                                                                              #
   5 # Vertex Color to Vertex Group allow you to convert colors channels to weight  #
   6 # maps.                                                                        #
   7 # The main purpose is to use vertex colors to store information when importing #
   8 # files from other software. The script works with the active vertex color     #
   9 # slot.                                                                        #
  10 # For use the command "Vertex Clors to Vertex Groups" use the search bar       #
  11 # (space bar).                                                                 #
  12 #                                                                              #
  13 #                          (c)  Alessandro Zomparelli                          #
  14 #                                     (2017)                                   #
  15 #                                                                              #
  16 # http://www.co-de-it.com/                                                     #
  17 #                                                                              #
  18 ################################################################################
  19
  20 import bpy, bmesh, os
  21 import numpy as np
  22 import math, timeit, time
  23 from math import pi
  24 from statistics import mean, stdev
  25 from mathutils import Vector
  26 from mathutils.kdtree import KDTree
  27 from numpy import *
  28 try: from .numba_functions import numba_reaction_diffusion, numba_reaction_diffusion_anisotropic, integrate_field
  29 except: pass
  30 #from .numba_functions import integrate_field
  31 #from .numba_functions import numba_reaction_diffusion
  32 try: import numexpr as ne
  33 except: pass
  34
  35 # Reaction-Diffusion cache
  36 from pathlib import Path
  37 import random as rnd
  38 import string
  39
  40 from bpy.types import (
  41         Operator,
  42         Panel,
  43         PropertyGroup,
  44         )
  45
  46 from bpy.props import (
  47     BoolProperty,
  48     EnumProperty,
  49     FloatProperty,
  50     IntProperty,
  51     StringProperty,
  52     FloatVectorProperty,
  53     IntVectorProperty
  54 )
  55
  56 from .utils import *
  57
  58 def reaction_diffusion_add_handler(self, context):
  59     # remove existing handlers
  60     reaction_diffusion_remove_handler(self, context)
  61     # add new handler
  62     bpy.app.handlers.frame_change_post.append(reaction_diffusion_scene)
  63
  64 def reaction_diffusion_remove_handler(self, context):
  65     # remove existing handlers
  66     old_handlers = []
  67     for h in bpy.app.handlers.frame_change_post:
  68         if "reaction_diffusion" in str(h):
  69             old_handlers.append(h)
  70     for h in old_handlers: bpy.app.handlers.frame_change_post.remove(h)
  71
  72 class formula_prop(PropertyGroup):
  73     name : StringProperty()
  74     formula : StringProperty()
  75     float_var : FloatVectorProperty(name="", description="", default=(0, 0, 0, 0, 0), size=5)
  76     int_var : IntVectorProperty(name="", description="", default=(0, 0, 0, 0, 0), size=5)
  77
  78 class reaction_diffusion_prop(PropertyGroup):
  79     run : BoolProperty(default=False, update = reaction_diffusion_add_handler,
  80         description='Compute a new iteration on frame changes. Currently is not working during  Render Animation')
  81
  82     time_steps : IntProperty(
  83         name="Steps", default=10, min=0, soft_max=50,
  84         description="Number of Steps")
  85
  86     dt : FloatProperty(
  87         name="dt", default=1, min=0, soft_max=0.2,
  88         description="Time Step")
  89
  90     diff_a : FloatProperty(
  91         name="Diff A", default=0.1, min=0, soft_max=2, precision=3,
  92         description="Diffusion A")
  93
  94     diff_b : FloatProperty(
  95         name="Diff B", default=0.05, min=0, soft_max=2, precision=3,
  96         description="Diffusion B")
  97
  98     f : FloatProperty(
  99         name="f", default=0.055, soft_min=0.01, soft_max=0.06, precision=4, step=0.05,
 100         description="Feed Rate")
 101
 102     k : FloatProperty(
 103         name="k", default=0.062, soft_min=0.035, soft_max=0.065, precision=4, step=0.05,
 104         description="Kill Rate")
 105
 106     diff_mult : FloatProperty(
 107         name="Scale", default=1, min=0, soft_max=1, max=10, precision=2,
 108         description="Multiplier for the diffusion of both substances")
 109
 110     vertex_group_diff_a : StringProperty(
 111         name="Diff A", default='',
 112         description="Vertex Group used for A diffusion")
 113
 114     vertex_group_diff_b : StringProperty(
 115         name="Diff B", default='',
 116         description="Vertex Group used for B diffusion")
 117
 118     vertex_group_scale : StringProperty(
 119         name="Scale", default='',
 120         description="Vertex Group used for Scale value")
 121
 122     vertex_group_f : StringProperty(
 123         name="f", default='',
 124         description="Vertex Group used for Feed value (f)")
 125
 126     vertex_group_k : StringProperty(
 127         name="k", default='',
 128         description="Vertex Group used for Kill value (k)")
 129
 130     vertex_group_brush : StringProperty(
 131         name="Brush", default='',
 132         description="Vertex Group used for adding/removing B")
 133
 134     invert_vertex_group_diff_a : BoolProperty(default=False,
 135         description='Inverte the value of the Vertex Group Diff A')
 136
 137     invert_vertex_group_diff_b : BoolProperty(default=False,
 138         description='Inverte the value of the Vertex Group Diff B')
 139
 140     invert_vertex_group_scale : BoolProperty(default=False,
 141         description='Inverte the value of the Vertex Group Scale')
 142
 143     invert_vertex_group_f : BoolProperty(default=False,
 144         description='Inverte the value of the Vertex Group f')
 145
 146     invert_vertex_group_k : BoolProperty(default=False,
 147         description='Inverte the value of the Vertex Group k')
 148
 149     min_diff_a : FloatProperty(
 150         name="Min Diff A", default=0.1, min=0, soft_max=2, precision=3,
 151         description="Min Diff A")
 152
 153     max_diff_a : FloatProperty(
 154         name="Max Diff A", default=0.1, min=0, soft_max=2, precision=3,
 155         description="Max Diff A")
 156
 157     min_diff_b : FloatProperty(
 158         name="Min Diff B", default=0.1, min=0, soft_max=2, precision=3,
 159         description="Min Diff B")
 160
 161     max_diff_b : FloatProperty(
 162         name="Max Diff B", default=0.1, min=0, soft_max=2, precision=3,
 163         description="Max Diff B")
 164
 165     min_scale : FloatProperty(
 166         name="Scale", default=0.35, min=0, soft_max=1, max=10, precision=2,
 167         description="Min Scale Value")
 168
 169     max_scale : FloatProperty(
 170         name="Scale", default=1, min=0, soft_max=1, max=10, precision=2,
 171         description="Max Scale value")
 172
 173     min_f : FloatProperty(
 174         name="Min f", default=0.02, min=0, soft_min=0.01, soft_max=0.06, max=0.1, precision=4, step=0.05,
 175         description="Min Feed Rate")
 176
 177     max_f : FloatProperty(
 178         name="Max f", default=0.055, min=0, soft_min=0.01, soft_max=0.06, max=0.1, precision=4, step=0.05,
 179         description="Max Feed Rate")
 180
 181     min_k : FloatProperty(
 182         name="Min k", default=0.035, min=0, soft_min=0.035, soft_max=0.065, max=0.1, precision=4, step=0.05,
 183         description="Min Kill Rate")
 184
 185     max_k : FloatProperty(
 186         name="Max k", default=0.062, min=0, soft_min=0.035, soft_max=0.065, max=0.1, precision=4, step=0.05,
 187         description="Max Kill Rate")
 188
 189     brush_mult : FloatProperty(
 190         name="Mult", default=0.5, min=-1, max=1, precision=3, step=0.05,
 191         description="Multiplier for brush value")
 192
 193     bool_mod : BoolProperty(
 194         name="Use Modifiers", default=False,
 195         description="Read modifiers affect the vertex groups")
 196
 197     bool_cache : BoolProperty(
 198         name="Use Cache", default=False,
 199         description="Read modifiers affect the vertex groups")
 200
 201     cache_frame_start : IntProperty(
 202         name="Start", default=1,
 203         description="Frame on which the simulation starts")
 204
 205     cache_frame_end : IntProperty(
 206         name="End", default=250,
 207         description="Frame on which the simulation ends")
 208
 209     cache_dir : StringProperty(
 210         name="Cache directory", default="", subtype='FILE_PATH',
 211         description = 'Directory that contains Reaction-Diffusion cache files'
 212         )
 213
 214     update_weight_a : BoolProperty(
 215         name="Update Vertex Group A", default=True,
 216         description="Transfer Cache to the Vertex Groups named A")
 217
 218     update_weight_b : BoolProperty(
 219         name="Update Vertex Group B", default=True,
 220         description="Transfer Cache to the Vertex Groups named B")
 221
 222     update_colors_a : BoolProperty(
 223         name="Update Vertex Color A", default=False,
 224         description="Transfer Cache to the Vertex Color named A")
 225
 226     update_colors_b : BoolProperty(
 227         name="Update Vertex Color B", default=False,
 228         description="Transfer Cache to the Vertex Color named B")
 229
 230     update_colors : BoolProperty(
 231         name="Update Vertex Color AB", default=False,
 232         description="Transfer Cache to the Vertex Color named AB")
 233
 234     update_uv : BoolProperty(
 235         name="Update UV", default=False,
 236         description="Transfer Cache to the UV Map Layer named AB")
 237
 238     normalize : BoolProperty(
 239         name="Normalize values", default=False,
 240         description="Normalize values from 0 to 1")
 241
 242     fast_bake : BoolProperty(
 243         name="Fast Bake", default=True,
 244         description="Do not update modifiers or vertex groups while baking. Much faster!")
 245
 246
 247 from numpy import *
 248 def compute_formula(ob=None, formula="rx", float_var=(0,0,0,0,0), int_var=(0,0,0,0,0)):
 249     verts = ob.data.vertices
 250     n_verts = len(verts)
 251
 252     f1,f2,f3,f4,f5 = float_var
 253     i1,i2,i3,i4,i5 = int_var
 254
 255     do_groups = "w[" in formula
 256     do_local = "lx" in formula or "ly" in formula or "lz" in formula
 257     do_global = "gx" in formula or "gy" in formula or "gz" in formula
 258     do_relative = "rx" in formula or "ry" in formula or "rz" in formula
 259     do_normal = "nx" in formula or "ny" in formula or "nz" in formula
 260     mat = ob.matrix_world
 261
 262     for i in range(1000):
 263         if "w["+str(i)+"]" in formula and i > len(ob.vertex_groups)-1:
 264             return "w["+str(i)+"] not found"
 265
 266     w = []
 267     for i in range(len(ob.vertex_groups)):
 268         w.append([])
 269         if "w["+str(i)+"]" in formula:
 270             vg = ob.vertex_groups[i]
 271             for v in verts:
 272                 try:
 273                     w[i].append(vg.weight(v.index))
 274                 except:
 275                     w[i].append(0)
 276             w[i] = array(w[i])
 277
 278     start_time = timeit.default_timer()
 279     # compute vertex coordinates
 280     if do_local or do_relative or do_global:
 281         co = [0]*n_verts*3
 282         verts.foreach_get('co', co)
 283         np_co = array(co).reshape((n_verts, 3))
 284         lx, ly, lz = array(np_co).transpose()
 285         if do_relative:
 286             rx = np.interp(lx, (lx.min(), lx.max()), (0, +1))
 287             ry = np.interp(ly, (ly.min(), ly.max()), (0, +1))
 288             rz = np.interp(lz, (lz.min(), lz.max()), (0, +1))
 289         if do_global:
 290             co = [v.co for v in verts]
 291             global_co = []
 292             for v in co:
 293                 global_co.append(mat @ v)
 294             global_co = array(global_co).reshape((n_verts, 3))
 295             gx, gy, gz = array(global_co).transpose()
 296     # compute vertex normals
 297     if do_normal:
 298         normal = [0]*n_verts*3
 299         verts.foreach_get('normal', normal)
 300         normal = array(normal).reshape((n_verts, 3))
 301         nx, ny, nz = array(normal).transpose()
 302
 303     try:
 304         weight = eval(formula)
 305         return weight
 306     except:
 307         return "There is something wrong"
 308     print("Weight Formula: " + str(timeit.default_timer() - start_time))
 309
 310 class weight_formula_wiki(Operator):
 311     bl_idname = "scene.weight_formula_wiki"
 312     bl_label = "Online Documentation"
 313     bl_options = {'REGISTER', 'UNDO'}
 314
 315     def execute(self, context):
 316         bpy.ops.wm.url_open(url="https://github.com/alessandro-zomparelli/tissue/wiki/Weight-Tools#weight-formula")
 317         return {'FINISHED'}
 318
 319 class weight_formula(Operator):
 320     bl_idname = "object.weight_formula"
 321     bl_label = "Weight Formula"
 322     bl_description = "Generate a Vertex Group according to a mathematical formula"
 323     bl_options = {'REGISTER', 'UNDO'}
 324
 325     ex_items = [
 326         ('cos(arctan(nx/ny)*i1*2 + sin(rz*i3))/i2 + cos(arctan(nx/ny)*i1*2 - sin(rz*i3))/i2 + 0.5','Vertical Spots'),
 327         ('cos(arctan(nx/ny)*i1*2 + sin(rz*i2))/2 + cos(arctan(nx/ny)*i1*2 - sin(rz*i2))/2','Vertical Spots'),
 328         ('(sin(arctan(nx/ny)*i1*2)*sin(nz*i1*2)+1)/2','Grid Spots'),
 329         ('cos(arctan(nx/ny)*f1)','Vertical Stripes'),
 330         ('cos(arctan(lx/ly)*f1 + sin(rz*f2)*f3)','Curly Stripes'),
 331         ('sin(rz*pi*i1+arctan2(nx,ny))/2+0.5', 'Vertical Spiral'),
 332         ('sin(nx*15)<sin(ny*15)','Chess'),
 333         ('cos(ny*rz**2*i1)','Hyperbolic'),
 334         ('sin(rx*30) > 0','Step Stripes'),
 335         ('sin(nz*i1)','Normal Stripes'),
 336         ('w[0]**2','Vertex Group square'),
 337         ('abs(0.5-rz)*2','Double vertical gradient'),
 338         ('rz', 'Vertical Gradient')
 339     ]
 340     _ex_items = list((str(i),'{}   ( {} )'.format(s[0],s[1]),s[1]) for i,s in enumerate(ex_items))
 341     _ex_items.append(('CUSTOM', "User Formula", ""))
 342
 343     examples : EnumProperty(
 344         items = _ex_items, default='CUSTOM', name="Examples")
 345
 346     old_ex = ""
 347
 348     formula : StringProperty(
 349         name="Formula", default="", description="Formula to Evaluate")
 350
 351     slider_f01 : FloatProperty(
 352         name="f1", default=1, description="Slider Float 1")
 353     slider_f02 : FloatProperty(
 354         name="f2", default=1, description="Slider Float 2")
 355     slider_f03 : FloatProperty(
 356         name="f3", default=1, description="Slider Float 3")
 357     slider_f04 : FloatProperty(
 358         name="f4", default=1, description="Slider Float 4")
 359     slider_f05 : FloatProperty(
 360         name="f5", default=1, description="Slider Float 5")
 361     slider_i01 : IntProperty(
 362         name="i1", default=1, description="Slider Integer 1")
 363     slider_i02 : IntProperty(
 364         name="i2", default=1, description="Slider Integer 2")
 365     slider_i03 : IntProperty(
 366         name="i3", default=1, description="Slider Integer 3")
 367     slider_i04 : IntProperty(
 368         name="i4", default=1, description="Slider Integer 4")
 369     slider_i05 : IntProperty(
 370         name="i5", default=1, description="Slider Integer 5")
 371
 372     def invoke(self, context, event):
 373         return context.window_manager.invoke_props_dialog(self, width=350)
 374
 375     def draw(self, context):
 376         layout = self.layout
 377         #layout.label(text="Examples")
 378         layout.prop(self, "examples", text="Examples")
 379         #if self.examples == 'CUSTOM':
 380         layout.label(text="Formula")
 381         layout.prop(self, "formula", text="")
 382         #try: self.examples = self.formula
 383         #except: pass
 384
 385         if self.examples != 'CUSTOM':
 386             example = self.ex_items[int(self.examples)][0]
 387             if example != self.old_ex:
 388                 self.formula = example
 389                 self.old_ex = example
 390             elif self.formula != example:
 391                 self.examples = 'CUSTOM'
 392         formula = self.formula
 393
 394         layout.separator()
 395         if "f1" in formula: layout.prop(self, "slider_f01")
 396         if "f2" in formula: layout.prop(self, "slider_f02")
 397         if "f3" in formula: layout.prop(self, "slider_f03")
 398         if "f4" in formula: layout.prop(self, "slider_f04")
 399         if "f5" in formula: layout.prop(self, "slider_f05")
 400         if "i1" in formula: layout.prop(self, "slider_i01")
 401         if "i2" in formula: layout.prop(self, "slider_i02")
 402         if "i3" in formula: layout.prop(self, "slider_i03")
 403         if "i4" in formula: layout.prop(self, "slider_i04")
 404         if "i5" in formula: layout.prop(self, "slider_i05")
 405
 406         layout.label(text="Variables (for each vertex):")
 407         layout.label(text="lx, ly, lz: Local Coordinates", icon='ORIENTATION_LOCAL')
 408         layout.label(text="gx, gy, gz: Global Coordinates", icon='WORLD')
 409         layout.label(text="rx, ry, rz: Local Coordinates (0 to 1)", icon='NORMALIZE_FCURVES')
 410         layout.label(text="nx, ny, nz: Normal Coordinates", icon='SNAP_NORMAL')
 411         layout.label(text="w[0], w[1], w[2], ... : Vertex Groups", icon="GROUP_VERTEX")
 412         layout.separator()
 413         layout.label(text="f1, f2, f3, f4, f5: Float Sliders", icon='MOD_HUE_SATURATION')#PROPERTIES
 414         layout.label(text="i1, i2, i3, i4, i5: Integer Sliders", icon='MOD_HUE_SATURATION')
 415         layout.separator()
 416         #layout.label(text="All mathematical functions are based on Numpy", icon='INFO')
 417         #layout.label(text="https://docs.scipy.org/doc/numpy-1.13.0/reference/routines.math.html", icon='INFO')
 418         layout.operator("scene.weight_formula_wiki", icon="HELP")
 419         #layout.label(text="(where 'i' is the index of the Vertex Group)")
 420
 421     def execute(self, context):
 422         ob = context.active_object
 423         n_verts = len(ob.data.vertices)
 424         #if self.examples == 'CUSTOM':
 425         #    formula = self.formula
 426         #else:
 427         #self.formula = self.examples
 428         #    formula = self.examples
 429
 430         #f1, f2, f3, f4, f5 = self.slider_f01, self.slider_f02, self.slider_f03, self.slider_f04, self.slider_f05
 431         #i1, i2, i3, i4, i5 = self.slider_i01, self.slider_i02, self.slider_i03, self.slider_i04, self.slider_i05
 432         f_sliders = self.slider_f01, self.slider_f02, self.slider_f03, self.slider_f04, self.slider_f05
 433         i_sliders = self.slider_i01, self.slider_i02, self.slider_i03, self.slider_i04, self.slider_i05
 434
 435         if self.examples != 'CUSTOM':
 436             example = self.ex_items[int(self.examples)][0]
 437             if example != self.old_ex:
 438                 self.formula = example
 439                 self.old_ex = example
 440             elif self.formula != example:
 441                 self.examples = 'CUSTOM'
 442         formula = self.formula
 443
 444         if formula == "": return {'FINISHED'}
 445         # replace numeric sliders value
 446         for i, slider in enumerate(f_sliders):
 447             formula = formula.replace('f'+str(i+1),"{0:.2f}".format(slider))
 448         for i, slider in enumerate(i_sliders):
 449             formula =formula.replace('i'+str(i+1),str(slider))
 450         vertex_group_name = "" + formula
 451         ob.vertex_groups.new(name=vertex_group_name)
 452
 453         weight = compute_formula(ob, formula=formula, float_var=f_sliders, int_var=i_sliders)
 454         if type(weight) == str:
 455             self.report({'ERROR'}, weight)
 456             return {'CANCELLED'}
 457
 458         #start_time = timeit.default_timer()
 459         weight = nan_to_num(weight)
 460         vg = ob.vertex_groups[-1]
 461         if type(weight) == int or type(weight) == float:
 462             for i in range(n_verts):
 463                 vg.add([i], weight, 'REPLACE')
 464         elif type(weight) == ndarray:
 465             for i in range(n_verts):
 466                 vg.add([i], weight[i], 'REPLACE')
 467         ob.data.update()
 468         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
 469
 470         # Store formula settings
 471         new_formula = ob.formula_settings.add()
 472         new_formula.name = ob.vertex_groups[-1].name
 473         new_formula.formula = formula
 474         new_formula.int_var = i_sliders
 475         new_formula.float_var = f_sliders
 476
 477         #for f in ob.formula_settings:
 478         #    print(f.name, f.formula, f.int_var, f.float_var)
 479         return {'FINISHED'}
 480
 481
 482 class update_weight_formula(Operator):
 483     bl_idname = "object.update_weight_formula"
 484     bl_label = "Update Weight Formula"
 485     bl_description = "Update an existing Vertex Group. Make sure that the name\nof the active Vertex Group is a valid formula"
 486     bl_options = {'REGISTER', 'UNDO'}
 487
 488     @classmethod
 489     def poll(cls, context):
 490         return len(context.object.vertex_groups) > 0
 491
 492     def execute(self, context):
 493         ob = context.active_object
 494         n_verts = len(ob.data.vertices)
 495
 496         vg = ob.vertex_groups.active
 497         formula = vg.name
 498         weight = compute_formula(ob, formula=formula)
 499         if type(weight) == str:
 500             self.report({'ERROR'}, "The name of the active Vertex Group\nis not a valid Formula")
 501             return {'CANCELLED'}
 502
 503         #start_time = timeit.default_timer()
 504         weight = nan_to_num(weight)
 505         if type(weight) == int or type(weight) == float:
 506             for i in range(n_verts):
 507                 vg.add([i], weight, 'REPLACE')
 508         elif type(weight) == ndarray:
 509             for i in range(n_verts):
 510                 vg.add([i], weight[i], 'REPLACE')
 511         ob.data.update()
 512         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
 513         return {'FINISHED'}
 514
 515
 516 class _weight_laplacian(Operator):
 517     bl_idname = "object._weight_laplacian"
 518     bl_label = "Weight Laplacian"
 519     bl_description = ("Compute the Vertex Group Laplacian")
 520     bl_options = {'REGISTER', 'UNDO'}
 521
 522     bounds : EnumProperty(
 523         items=(('MANUAL', "Manual Bounds", ""),
 524             ('POSITIVE', "Positive Only", ""),
 525             ('NEGATIVE', "Negative Only", ""),
 526             ('AUTOMATIC', "Automatic Bounds", "")),
 527         default='AUTOMATIC', name="Bounds")
 528
 529     mode : EnumProperty(
 530         items=(('LENGTH', "Length Weight", ""),
 531             ('SIMPLE', "Simple", "")),
 532         default='SIMPLE', name="Evaluation Mode")
 533
 534     min_def : FloatProperty(
 535         name="Min", default=0, soft_min=-1, soft_max=0,
 536         description="Laplacian value with 0 weight")
 537
 538     max_def : FloatProperty(
 539         name="Max", default=0.5, soft_min=0, soft_max=5,
 540         description="Laplacian value with 1 weight")
 541
 542     bounds_string = ""
 543
 544     frame = None
 545
 546     @classmethod
 547     def poll(cls, context):
 548         return len(context.object.vertex_groups) > 0
 549
 550     def draw(self, context):
 551         layout = self.layout
 552         col = layout.column(align=True)
 553         col.label(text="Evaluation Mode")
 554         col.prop(self, "mode", text="")
 555         col.label(text="Bounds")
 556         col.prop(self, "bounds", text="")
 557         if self.bounds == 'MANUAL':
 558             col.label(text="Strain Rate \u03B5:")
 559             col.prop(self, "min_def")
 560             col.prop(self, "max_def")
 561         col.label(text="\u03B5" + ": from " + self.bounds_string)
 562
 563
 564     def execute(self, context):
 565         try: ob = context.object
 566         except:
 567             self.report({'ERROR'}, "Please select an Object")
 568             return {'CANCELLED'}
 569
 570         group_id = ob.vertex_groups.active_index
 571         input_group = ob.vertex_groups[group_id].name
 572
 573         group_name = "Laplacian"
 574         ob.vertex_groups.new(name=group_name)
 575         me = ob.data
 576         bm = bmesh.new()
 577         bm.from_mesh(me)
 578         bm.edges.ensure_lookup_table()
 579
 580         # store weight values
 581         weight = []
 582         for v in me.vertices:
 583             try:
 584                 weight.append(ob.vertex_groups[input_group].weight(v.index))
 585             except:
 586                 weight.append(0)
 587
 588         n_verts = len(bm.verts)
 589         lap = [0]*n_verts
 590         for e in bm.edges:
 591             if self.mode == 'LENGTH':
 592                 length = e.calc_length()
 593                 if length == 0: continue
 594                 id0 = e.verts[0].index
 595                 id1 = e.verts[1].index
 596                 lap[id0] += weight[id1]/length - weight[id0]/length
 597                 lap[id1] += weight[id0]/length - weight[id1]/length
 598             else:
 599                 id0 = e.verts[0].index
 600                 id1 = e.verts[1].index
 601                 lap[id0] += weight[id1] - weight[id0]
 602                 lap[id1] += weight[id0] - weight[id1]
 603
 604         mean_lap = mean(lap)
 605         stdev_lap = stdev(lap)
 606         filter_lap = [i for i in lap if mean_lap-2*stdev_lap < i < mean_lap+2*stdev_lap]
 607         if self.bounds == 'MANUAL':
 608             min_def = self.min_def
 609             max_def = self.max_def
 610         elif self.bounds == 'AUTOMATIC':
 611             min_def = min(filter_lap)
 612             max_def = max(filter_lap)
 613             self.min_def = min_def
 614             self.max_def = max_def
 615         elif self.bounds == 'NEGATIVE':
 616             min_def = 0
 617             max_def = min(filter_lap)
 618             self.min_def = min_def
 619             self.max_def = max_def
 620         elif self.bounds == 'POSITIVE':
 621             min_def = 0
 622             max_def = max(filter_lap)
 623             self.min_def = min_def
 624             self.max_def = max_def
 625         delta_def = max_def - min_def
 626
 627         # check undeformed errors
 628         if delta_def == 0: delta_def = 0.0001
 629
 630         for i in range(len(lap)):
 631             val = (lap[i]-min_def)/delta_def
 632             if val > 0.7: print(str(val) + " " + str(lap[i]))
 633             #val = weight[i] + 0.2*lap[i]
 634             ob.vertex_groups[-1].add([i], val, 'REPLACE')
 635         self.bounds_string = str(round(min_def,2)) + " to " + str(round(max_def,2))
 636         ob.vertex_groups[-1].name = group_name + " " + self.bounds_string
 637         ob.vertex_groups.update()
 638         ob.data.update()
 639         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
 640         bm.free()
 641         return {'FINISHED'}
 642
 643 class ok_weight_laplacian(Operator):
 644     bl_idname = "object.weight_laplacian"
 645     bl_label = "Weight Laplacian"
 646     bl_description = ("Compute the Vertex Group Laplacian")
 647     bl_options = {'REGISTER', 'UNDO'}
 648
 649     bounds_string = ""
 650
 651     frame = None
 652
 653     @classmethod
 654     def poll(cls, context):
 655         return len(context.object.vertex_groups) > 0
 656
 657
 658     def execute(self, context):
 659         try: ob = context.object
 660         except:
 661             self.report({'ERROR'}, "Please select an Object")
 662             return {'CANCELLED'}
 663
 664         me = ob.data
 665         bm = bmesh.new()
 666         bm.from_mesh(me)
 667         bm.edges.ensure_lookup_table()
 668
 669         group_id = ob.vertex_groups.active_index
 670         input_group = ob.vertex_groups[group_id].name
 671
 672         group_name = "Laplacian"
 673         ob.vertex_groups.new(name=group_name)
 674
 675         # store weight values
 676         a = []
 677         for v in me.vertices:
 678             try:
 679                 a.append(ob.vertex_groups[input_group].weight(v.index))
 680             except:
 681                 a.append(0)
 682
 683         a = array(a)
 684
 685
 686         # initialize
 687         n_verts = len(bm.verts)
 688         # find max number of edges for vertex
 689         max_edges = 0
 690         n_neighbors = []
 691         id_neighbors = []
 692         for v in bm.verts:
 693             n_edges = len(v.link_edges)
 694             max_edges = max(max_edges, n_edges)
 695             n_neighbors.append(n_edges)
 696             neighbors = []
 697             for e in v.link_edges:
 698                 for v1 in e.verts:
 699                     if v != v1: neighbors.append(v1.index)
 700             id_neighbors.append(neighbors)
 701         n_neighbors = array(n_neighbors)
 702
 703
 704         lap_map = [[] for i in range(n_verts)]
 705         #lap_map = []
 706         '''
 707         for e in bm.edges:
 708             id0 = e.verts[0].index
 709             id1 = e.verts[1].index
 710             lap_map[id0].append(id1)
 711             lap_map[id1].append(id0)
 712         '''
 713         lap = zeros((n_verts))#[0]*n_verts
 714         n_records = zeros((n_verts))
 715         for e in bm.edges:
 716             id0 = e.verts[0].index
 717             id1 = e.verts[1].index
 718             length = e.calc_length()
 719             if length == 0: continue
 720             #lap[id0] += abs(a[id1] - a[id0])/length
 721             #lap[id1] += abs(a[id0] - a[id1])/length
 722             lap[id0] += (a[id1] - a[id0])/length
 723             lap[id1] += (a[id0] - a[id1])/length
 724             n_records[id0]+=1
 725             n_records[id1]+=1
 726         lap /= n_records
 727         lap /= max(lap)
 728
 729         for i in range(n_verts):
 730             ob.vertex_groups['Laplacian'].add([i], lap[i], 'REPLACE')
 731         ob.vertex_groups.update()
 732         ob.data.update()
 733         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
 734         bm.free()
 735         return {'FINISHED'}
 736
 737 class weight_laplacian(Operator):
 738     bl_idname = "object.weight_laplacian"
 739     bl_label = "Weight Laplacian"
 740     bl_description = ("Compute the Vertex Group Laplacian")
 741     bl_options = {'REGISTER', 'UNDO'}
 742
 743     bounds_string = ""
 744
 745     frame = None
 746
 747     @classmethod
 748     def poll(cls, context):
 749         return len(context.object.vertex_groups) > 0
 750
 751
 752     def execute(self, context):
 753         try: ob = context.object
 754         except:
 755             self.report({'ERROR'}, "Please select an Object")
 756             return {'CANCELLED'}
 757
 758         me = ob.data
 759         bm = bmesh.new()
 760         bm.from_mesh(me)
 761         bm.edges.ensure_lookup_table()
 762         n_verts = len(me.vertices)
 763
 764         group_id = ob.vertex_groups.active_index
 765         input_group = ob.vertex_groups[group_id].name
 766
 767         group_name = "Laplacian"
 768         vg = ob.vertex_groups.new(name=group_name)
 769
 770         # store weight values
 771         dvert_lay = bm.verts.layers.deform.active
 772         weight = bmesh_get_weight_numpy(group_id, dvert_lay, bm.verts)
 773
 774         #verts, normals = get_vertices_and_normals_numpy(me)
 775
 776         #lap = zeros((n_verts))#[0]*n_verts
 777         lap = [Vector((0,0,0)) for i in range(n_verts)]
 778         n_records = zeros((n_verts))
 779         for e in bm.edges:
 780             vert0 = e.verts[0]
 781             vert1 = e.verts[1]
 782             id0 = vert0.index
 783             id1 = vert1.index
 784             v0 = vert0.co
 785             v1 = vert1.co
 786             v01 = v1-v0
 787             v10 = -v01
 788             v01 -= v01.project(vert0.normal)
 789             v10 -= v10.project(vert1.normal)
 790             length = e.calc_length()
 791             if length == 0: continue
 792             dw = (weight[id1] - weight[id0])/length
 793             lap[id0] += v01.normalized() * dw
 794             lap[id1] -= v10.normalized() * dw
 795             n_records[id0]+=1
 796             n_records[id1]+=1
 797         #lap /= n_records[:,np.newaxis]
 798         lap = [l.length/r for r,l in zip(n_records,lap)]
 799
 800         lap = np.array(lap)
 801         lap /= np.max(lap)
 802         lap = list(lap)
 803         print(lap)
 804
 805         for i in range(n_verts):
 806             vg.add([i], lap[i], 'REPLACE')
 807         ob.vertex_groups.update()
 808         ob.data.update()
 809         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
 810         bm.free()
 811         return {'FINISHED'}
 812
 813
 814 class reaction_diffusion(Operator):
 815     bl_idname = "object.reaction_diffusion"
 816     bl_label = "Reaction Diffusion"
 817     bl_description = ("Run a Reaction-Diffusion based on existing Vertex Groups: A and B")
 818     bl_options = {'REGISTER', 'UNDO'}
 819
 820     steps : IntProperty(
 821         name="Steps", default=10, min=0, soft_max=50,
 822         description="Number of Steps")
 823
 824     dt : FloatProperty(
 825         name="dt", default=0.2, min=0, soft_max=0.2,
 826         description="Time Step")
 827
 828     diff_a : FloatProperty(
 829         name="Diff A", default=1, min=0, soft_max=2,
 830         description="Diffusion A")
 831
 832     diff_b : FloatProperty(
 833         name="Diff B", default=0.5, min=0, soft_max=2,
 834         description="Diffusion B")
 835
 836     f : FloatProperty(
 837         name="f", default=0.055, min=0, soft_min=0.01, soft_max=0.06, max=0.1, precision=4,
 838         description="Feed Rate")
 839
 840     k : FloatProperty(
 841         name="k", default=0.062, min=0, soft_min=0.035, soft_max=0.065, max=0.1, precision=4,
 842         description="Kill Rate")
 843
 844     bounds_string = ""
 845
 846     frame = None
 847
 848     @classmethod
 849     def poll(cls, context):
 850         return len(context.object.vertex_groups) > 0
 851
 852
 853     def execute(self, context):
 854         #bpy.app.handlers.frame_change_post.remove(reaction_diffusion_def)
 855         reaction_diffusion_add_handler(self, context)
 856         set_animatable_fix_handler(self, context)
 857         try: ob = context.object
 858         except:
 859             self.report({'ERROR'}, "Please select an Object")
 860             return {'CANCELLED'}
 861
 862         me = ob.data
 863         bm = bmesh.new()
 864         bm.from_mesh(me)
 865         bm.edges.ensure_lookup_table()
 866
 867         # store weight values
 868         a = []
 869         b = []
 870         for v in me.vertices:
 871             try:
 872                 a.append(ob.vertex_groups["A"].weight(v.index))
 873             except:
 874                 a.append(0)
 875             try:
 876                 b.append(ob.vertex_groups["B"].weight(v.index))
 877             except:
 878                 b.append(0)
 879
 880         a = array(a)
 881         b = array(b)
 882         f = self.f
 883         k = self.k
 884         diff_a = self.diff_a
 885         diff_b = self.diff_b
 886         dt = self.dt
 887         n_verts = len(bm.verts)
 888
 889         for i in range(self.steps):
 890
 891             lap_a = zeros((n_verts))#[0]*n_verts
 892             lap_b = zeros((n_verts))#[0]*n_verts
 893             for e in bm.edges:
 894                 id0 = e.verts[0].index
 895                 id1 = e.verts[1].index
 896                 lap_a[id0] += a[id1] - a[id0]
 897                 lap_a[id1] += a[id0] - a[id1]
 898                 lap_b[id0] += b[id1] - b[id0]
 899                 lap_b[id1] += b[id0] - b[id1]
 900             ab2 = a*b**2
 901             a += (diff_a*lap_a - ab2 + f*(1-a))*dt
 902             b += (diff_b*lap_b + ab2 - (k+f)*b)*dt
 903
 904             for i in range(n_verts):
 905                 ob.vertex_groups['A'].add([i], a[i], 'REPLACE')
 906                 ob.vertex_groups['B'].add([i], b[i], 'REPLACE')
 907             ob.vertex_groups.update()
 908             ob.data.update()
 909
 910             bpy.ops.wm.redraw_timer(type='DRAW_WIN_SWAP', iterations=1)
 911
 912         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
 913         bm.free()
 914         return {'FINISHED'}
 915
 916
 917 class edges_deformation(Operator):
 918     bl_idname = "object.edges_deformation"
 919     bl_label = "Edges Deformation"
 920     bl_description = (
 921         "Compute Weight based on the deformation of edges "
 922         "according to visible modifiers"
 923     )
 924     bl_options = {'REGISTER', 'UNDO'}
 925
 926     bounds : EnumProperty(
 927         items=(('MANUAL', "Manual Bounds", ""),
 928             ('COMPRESSION', "Compressed Only", ""),
 929             ('TENSION', "Extended Only", ""),
 930             ('AUTOMATIC', "Automatic Bounds", "")),
 931         default='AUTOMATIC', name="Bounds")
 932
 933     mode : EnumProperty(
 934         items=(('MAX', "Max Deformation", ""),
 935             ('MEAN', "Average Deformation", "")),
 936         default='MEAN', name="Evaluation Mode")
 937
 938     min_def : FloatProperty(
 939         name="Min", default=0, soft_min=-1, soft_max=0,
 940         description="Deformations with 0 weight")
 941
 942     max_def : FloatProperty(
 943         name="Max", default=0.5, soft_min=0, soft_max=5,
 944         description="Deformations with 1 weight")
 945
 946     bounds_string = ""
 947
 948     frame = None
 949
 950     @classmethod
 951     def poll(cls, context):
 952         return len(context.object.modifiers) > 0
 953
 954     def draw(self, context):
 955         layout = self.layout
 956         col = layout.column(align=True)
 957         col.label(text="Evaluation Mode")
 958         col.prop(self, "mode", text="")
 959         col.label(text="Bounds")
 960         col.prop(self, "bounds", text="")
 961         if self.bounds == 'MANUAL':
 962             col.label(text="Strain Rate \u03B5:")
 963             col.prop(self, "min_def")
 964             col.prop(self, "max_def")
 965         col.label(text="\u03B5" + ": from " + self.bounds_string)
 966
 967     def execute(self, context):
 968         try: ob = context.object
 969         except:
 970             self.report({'ERROR'}, "Please select an Object")
 971             return {'CANCELLED'}
 972
 973         # check if the object is Cloth or Softbody
 974         physics = False
 975         for m in ob.modifiers:
 976             if m.type == 'CLOTH' or m.type == 'SOFT_BODY':
 977                 physics = True
 978                 if context.scene.frame_current == 1 and self.frame != None:
 979                     context.scene.frame_current = self.frame
 980                 break
 981         if not physics: self.frame = None
 982
 983         if self.mode == 'MEAN': group_name = "Average Deformation"
 984         elif self.mode == 'MAX': group_name = "Max Deformation"
 985         ob.vertex_groups.new(name=group_name)
 986         me0 = ob.data
 987
 988         me = simple_to_mesh(ob) #ob.to_mesh(preserve_all_data_layers=True, depsgraph=bpy.context.evaluated_depsgraph_get()).copy()
 989         if len(me.vertices) != len(me0.vertices) or len(me.edges) != len(me0.edges):
 990             self.report({'ERROR'}, "The topology of the object should be" +
 991                 "unaltered")
 992             return {'CANCELLED'}
 993
 994         bm0 = bmesh.new()
 995         bm0.from_mesh(me0)
 996         bm = bmesh.new()
 997         bm.from_mesh(me)
 998         deformations = []
 999         for e0, e in zip(bm0.edges, bm.edges):
1000             try:
1001                 l0 = e0.calc_length()
1002                 l1 = e.calc_length()
1003                 epsilon = (l1 - l0)/l0
1004                 deformations.append(epsilon)
1005             except: deformations.append(1)
1006         v_deformations = []
1007         for v in bm.verts:
1008             vdef = []
1009             for e in v.link_edges:
1010                 vdef.append(deformations[e.index])
1011             if self.mode == 'MEAN': v_deformations.append(mean(vdef))
1012             elif self.mode == 'MAX': v_deformations.append(max(vdef, key=abs))
1013             #elif self.mode == 'MIN': v_deformations.append(min(vdef, key=abs))
1014
1015         if self.bounds == 'MANUAL':
1016             min_def = self.min_def
1017             max_def = self.max_def
1018         elif self.bounds == 'AUTOMATIC':
1019             min_def = min(v_deformations)
1020             max_def = max(v_deformations)
1021             self.min_def = min_def
1022             self.max_def = max_def
1023         elif self.bounds == 'COMPRESSION':
1024             min_def = 0
1025             max_def = min(v_deformations)
1026             self.min_def = min_def
1027             self.max_def = max_def
1028         elif self.bounds == 'TENSION':
1029             min_def = 0
1030             max_def = max(v_deformations)
1031             self.min_def = min_def
1032             self.max_def = max_def
1033         delta_def = max_def - min_def
1034
1035         # check undeformed errors
1036         if delta_def == 0:
1037             if self.bounds == 'MANUAL':
1038                 delta_def = 0.0001
1039             else:
1040                 message = "The object doesn't have deformations."
1041                 if physics:
1042                     message = message + ("\nIf you are using Physics try to " +
1043                         "save it in the cache before.")
1044                 self.report({'ERROR'}, message)
1045                 return {'CANCELLED'}
1046         else:
1047             if physics:
1048                 self.frame = context.scene.frame_current
1049
1050         for i in range(len(v_deformations)):
1051             weight = (v_deformations[i] - min_def)/delta_def
1052             ob.vertex_groups[-1].add([i], weight, 'REPLACE')
1053         self.bounds_string = str(round(min_def,2)) + " to " + str(round(max_def,2))
1054         ob.vertex_groups[-1].name = group_name + " " + self.bounds_string
1055         ob.vertex_groups.update()
1056         ob.data.update()
1057         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
1058         bpy.data.meshes.remove(me)
1059         bm.free()
1060         bm0.free()
1061         return {'FINISHED'}
1062
1063 class edges_bending(Operator):
1064     bl_idname = "object.edges_bending"
1065     bl_label = "Edges Bending"
1066     bl_description = (
1067         "Compute Weight based on the bending of edges"
1068         "according to visible modifiers"
1069     )
1070     bl_options = {'REGISTER', 'UNDO'}
1071
1072     bounds : EnumProperty(
1073         items=(('MANUAL', "Manual Bounds", ""),
1074             ('POSITIVE', "Positive Only", ""),
1075             ('NEGATIVE', "Negative Only", ""),
1076             ('UNSIGNED', "Absolute Bending", ""),
1077             ('AUTOMATIC', "Signed Bending", "")),
1078         default='AUTOMATIC', name="Bounds")
1079
1080     min_def : FloatProperty(
1081         name="Min", default=-10, soft_min=-45, soft_max=45,
1082         description="Deformations with 0 weight")
1083
1084     max_def : FloatProperty(
1085         name="Max", default=10, soft_min=-45, soft_max=45,
1086         description="Deformations with 1 weight")
1087
1088     bounds_string = ""
1089     frame = None
1090
1091     @classmethod
1092     def poll(cls, context):
1093         return len(context.object.modifiers) > 0
1094
1095     def draw(self, context):
1096         layout = self.layout
1097         layout.label(text="Bounds")
1098         layout.prop(self, "bounds", text="")
1099         if self.bounds == 'MANUAL':
1100             layout.prop(self, "min_def")
1101             layout.prop(self, "max_def")
1102
1103     def execute(self, context):
1104         try: ob = context.object
1105         except:
1106             self.report({'ERROR'}, "Please select an Object")
1107             return {'CANCELLED'}
1108
1109         group_name = "Edges Bending"
1110         ob.vertex_groups.new(name=group_name)
1111
1112         # check if the object is Cloth or Softbody
1113         physics = False
1114         for m in ob.modifiers:
1115             if m.type == 'CLOTH' or m.type == 'SOFT_BODY':
1116                 physics = True
1117                 if context.scene.frame_current == 1 and self.frame != None:
1118                     context.scene.frame_current = self.frame
1119                 break
1120         if not physics: self.frame = None
1121
1122         #ob.data.update()
1123         #context.scene.update()
1124         me0 = ob.data
1125         me = simple_to_mesh(ob) #ob.to_mesh(preserve_all_data_layers=True, depsgraph=bpy.context.evaluated_depsgraph_get()).copy()
1126         if len(me.vertices) != len(me0.vertices) or len(me.edges) != len(me0.edges):
1127             self.report({'ERROR'}, "The topology of the object should be" +
1128                 "unaltered")
1129         bm0 = bmesh.new()
1130         bm0.from_mesh(me0)
1131         bm = bmesh.new()
1132         bm.from_mesh(me)
1133         deformations = []
1134         for e0, e in zip(bm0.edges, bm.edges):
1135             try:
1136                 ang = e.calc_face_angle_signed()
1137                 ang0 = e0.calc_face_angle_signed()
1138                 if self.bounds == 'UNSIGNED':
1139                     deformations.append(abs(ang-ang0))
1140                 else:
1141                     deformations.append(ang-ang0)
1142             except: deformations.append(0)
1143         v_deformations = []
1144         for v in bm.verts:
1145             vdef = []
1146             for e in v.link_edges:
1147                 vdef.append(deformations[e.index])
1148             v_deformations.append(mean(vdef))
1149         if self.bounds == 'MANUAL':
1150             min_def = radians(self.min_def)
1151             max_def = radians(self.max_def)
1152         elif self.bounds == 'AUTOMATIC':
1153             min_def = min(v_deformations)
1154             max_def = max(v_deformations)
1155         elif self.bounds == 'POSITIVE':
1156             min_def = 0
1157             max_def = min(v_deformations)
1158         elif self.bounds == 'NEGATIVE':
1159             min_def = 0
1160             max_def = max(v_deformations)
1161         elif self.bounds == 'UNSIGNED':
1162             min_def = 0
1163             max_def = max(v_deformations)
1164         delta_def = max_def - min_def
1165
1166         # check undeformed errors
1167         if delta_def == 0:
1168             if self.bounds == 'MANUAL':
1169                 delta_def = 0.0001
1170             else:
1171                 message = "The object doesn't have deformations."
1172                 if physics:
1173                     message = message + ("\nIf you are using Physics try to " +
1174                         "save it in the cache before.")
1175                 self.report({'ERROR'}, message)
1176                 return {'CANCELLED'}
1177         else:
1178             if physics:
1179                 self.frame = context.scene.frame_current
1180
1181         for i in range(len(v_deformations)):
1182             weight = (v_deformations[i] - min_def)/delta_def
1183             ob.vertex_groups[-1].add([i], weight, 'REPLACE')
1184         self.bounds_string = str(round(min_def,2)) + " to " + str(round(max_def,2))
1185         ob.vertex_groups[-1].name = group_name + " " + self.bounds_string
1186         ob.vertex_groups.update()
1187         ob.data.update()
1188         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
1189         bpy.data.meshes.remove(me)
1190         bm0.free()
1191         bm.free()
1192         return {'FINISHED'}
1193
1194 class weight_contour_displace(Operator):
1195     bl_idname = "object.weight_contour_displace"
1196     bl_label = "Contour Displace"
1197     bl_description = ("")
1198     bl_options = {'REGISTER', 'UNDO'}
1199
1200     use_modifiers : BoolProperty(
1201         name="Use Modifiers", default=True,
1202         description="Apply all the modifiers")
1203     min_iso : FloatProperty(
1204         name="Min Iso Value", default=0.49, min=0, max=1,
1205         description="Threshold value")
1206     max_iso : FloatProperty(
1207         name="Max Iso Value", default=0.51, min=0, max=1,
1208         description="Threshold value")
1209     n_cuts : IntProperty(
1210         name="Cuts", default=2, min=1, soft_max=10,
1211         description="Number of cuts in the selected range of values")
1212     bool_displace : BoolProperty(
1213         name="Add Displace", default=True, description="Add Displace Modifier")
1214     bool_flip : BoolProperty(
1215         name="Flip", default=False, description="Flip Output Weight")
1216
1217     weight_mode : EnumProperty(
1218         items=[('Remapped', 'Remapped', 'Remap values'),
1219                ('Alternate', 'Alternate', 'Alternate 0 and 1'),
1220                ('Original', 'Original', 'Keep original Vertex Group')],
1221         name="Weight", description="Choose how to convert vertex group",
1222         default="Remapped", options={'LIBRARY_EDITABLE'})
1223
1224     @classmethod
1225     def poll(cls, context):
1226         return len(context.object.vertex_groups) > 0
1227
1228     def invoke(self, context, event):
1229         return context.window_manager.invoke_props_dialog(self, width=350)
1230
1231     def execute(self, context):
1232         start_time = timeit.default_timer()
1233         try:
1234             check = context.object.vertex_groups[0]
1235         except:
1236             self.report({'ERROR'}, "The object doesn't have Vertex Groups")
1237             return {'CANCELLED'}
1238
1239         ob0 = context.object
1240
1241         group_id = ob0.vertex_groups.active_index
1242         vertex_group_name = ob0.vertex_groups[group_id].name
1243
1244         bpy.ops.object.mode_set(mode='EDIT')
1245         bpy.ops.mesh.select_all(action='SELECT')
1246         bpy.ops.object.mode_set(mode='OBJECT')
1247         if self.use_modifiers:
1248             #me0 = ob0.to_mesh(preserve_all_data_layers=True, depsgraph=bpy.context.evaluated_depsgraph_get()).copy()
1249             me0 = simple_to_mesh(ob0)
1250         else:
1251             me0 = ob0.data.copy()
1252
1253         # generate new bmesh
1254         bm = bmesh.new()
1255         bm.from_mesh(me0)
1256         bm.verts.ensure_lookup_table()
1257         bm.edges.ensure_lookup_table()
1258         bm.faces.ensure_lookup_table()
1259
1260         # store weight values
1261         weight = []
1262         ob = bpy.data.objects.new("temp", me0)
1263         for g in ob0.vertex_groups:
1264             ob.vertex_groups.new(name=g.name)
1265         for v in me0.vertices:
1266             try:
1267                 weight.append(ob.vertex_groups[vertex_group_name].weight(v.index))
1268             except:
1269                 weight.append(0)
1270
1271         # define iso values
1272         iso_values = []
1273         for i_cut in range(self.n_cuts):
1274             delta_iso = abs(self.max_iso - self.min_iso)
1275             min_iso = min(self.min_iso, self.max_iso)
1276             max_iso = max(self.min_iso, self.max_iso)
1277             if delta_iso == 0: iso_val = min_iso
1278             elif self.n_cuts > 1: iso_val = i_cut/(self.n_cuts-1)*delta_iso + min_iso
1279             else: iso_val = (self.max_iso + self.min_iso)/2
1280             iso_values.append(iso_val)
1281
1282         # Start Cuts Iterations
1283         filtered_edges = bm.edges
1284         for iso_val in iso_values:
1285             delete_edges = []
1286
1287             faces_mask = []
1288             for f in bm.faces:
1289                 w_min = 2
1290                 w_max = 2
1291                 for v in f.verts:
1292                     w = weight[v.index]
1293                     if w_min == 2:
1294                         w_max = w_min = w
1295                     if w > w_max: w_max = w
1296                     if w < w_min: w_min = w
1297                     if w_min < iso_val and w_max > iso_val:
1298                         faces_mask.append(f)
1299                         break
1300
1301             #link_faces = [[f for f in e.link_faces] for e in bm.edges]
1302
1303             #faces_todo = [f.select for f in bm.faces]
1304             #faces_todo = [True for f in bm.faces]
1305             verts = []
1306             edges = []
1307             edges_id = {}
1308             _filtered_edges = []
1309             n_verts = len(bm.verts)
1310             count = n_verts
1311             for e in filtered_edges:
1312                 #id0 = e.vertices[0]
1313                 #id1 = e.vertices[1]
1314                 id0 = e.verts[0].index
1315                 id1 = e.verts[1].index
1316                 w0 = weight[id0]
1317                 w1 = weight[id1]
1318
1319                 if w0 == w1: continue
1320                 elif w0 > iso_val and w1 > iso_val:
1321                     _filtered_edges.append(e)
1322                     continue
1323                 elif w0 < iso_val and w1 < iso_val: continue
1324                 elif w0 == iso_val or w1 == iso_val:
1325                     _filtered_edges.append(e)
1326                     continue
1327                 else:
1328                     v0 = bm.verts[id0].co
1329                     v1 = bm.verts[id1].co
1330                     v = v0.lerp(v1, (iso_val-w0)/(w1-w0))
1331                     if e not in delete_edges:
1332                         delete_edges.append(e)
1333                     verts.append(v)
1334                     edges_id[str(id0)+"_"+str(id1)] = count
1335                     edges_id[str(id1)+"_"+str(id0)] = count
1336                     count += 1
1337                     _filtered_edges.append(e)
1338             filtered_edges = _filtered_edges
1339             splitted_faces = []
1340
1341             switch = False
1342             # splitting faces
1343             for f in faces_mask:
1344                 # create sub-faces slots. Once a new vertex is reached it will
1345                 # change slot, storing the next vertices for a new face.
1346                 build_faces = [[],[]]
1347                 #switch = False
1348                 verts0 = [v.index for v in f.verts]
1349                 verts1 = list(verts0)
1350                 verts1.append(verts1.pop(0)) # shift list
1351                 for id0, id1 in zip(verts0, verts1):
1352
1353                     # add first vertex to active slot
1354                     build_faces[switch].append(id0)
1355
1356                     # try to split edge
1357                     try:
1358                         # check if the edge must be splitted
1359                         new_vert = edges_id[str(id0)+"_"+str(id1)]
1360                         # add new vertex
1361                         build_faces[switch].append(new_vert)
1362                         # if there is an open face on the other slot
1363                         if len(build_faces[not switch]) > 0:
1364                             # store actual face
1365                             splitted_faces.append(build_faces[switch])
1366                             # reset actual faces and switch
1367                             build_faces[switch] = []
1368                             # change face slot
1369                         switch = not switch
1370                         # continue previous face
1371                         build_faces[switch].append(new_vert)
1372                     except: pass
1373                 if len(build_faces[not switch]) == 2:
1374                     build_faces[not switch].append(id0)
1375                 if len(build_faces[not switch]) > 2:
1376                     splitted_faces.append(build_faces[not switch])
1377                 # add last face
1378                 splitted_faces.append(build_faces[switch])
1379                 #del_faces.append(f.index)
1380
1381             # adding new vertices
1382             _new_vert = bm.verts.new
1383             for v in verts: new_vert = _new_vert(v)
1384             bm.verts.index_update()
1385             bm.verts.ensure_lookup_table()
1386             # adding new faces
1387             _new_face = bm.faces.new
1388             missed_faces = []
1389             added_faces = []
1390             for f in splitted_faces:
1391                 try:
1392                     face_verts = [bm.verts[i] for i in f]
1393                     new_face = _new_face(face_verts)
1394                     for e in new_face.edges:
1395                         filtered_edges.append(e)
1396                 except:
1397                     missed_faces.append(f)
1398
1399             bm.faces.ensure_lookup_table()
1400             # updating weight values
1401             weight = weight + [iso_val]*len(verts)
1402
1403             # deleting old edges/faces
1404             _remove_edge = bm.edges.remove
1405             bm.edges.ensure_lookup_table()
1406             for e in delete_edges:
1407                 _remove_edge(e)
1408             _filtered_edges = []
1409             for e in filtered_edges:
1410                 if e not in delete_edges: _filtered_edges.append(e)
1411             filtered_edges = _filtered_edges
1412
1413         name = ob0.name + '_ContourDisp'
1414         me = bpy.data.meshes.new(name)
1415         bm.to_mesh(me)
1416         bm.free()
1417         ob = bpy.data.objects.new(name, me)
1418
1419         # Link object to scene and make active
1420         scn = context.scene
1421         context.collection.objects.link(ob)
1422         context.view_layer.objects.active = ob
1423         ob.select_set(True)
1424         ob0.select_set(False)
1425
1426         # generate new vertex group
1427         for g in ob0.vertex_groups:
1428             ob.vertex_groups.new(name=g.name)
1429         #ob.vertex_groups.new(name=vertex_group_name)
1430
1431         all_weight = weight + [iso_val]*len(verts)
1432         #mult = 1/(1-iso_val)
1433         for id in range(len(all_weight)):
1434             #if False: w = (all_weight[id]-iso_val)*mult
1435             w = all_weight[id]
1436             if self.weight_mode == 'Alternate':
1437                 direction = self.bool_flip
1438                 for i in range(len(iso_values)-1):
1439                     val0, val1 = iso_values[i], iso_values[i+1]
1440                     if val0 < w <= val1:
1441                         if direction: w1 = (w-val0)/(val1-val0)
1442                         else: w1 = (val1-w)/(val1-val0)
1443                     direction = not direction
1444                 if w < iso_values[0]: w1 = not self.bool_flip
1445                 if w > iso_values[-1]: w1 = not direction
1446             elif self.weight_mode == 'Remapped':
1447                 if w < min_iso: w1 = 0
1448                 elif w > max_iso: w1 = 1
1449                 else: w1 = (w - min_iso)/delta_iso
1450             else:
1451                 if self.bool_flip: w1 = 1-w
1452                 else: w1 = w
1453             ob.vertex_groups[vertex_group_name].add([id], w1, 'REPLACE')
1454
1455         ob.vertex_groups.active_index = group_id
1456
1457         # align new object
1458         ob.matrix_world = ob0.matrix_world
1459
1460         # Displace Modifier
1461         if self.bool_displace:
1462             displace_modifier = ob.modifiers.new(type='DISPLACE', name='Displace')
1463             displace_modifier.mid_level = 0
1464             displace_modifier.strength = 0.1
1465             displace_modifier.vertex_group = vertex_group_name
1466
1467         bpy.ops.object.mode_set(mode='EDIT')
1468         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
1469         print("Contour Displace time: " + str(timeit.default_timer() - start_time) + " sec")
1470
1471         bpy.data.meshes.remove(me0)
1472
1473         return {'FINISHED'}
1474
1475 class weight_contour_mask(Operator):
1476     bl_idname = "object.weight_contour_mask"
1477     bl_label = "Contour Mask"
1478     bl_description = ("")
1479     bl_options = {'REGISTER', 'UNDO'}
1480
1481     use_modifiers : BoolProperty(
1482         name="Use Modifiers", default=True,
1483         description="Apply all the modifiers")
1484     iso : FloatProperty(
1485         name="Iso Value", default=0.5, soft_min=0, soft_max=1,
1486         description="Threshold value")
1487     bool_solidify : BoolProperty(
1488         name="Solidify", default=True, description="Add Solidify Modifier")
1489     offset : FloatProperty(
1490         name="Offset", default=1, min=0, max=1,
1491         description="Offset")
1492     thickness : FloatProperty(
1493         name="Thickness", default=0.5, soft_min=0, soft_max=1,
1494         description="Thickness")
1495     normalize_weight : BoolProperty(
1496         name="Normalize Weight", default=True,
1497         description="Normalize weight of remaining vertices")
1498
1499     @classmethod
1500     def poll(cls, context):
1501         return len(context.object.vertex_groups) > 0
1502
1503     def invoke(self, context, event):
1504         return context.window_manager.invoke_props_dialog(self, width=350)
1505
1506     def execute(self, context):
1507         start_time = timeit.default_timer()
1508         try:
1509             check = context.object.vertex_groups[0]
1510         except:
1511             self.report({'ERROR'}, "The object doesn't have Vertex Groups")
1512             return {'CANCELLED'}
1513
1514         ob0 = bpy.context.object
1515
1516         iso_val = self.iso
1517         group_id = ob0.vertex_groups.active_index
1518         vertex_group_name = ob0.vertex_groups[group_id].name
1519
1520         bpy.ops.object.mode_set(mode='EDIT')
1521         bpy.ops.mesh.select_all(action='SELECT')
1522         bpy.ops.object.mode_set(mode='OBJECT')
1523         if self.use_modifiers:
1524             me0 = simple_to_mesh(ob0)#ob0.to_mesh(preserve_all_data_layers=True, depsgraph=bpy.context.evaluated_depsgraph_get()).copy()
1525         else:
1526             me0 = ob0.data.copy()
1527
1528         # generate new bmesh
1529         bm = bmesh.new()
1530         bm.from_mesh(me0)
1531         bm.verts.ensure_lookup_table()
1532         bm.edges.ensure_lookup_table()
1533         bm.faces.ensure_lookup_table()
1534
1535         # store weight values
1536         weight = []
1537         ob = bpy.data.objects.new("temp", me0)
1538         for g in ob0.vertex_groups:
1539             ob.vertex_groups.new(name=g.name)
1540         for v in me0.vertices:
1541             try:
1542                 #weight.append(v.groups[vertex_group_name].weight)
1543                 weight.append(ob.vertex_groups[vertex_group_name].weight(v.index))
1544             except:
1545                 weight.append(0)
1546
1547         faces_mask = []
1548         for f in bm.faces:
1549             w_min = 2
1550             w_max = 2
1551             for v in f.verts:
1552                 w = weight[v.index]
1553                 if w_min == 2:
1554                     w_max = w_min = w
1555                 if w > w_max: w_max = w
1556                 if w < w_min: w_min = w
1557                 if w_min < iso_val and w_max > iso_val:
1558                     faces_mask.append(f)
1559                     break
1560
1561         filtered_edges = bm.edges# me0.edges
1562         faces_todo = [f.select for f in bm.faces]
1563         verts = []
1564         edges = []
1565         delete_edges = []
1566         edges_id = {}
1567         _filtered_edges = []
1568         n_verts = len(bm.verts)
1569         count = n_verts
1570         for e in filtered_edges:
1571             id0 = e.verts[0].index
1572             id1 = e.verts[1].index
1573             w0 = weight[id0]
1574             w1 = weight[id1]
1575
1576             if w0 == w1: continue
1577             elif w0 > iso_val and w1 > iso_val:
1578                 continue
1579             elif w0 < iso_val and w1 < iso_val: continue
1580             elif w0 == iso_val or w1 == iso_val: continue
1581             else:
1582                 v0 = me0.vertices[id0].co
1583                 v1 = me0.vertices[id1].co
1584                 v = v0.lerp(v1, (iso_val-w0)/(w1-w0))
1585                 delete_edges.append(e)
1586                 verts.append(v)
1587                 edges_id[str(id0)+"_"+str(id1)] = count
1588                 edges_id[str(id1)+"_"+str(id0)] = count
1589                 count += 1
1590
1591         splitted_faces = []
1592
1593         switch = False
1594         # splitting faces
1595         for f in faces_mask:
1596             # create sub-faces slots. Once a new vertex is reached it will
1597             # change slot, storing the next vertices for a new face.
1598             build_faces = [[],[]]
1599             #switch = False
1600             verts0 = list(me0.polygons[f.index].vertices)
1601             verts1 = list(verts0)
1602             verts1.append(verts1.pop(0)) # shift list
1603             for id0, id1 in zip(verts0, verts1):
1604
1605                 # add first vertex to active slot
1606                 build_faces[switch].append(id0)
1607
1608                 # try to split edge
1609                 try:
1610                     # check if the edge must be splitted
1611                     new_vert = edges_id[str(id0)+"_"+str(id1)]
1612                     # add new vertex
1613                     build_faces[switch].append(new_vert)
1614                     # if there is an open face on the other slot
1615                     if len(build_faces[not switch]) > 0:
1616                         # store actual face
1617                         splitted_faces.append(build_faces[switch])
1618                         # reset actual faces and switch
1619                         build_faces[switch] = []
1620                         # change face slot
1621                     switch = not switch
1622                     # continue previous face
1623                     build_faces[switch].append(new_vert)
1624                 except: pass
1625             if len(build_faces[not switch]) == 2:
1626                 build_faces[not switch].append(id0)
1627             if len(build_faces[not switch]) > 2:
1628                 splitted_faces.append(build_faces[not switch])
1629             # add last face
1630             splitted_faces.append(build_faces[switch])
1631
1632         # adding new vertices
1633         _new_vert = bm.verts.new
1634         for v in verts: _new_vert(v)
1635         bm.verts.ensure_lookup_table()
1636
1637         # deleting old edges/faces
1638         _remove_edge = bm.edges.remove
1639         bm.edges.ensure_lookup_table()
1640         remove_edges = []
1641         for e in delete_edges: _remove_edge(e)
1642
1643         bm.verts.ensure_lookup_table()
1644         # adding new faces
1645         _new_face = bm.faces.new
1646         missed_faces = []
1647         for f in splitted_faces:
1648             try:
1649                 face_verts = [bm.verts[i] for i in f]
1650                 _new_face(face_verts)
1651             except:
1652                 missed_faces.append(f)
1653
1654         # Mask geometry
1655         if(True):
1656             _remove_vert = bm.verts.remove
1657             all_weight = weight + [iso_val+0.0001]*len(verts)
1658             weight = []
1659             for w, v in zip(all_weight, bm.verts):
1660                 if w < iso_val: _remove_vert(v)
1661                 else: weight.append(w)
1662
1663         # Create mesh and object
1664         name = ob0.name + '_ContourMask_{:.3f}'.format(iso_val)
1665         me = bpy.data.meshes.new(name)
1666         bm.to_mesh(me)
1667         bm.free()
1668         ob = bpy.data.objects.new(name, me)
1669
1670         # Link object to scene and make active
1671         scn = context.scene
1672         context.collection.objects.link(ob)
1673         context.view_layer.objects.active = ob
1674         ob.select_set(True)
1675         ob0.select_set(False)
1676
1677         # generate new vertex group
1678         for g in ob0.vertex_groups:
1679             ob.vertex_groups.new(name=g.name)
1680
1681         if iso_val != 1: mult = 1/(1-iso_val)
1682         else: mult = 1
1683         for id in range(len(weight)):
1684             if self.normalize_weight: w = (weight[id]-iso_val)*mult
1685             else: w = weight[id]
1686             ob.vertex_groups[vertex_group_name].add([id], w, 'REPLACE')
1687         ob.vertex_groups.active_index = group_id
1688
1689         # align new object
1690         ob.matrix_world = ob0.matrix_world
1691
1692         # Add Solidify
1693         if self.bool_solidify and True:
1694             ob.modifiers.new(type='SOLIDIFY', name='Solidify')
1695             ob.modifiers['Solidify'].thickness = self.thickness
1696             ob.modifiers['Solidify'].offset = self.offset
1697             ob.modifiers['Solidify'].vertex_group = vertex_group_name
1698
1699         bpy.ops.object.mode_set(mode='EDIT')
1700         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
1701         print("Contour Mask time: " + str(timeit.default_timer() - start_time) + " sec")
1702
1703         bpy.data.meshes.remove(me0)
1704
1705         return {'FINISHED'}
1706
1707
1708 class weight_contour_mask_wip(Operator):
1709     bl_idname = "object.weight_contour_mask"
1710     bl_label = "Contour Mask"
1711     bl_description = ("")
1712     bl_options = {'REGISTER', 'UNDO'}
1713
1714     use_modifiers : BoolProperty(
1715         name="Use Modifiers", default=True,
1716         description="Apply all the modifiers")
1717     iso : FloatProperty(
1718         name="Iso Value", default=0.5, soft_min=0, soft_max=1,
1719         description="Threshold value")
1720     bool_solidify : BoolProperty(
1721         name="Solidify", default=True, description="Add Solidify Modifier")
1722     normalize_weight : BoolProperty(
1723         name="Normalize Weight", default=True,
1724         description="Normalize weight of remaining vertices")
1725
1726     @classmethod
1727     def poll(cls, context):
1728         return len(context.object.vertex_groups) > 0
1729
1730     def execute(self, context):
1731         start_time = timeit.default_timer()
1732         try:
1733             check = context.object.vertex_groups[0]
1734         except:
1735             self.report({'ERROR'}, "The object doesn't have Vertex Groups")
1736             return {'CANCELLED'}
1737
1738         ob0 = bpy.context.object
1739
1740         iso_val = self.iso
1741         group_id = ob0.vertex_groups.active_index
1742         vertex_group_name = ob0.vertex_groups[group_id].name
1743
1744         #bpy.ops.object.mode_set(mode='EDIT')
1745         #bpy.ops.mesh.select_all(action='SELECT')
1746         #bpy.ops.object.mode_set(mode='OBJECT')
1747         if self.use_modifiers:
1748             me0 = simple_to_mesh(ob0)#ob0.to_mesh(preserve_all_data_layers=True, depsgraph=bpy.context.evaluated_depsgraph_get()).copy()
1749         else:
1750             me0 = ob0.data.copy()
1751
1752         # generate new bmesh
1753         bm = bmesh.new()
1754         bm.from_mesh(me0)
1755
1756         # store weight values
1757         weight = []
1758         ob = bpy.data.objects.new("temp", me0)
1759         for g in ob0.vertex_groups:
1760             ob.vertex_groups.new(name=g.name)
1761         weight = get_weight_numpy(ob.vertex_groups[vertex_group_name], len(me0.vertices))
1762
1763         me0, bm, weight = contour_bmesh(me0, bm, weight, iso_val)
1764
1765         # Mask geometry
1766         mask = weight >= iso_val
1767         weight = weight[mask]
1768         mask = np.logical_not(mask)
1769         delete_verts = np.array(bm.verts)[mask]
1770         #for v in delete_verts: bm.verts.remove(v)
1771
1772         # Create mesh and object
1773         name = ob0.name + '_ContourMask_{:.3f}'.format(iso_val)
1774         me = bpy.data.meshes.new(name)
1775         bm.to_mesh(me)
1776         bm.free()
1777         ob = bpy.data.objects.new(name, me)
1778
1779         # Link object to scene and make active
1780         scn = context.scene
1781         context.collection.objects.link(ob)
1782         context.view_layer.objects.active = ob
1783         ob.select_set(True)
1784         ob0.select_set(False)
1785
1786         # generate new vertex group
1787         for g in ob0.vertex_groups:
1788             ob.vertex_groups.new(name=g.name)
1789
1790         if iso_val != 1: mult = 1/(1-iso_val)
1791         else: mult = 1
1792         for id in range(len(weight)):
1793             if self.normalize_weight: w = (weight[id]-iso_val)*mult
1794             else: w = weight[id]
1795             ob.vertex_groups[vertex_group_name].add([id], w, 'REPLACE')
1796         ob.vertex_groups.active_index = group_id
1797
1798         # align new object
1799         ob.matrix_world = ob0.matrix_world
1800
1801         # Add Solidify
1802         if self.bool_solidify and True:
1803             ob.modifiers.new(type='SOLIDIFY', name='Solidify')
1804             ob.modifiers['Solidify'].thickness = 0.05
1805             ob.modifiers['Solidify'].offset = 0
1806             ob.modifiers['Solidify'].vertex_group = vertex_group_name
1807
1808         bpy.ops.object.mode_set(mode='EDIT')
1809         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
1810         print("Contour Mask time: " + str(timeit.default_timer() - start_time) + " sec")
1811
1812         bpy.data.meshes.remove(me0)
1813
1814         return {'FINISHED'}
1815
1816
1817 class weight_contour_curves(Operator):
1818     bl_idname = "object.weight_contour_curves"
1819     bl_label = "Contour Curves"
1820     bl_description = ("")
1821     bl_options = {'REGISTER', 'UNDO'}
1822
1823     use_modifiers : BoolProperty(
1824         name="Use Modifiers", default=True,
1825         description="Apply all the modifiers")
1826
1827     min_iso : FloatProperty(
1828         name="Min Value", default=0., soft_min=0, soft_max=1,
1829         description="Minimum weight value")
1830     max_iso : FloatProperty(
1831         name="Max Value", default=1, soft_min=0, soft_max=1,
1832         description="Maximum weight value")
1833     n_curves : IntProperty(
1834         name="Curves", default=3, soft_min=1, soft_max=10,
1835         description="Number of Contour Curves")
1836
1837     min_rad : FloatProperty(
1838         name="Min Radius", default=1, soft_min=0, soft_max=1,
1839         description="Change radius according to Iso Value")
1840     max_rad : FloatProperty(
1841         name="Max Radius", default=1, soft_min=0, soft_max=1,
1842         description="Change radius according to Iso Value")
1843
1844     @classmethod
1845     def poll(cls, context):
1846         ob = context.object
1847         return len(ob.vertex_groups) > 0 or ob.type == 'CURVE'
1848
1849     def invoke(self, context, event):
1850         return context.window_manager.invoke_props_dialog(self, width=350)
1851
1852     def execute(self, context):
1853         start_time = timeit.default_timer()
1854         try:
1855             check = context.object.vertex_groups[0]
1856         except:
1857             self.report({'ERROR'}, "The object doesn't have Vertex Groups")
1858             return {'CANCELLED'}
1859         ob0 = context.object
1860
1861         group_id = ob0.vertex_groups.active_index
1862         vertex_group_name = ob0.vertex_groups[group_id].name
1863
1864         bpy.ops.object.mode_set(mode='EDIT')
1865         bpy.ops.mesh.select_all(action='SELECT')
1866         bpy.ops.object.mode_set(mode='OBJECT')
1867         if self.use_modifiers:
1868             me0 = simple_to_mesh(ob0) #ob0.to_mesh(preserve_all_data_layers=True, depsgraph=bpy.context.evaluated_depsgraph_get()).copy()
1869         else:
1870             me0 = ob0.data.copy()
1871
1872         # generate new bmesh
1873         bm = bmesh.new()
1874         bm.from_mesh(me0)
1875         bm.verts.ensure_lookup_table()
1876         bm.edges.ensure_lookup_table()
1877         bm.faces.ensure_lookup_table()
1878
1879         # store weight values
1880         weight = []
1881         ob = bpy.data.objects.new("temp", me0)
1882         for g in ob0.vertex_groups:
1883             ob.vertex_groups.new(name=g.name)
1884         weight = get_weight_numpy(ob.vertex_groups[vertex_group_name], len(bm.verts))
1885
1886         #filtered_edges = bm.edges
1887         total_verts = np.zeros((0,3))
1888         total_segments = []
1889         radius = []
1890
1891         # start iterate contours levels
1892         vertices = get_vertices_numpy(me0)
1893         filtered_edges = get_edges_id_numpy(me0)
1894
1895         faces_weight = [np.array([weight[v] for v in p.vertices]) for p in me0.polygons]
1896         fw_min = np.array([np.min(fw) for fw in faces_weight])
1897         fw_max = np.array([np.max(fw) for fw in faces_weight])
1898
1899         bm_faces = np.array(bm.faces)
1900
1901         ### Spiral
1902         normals = np.array([v.normal for v in me0.vertices])
1903
1904         for c in range(self.n_curves):
1905             min_iso = min(self.min_iso, self.max_iso)
1906             max_iso = max(self.min_iso, self.max_iso)
1907             try:
1908                 delta_iso = (max_iso-min_iso)/(self.n_curves-1)
1909                 iso_val = c*delta_iso + min_iso
1910                 if iso_val < 0: iso_val = (min_iso + max_iso)/2
1911             except:
1912                 iso_val = (min_iso + max_iso)/2
1913
1914             # remove passed faces
1915             bool_mask = iso_val < fw_max
1916             bm_faces = bm_faces[bool_mask]
1917             fw_min = fw_min[bool_mask]
1918             fw_max = fw_max[bool_mask]
1919
1920             # mask faces
1921             bool_mask = fw_min < iso_val
1922             faces_mask = bm_faces[bool_mask]
1923
1924             n_verts = len(bm.verts)
1925             count = len(total_verts)
1926
1927             # vertices indexes
1928             id0 = filtered_edges[:,0]
1929             id1 = filtered_edges[:,1]
1930             # vertices weight
1931             w0 = weight[id0]
1932             w1 = weight[id1]
1933             # weight condition
1934             bool_w0 = w0 < iso_val
1935             bool_w1 = w1 < iso_val
1936
1937             # mask all edges that have one weight value below the iso value
1938             mask_new_verts = np.logical_xor(bool_w0, bool_w1)
1939
1940             id0 = id0[mask_new_verts]
1941             id1 = id1[mask_new_verts]
1942             # filter arrays
1943             v0 = vertices[id0]
1944             v1 = vertices[id1]
1945             w0 = w0[mask_new_verts]
1946             w1 = w1[mask_new_verts]
1947             div = (w1-w0)
1948             if div == 0: div = 0.000001
1949
1950             param = np.expand_dims((iso_val-w0)/div,axis=1)
1951             verts = v0 + (v1-v0)*param
1952
1953             # indexes of edges with new vertices
1954             edges_index = filtered_edges[mask_new_verts][:,2]
1955             edges_id = {}
1956             for i, id in enumerate(edges_index): edges_id[id] = i+len(total_verts)
1957
1958             # remove all edges completely below the iso value
1959             mask_edges = np.logical_not(np.logical_and(bool_w0, bool_w1))
1960             filtered_edges = filtered_edges[mask_edges]
1961             if len(verts) == 0: continue
1962
1963             # finding segments
1964             segments = []
1965             for f in faces_mask:
1966                 seg = []
1967                 for e in f.edges:
1968                     try:
1969                         seg.append(edges_id[e.index])
1970                         if len(seg) == 2:
1971                             segments.append(seg)
1972                             seg = []
1973                     except: pass
1974
1975
1976             #curves_points_indexes = find_curves(segments)
1977             total_segments = total_segments + segments
1978             total_verts = np.concatenate((total_verts,verts))
1979
1980             if self.min_rad != self.max_rad:
1981                 try:
1982                     iso_rad = c*(self.max_rad-self.min_rad)/(self.n_curves-1)+self.min_rad
1983                     if iso_rad < 0: iso_rad = (self.min_rad + self.max_rad)/2
1984                 except:
1985                     iso_rad = (self.min_rad + self.max_rad)/2
1986                 radius = radius + [iso_rad]*len(verts)
1987         print("Contour Curves, computing time: " + str(timeit.default_timer() - start_time) + " sec")
1988         bm.free()
1989         bm = bmesh.new()
1990         # adding new vertices  _local for fast access
1991         _new_vert = bm.verts.new
1992         for v in total_verts: _new_vert(v)
1993         bm.verts.ensure_lookup_table()
1994
1995         # adding new edges
1996         _new_edge = bm.edges.new
1997         for s in total_segments:
1998             try:
1999                 pts = [bm.verts[i] for i in s]
2000                 _new_edge(pts)
2001             except: pass
2002
2003
2004         try:
2005             name = ob0.name + '_ContourCurves'
2006             me = bpy.data.meshes.new(name)
2007             bm.to_mesh(me)
2008             bm.free()
2009             ob = bpy.data.objects.new(name, me)
2010             # Link object to scene and make active
2011             scn = context.scene
2012             context.collection.objects.link(ob)
2013             context.view_layer.objects.active = ob
2014             ob.select_set(True)
2015             ob0.select_set(False)
2016
2017             print("Contour Curves, bmesh time: " + str(timeit.default_timer() - start_time) + " sec")
2018             bpy.ops.object.convert(target='CURVE')
2019             ob = context.object
2020             if not (self.min_rad == 0 and self.max_rad == 0):
2021                 if self.min_rad != self.max_rad:
2022                     count = 0
2023                     for s in ob.data.splines:
2024                         for p in s.points:
2025                             p.radius = radius[count]
2026                             count += 1
2027                 else:
2028                     for s in ob.data.splines:
2029                         for p in s.points:
2030                             p.radius = self.min_rad
2031             ob.data.bevel_depth = 0.01
2032             ob.data.fill_mode = 'FULL'
2033             ob.data.bevel_resolution = 3
2034         except:
2035             self.report({'ERROR'}, "There are no values in the chosen range")
2036             return {'CANCELLED'}
2037
2038         # align new object
2039         ob.matrix_world = ob0.matrix_world
2040         print("Contour Curves time: " + str(timeit.default_timer() - start_time) + " sec")
2041
2042         bpy.data.meshes.remove(me0)
2043         bpy.data.meshes.remove(me)
2044
2045         return {'FINISHED'}
2046
2047 class tissue_weight_contour_curves_pattern(Operator):
2048     bl_idname = "object.tissue_weight_contour_curves_pattern"
2049     bl_label = "Contour Curves"
2050     bl_description = ("")
2051     bl_options = {'REGISTER', 'UNDO'}
2052
2053     use_modifiers : BoolProperty(
2054         name="Use Modifiers", default=True,
2055         description="Apply all the modifiers")
2056
2057     auto_bevel : BoolProperty(
2058         name="Automatic Bevel", default=False,
2059         description="Bevel depends on weight density")
2060
2061     min_iso : FloatProperty(
2062         name="Min Value", default=0., soft_min=0, soft_max=1,
2063         description="Minimum weight value")
2064     max_iso : FloatProperty(
2065         name="Max Value", default=1, soft_min=0, soft_max=1,
2066         description="Maximum weight value")
2067     n_curves : IntProperty(
2068         name="Curves", default=10, soft_min=1, soft_max=100,
2069         description="Number of Contour Curves")
2070     min_rad = 1
2071     max_rad = 1
2072
2073     in_displace : FloatProperty(
2074         name="Displace A", default=0, soft_min=-10, soft_max=10,
2075         description="Pattern displace strength")
2076     out_displace : FloatProperty(
2077         name="Displace B", default=2, soft_min=-10, soft_max=10,
2078         description="Pattern displace strength")
2079
2080     in_steps : IntProperty(
2081         name="Steps A", default=1, min=0, soft_max=10,
2082         description="Number of layers to move inwards")
2083     out_steps : IntProperty(
2084         name="Steps B", default=1, min=0, soft_max=10,
2085         description="Number of layers to move outwards")
2086     limit_z : BoolProperty(
2087         name="Limit Z", default=False,
2088         description="Limit Pattern in Z")
2089
2090     merge : BoolProperty(
2091         name="Merge Vertices", default=True,
2092         description="Merge points")
2093     merge_thres : FloatProperty(
2094         name="Merge Threshold", default=0.01, min=0, soft_max=1,
2095         description="Minimum Curve Radius")
2096
2097     bevel_depth : FloatProperty(
2098         name="Bevel Depth", default=0, min=0, soft_max=1,
2099         description="")
2100     min_bevel_depth : FloatProperty(
2101         name="Min Bevel Depth", default=0.1, min=0, soft_max=1,
2102         description="")
2103     max_bevel_depth : FloatProperty(
2104         name="Max Bevel Depth", default=1, min=0, soft_max=1,
2105         description="")
2106     remove_open_curves : BoolProperty(
2107         name="Remove Open Curves", default=False,
2108         description="Remove Open Curves")
2109
2110     vertex_group_pattern : StringProperty(
2111         name="Displace", default='',
2112         description="Vertex Group used for pattern displace")
2113
2114     vertex_group_bevel : StringProperty(
2115         name="Bevel", default='',
2116         description="Variable Bevel depth")
2117
2118     object_name : StringProperty(
2119         name="Active Object", default='',
2120         description="")
2121
2122     try: vg_name = bpy.context.object.vertex_groups.active.name
2123     except: vg_name = ''
2124
2125     vertex_group_contour : StringProperty(
2126         name="Contour", default=vg_name,
2127         description="Vertex Group used for contouring")
2128     clean_distance : FloatProperty(
2129         name="Clean Distance", default=0, min=0, soft_max=10,
2130         description="Remove short segments")
2131
2132
2133     spiralized: BoolProperty(
2134         name='Spiralized', default=False,
2135         description='Create a Spiral Contour. Works better with dense meshes'
2136     )
2137     spiral_axis: FloatVectorProperty(
2138         name="Spiral Axis", default=(0,0,1),
2139         description="Axis of the Spiral (in local coordinates)"
2140     )
2141     spiral_rotation : FloatProperty(
2142         name="Spiral Rotation", default=0, min=0, max=2*pi,
2143         description=""
2144     )
2145
2146     @classmethod
2147     def poll(cls, context):
2148         ob = context.object
2149         return ob and len(ob.vertex_groups) > 0 or ob.type == 'CURVE'
2150
2151     def invoke(self, context, event):
2152         return context.window_manager.invoke_props_dialog(self, width=250)
2153
2154     def draw(self, context):
2155         if not context.object.type == 'CURVE':
2156             self.object_name = context.object.name
2157         ob = bpy.data.objects[self.object_name]
2158         if self.vertex_group_contour not in [vg.name for vg in ob.vertex_groups]:
2159             self.vertex_group_contour = ob.vertex_groups.active.name
2160         layout = self.layout
2161         col = layout.column(align=True)
2162         col.prop(self, "use_modifiers")
2163         col.label(text="Contour Curves:")
2164         col.prop_search(self, 'vertex_group_contour', ob, "vertex_groups", text='')
2165         row = col.row(align=True)
2166         row.prop(self,'min_iso')
2167         row.prop(self,'max_iso')
2168         col.prop(self,'n_curves')
2169         col.separator()
2170         col.label(text='Curves Bevel:')
2171         col.prop(self,'auto_bevel')
2172         if not self.auto_bevel:
2173             col.prop_search(self, 'vertex_group_bevel', ob, "vertex_groups", text='')
2174         if self.vertex_group_bevel != '' or self.auto_bevel:
2175             row = col.row(align=True)
2176             row.prop(self,'min_bevel_depth')
2177             row.prop(self,'max_bevel_depth')
2178         else:
2179             col.prop(self,'bevel_depth')
2180         col.separator()
2181
2182         col.label(text="Displace Pattern:")
2183         col.prop_search(self, 'vertex_group_pattern', ob, "vertex_groups", text='')
2184         if self.vertex_group_pattern != '':
2185             row = col.row(align=True)
2186             row.prop(self,'in_steps')
2187             row.prop(self,'out_steps')
2188             row = col.row(align=True)
2189             row.prop(self,'in_displace')
2190             row.prop(self,'out_displace')
2191             col.prop(self,'limit_z')
2192         col.separator()
2193         row=col.row(align=True)
2194         row.prop(self,'spiralized')
2195         row.label(icon='MOD_SCREW')
2196         if self.spiralized:
2197             #row=col.row(align=True)
2198             #row.prop(self,'spiral_axis')
2199             #col.separator()
2200             col.prop(self,'spiral_rotation')
2201         col.separator()
2202
2203         col.label(text='Clean Curves:')
2204         col.prop(self,'clean_distance')
2205         col.prop(self,'remove_open_curves')
2206
2207     def execute(self, context):
2208         n_curves = self.n_curves
2209         start_time = timeit.default_timer()
2210         try:
2211             check = context.object.vertex_groups[0]
2212         except:
2213             self.report({'ERROR'}, "The object doesn't have Vertex Groups")
2214             return {'CANCELLED'}
2215         ob0 = bpy.data.objects[self.object_name]
2216
2217         dg = context.evaluated_depsgraph_get()
2218         ob = ob0.evaluated_get(dg)
2219         me0 = ob.data
2220
2221         # generate new bmesh
2222         bm = bmesh.new()
2223         bm.from_mesh(me0)
2224         n_verts = len(bm.verts)
2225
2226         # store weight values
2227         try:
2228             weight = get_weight_numpy(ob.vertex_groups[self.vertex_group_contour], len(me0.vertices))
2229         except:
2230             bm.free()
2231             self.report({'ERROR'}, "Please select a Vertex Group for contouring")
2232             return {'CANCELLED'}
2233
2234         try:
2235             pattern_weight = get_weight_numpy(ob.vertex_groups[self.vertex_group_pattern], len(me0.vertices))
2236         except:
2237             #self.report({'WARNING'}, "There is no Vertex Group assigned to the pattern displace")
2238             pattern_weight = np.zeros(len(me0.vertices))
2239
2240         variable_bevel = False
2241         try:
2242             bevel_weight = get_weight_numpy(ob.vertex_groups[self.vertex_group_bevel], len(me0.vertices))
2243             variable_bevel = True
2244         except:
2245             bevel_weight = np.ones(len(me0.vertices))
2246
2247         if self.auto_bevel:
2248             # calc weight density
2249             bevel_weight = np.ones(len(me0.vertices))*10000
2250             bevel_weight = np.zeros(len(me0.vertices))
2251             edges_length = np.array([e.calc_length() for e in bm.edges])
2252             edges_dw = np.array([max(abs(weight[e.verts[0].index]-weight[e.verts[1].index]),0.000001) for e in bm.edges])
2253             dens = edges_length/edges_dw
2254             n_records = np.zeros(len(me0.vertices))
2255             for i, e in enumerate(bm.edges):
2256                 for v in e.verts:
2257                     id = v.index
2258                     #bevel_weight[id] = min(bevel_weight[id], dens[i])
2259                     bevel_weight[id] += dens[i]
2260                     n_records[id] += 1
2261             bevel_weight = bevel_weight/n_records
2262             bevel_weight = (bevel_weight - min(bevel_weight))/(max(bevel_weight) - min(bevel_weight))
2263             #bevel_weight = 1-bevel_weight
2264             variable_bevel = True
2265
2266         #filtered_edges = bm.edges
2267         total_verts = np.zeros((0,3))
2268         total_radii = np.zeros((0,1))
2269         total_segments = []# np.array([])
2270         radius = []
2271
2272         # start iterate contours levels
2273         vertices, normals = get_vertices_and_normals_numpy(me0)
2274         filtered_edges = get_edges_id_numpy(me0)
2275
2276
2277         min_iso = min(self.min_iso, self.max_iso)
2278         max_iso = max(self.min_iso, self.max_iso)
2279
2280         # Spiral
2281         if self.spiralized:
2282             nx = normals[:,0]
2283             ny = normals[:,1]
2284             ang = self.spiral_rotation + weight*pi*n_curves+arctan2(nx,ny)
2285             weight = sin(ang)/2+0.5
2286             n_curves = 1
2287
2288         if n_curves > 1:
2289             delta_iso = (max_iso-min_iso)/(n_curves-1)
2290
2291         else:
2292             delta_iso = None
2293
2294         faces_weight = [np.array([weight[v] for v in p.vertices]) for p in me0.polygons]
2295         fw_min = np.array([np.min(fw) for fw in faces_weight])
2296         fw_max = np.array([np.max(fw) for fw in faces_weight])
2297
2298         bm_faces = np.array(bm.faces)
2299
2300         #print("Contour Curves, data loaded: " + str(timeit.default_timer() - start_time) + " sec")
2301         step_time = timeit.default_timer()
2302         for c in range(n_curves):
2303             if delta_iso:
2304                 iso_val = c*delta_iso + min_iso
2305                 if iso_val < 0: iso_val = (min_iso + max_iso)/2
2306             else:
2307                 iso_val = (min_iso + max_iso)/2
2308
2309             #if c == 0 and self.auto_bevel:
2310
2311
2312             # remove passed faces
2313             bool_mask = iso_val < fw_max
2314             bm_faces = bm_faces[bool_mask]
2315             fw_min = fw_min[bool_mask]
2316             fw_max = fw_max[bool_mask]
2317
2318             # mask faces
2319             bool_mask = fw_min < iso_val
2320             faces_mask = bm_faces[bool_mask]
2321
2322             count = len(total_verts)
2323
2324             new_filtered_edges, edges_index, verts, bevel = contour_edges_pattern(self, c, len(total_verts), iso_val, vertices, normals, filtered_edges, weight, pattern_weight, bevel_weight)
2325
2326             if len(edges_index) > 0:
2327                 if self.auto_bevel and False:
2328                     bevel = 1-dens[edges_index]
2329                     bevel = bevel[:,np.newaxis]
2330                 if self.max_bevel_depth != self.min_bevel_depth:
2331                     min_radius = self.min_bevel_depth / max(0.0001,self.max_bevel_depth)
2332                     radii = min_radius + bevel*(1 - min_radius)
2333                 else:
2334                     radii = bevel
2335             else:
2336                 continue
2337
2338             if verts[0,0] == None: continue
2339             else: filtered_edges = new_filtered_edges
2340             edges_id = {}
2341             for i, id in enumerate(edges_index): edges_id[id] = i + count
2342
2343             if len(verts) == 0: continue
2344
2345             # finding segments
2346             segments = []
2347             for f in faces_mask:
2348                 seg = []
2349                 for e in f.edges:
2350                     try:
2351                         #seg.append(new_ids[np.where(edges_index == e.index)[0][0]])
2352                         seg.append(edges_id[e.index])
2353                         if len(seg) == 2:
2354                             segments.append(seg)
2355                             seg = []
2356                     except: pass
2357
2358             total_segments = total_segments + segments
2359             total_verts = np.concatenate((total_verts, verts))
2360             total_radii = np.concatenate((total_radii, radii))
2361
2362             if self.min_rad != self.max_rad:
2363                 try:
2364                     iso_rad = c*(self.max_rad-self.min_rad)/(self.n_curves-1)+self.min_rad
2365                     if iso_rad < 0: iso_rad = (self.min_rad + self.max_rad)/2
2366                 except:
2367                     iso_rad = (self.min_rad + self.max_rad)/2
2368                 radius = radius + [iso_rad]*len(verts)
2369         #print("Contour Curves, points computing: " + str(timeit.default_timer() - step_time) + " sec")
2370         step_time = timeit.default_timer()
2371
2372         if len(total_segments) > 0:
2373             step_time = timeit.default_timer()
2374             ordered_points = find_curves(total_segments, len(total_verts))
2375
2376             #print("Contour Curves, point ordered in: " + str(timeit.default_timer() - step_time) + " sec")
2377             step_time = timeit.default_timer()
2378             crv = curve_from_pydata(total_verts, total_radii, ordered_points, ob0.name + '_ContourCurves', self.remove_open_curves, merge_distance=self.clean_distance)
2379             context.view_layer.objects.active = crv
2380             if variable_bevel: crv.data.bevel_depth = self.max_bevel_depth
2381             else: crv.data.bevel_depth = self.bevel_depth
2382
2383             crv.select_set(True)
2384             ob0.select_set(False)
2385             crv.matrix_world = ob0.matrix_world
2386             #print("Contour Curves, curves created in: " + str(timeit.default_timer() - step_time) + " sec")
2387         else:
2388             bm.free()
2389             self.report({'ERROR'}, "There are no values in the chosen range")
2390             return {'CANCELLED'}
2391         bm.free()
2392         print("Contour Curves, total time: " + str(timeit.default_timer() - start_time) + " sec")
2393         return {'FINISHED'}
2394
2395 class vertex_colors_to_vertex_groups(Operator):
2396     bl_idname = "object.vertex_colors_to_vertex_groups"
2397     bl_label = "Vertex Color"
2398     bl_options = {'REGISTER', 'UNDO'}
2399     bl_description = ("Convert the active Vertex Color into a Vertex Group")
2400
2401     red : BoolProperty(
2402         name="red channel", default=False, description="convert red channel")
2403     green : BoolProperty(
2404         name="green channel", default=False,
2405         description="convert green channel")
2406     blue : BoolProperty(
2407         name="blue channel", default=False, description="convert blue channel")
2408     value : BoolProperty(
2409         name="value channel", default=True, description="convert value channel")
2410     invert : BoolProperty(
2411          name="invert", default=False, description="invert all color channels")
2412
2413     @classmethod
2414     def poll(cls, context):
2415         try:
2416             return len(context.object.data.vertex_colors) > 0
2417         except: return False
2418
2419     def execute(self, context):
2420         obj = context.active_object
2421         id = len(obj.vertex_groups)
2422         id_red = id
2423         id_green = id
2424         id_blue = id
2425         id_value = id
2426
2427         boolCol = len(obj.data.vertex_colors)
2428         if(boolCol): col_name = obj.data.vertex_colors.active.name
2429         bpy.ops.object.mode_set(mode='EDIT')
2430         bpy.ops.mesh.select_all(action='SELECT')
2431
2432         if(self.red and boolCol):
2433             bpy.ops.object.vertex_group_add()
2434             bpy.ops.object.vertex_group_assign()
2435             id_red = id
2436             obj.vertex_groups[id_red].name = col_name + '_red'
2437             id+=1
2438         if(self.green and boolCol):
2439             bpy.ops.object.vertex_group_add()
2440             bpy.ops.object.vertex_group_assign()
2441             id_green = id
2442             obj.vertex_groups[id_green].name = col_name + '_green'
2443             id+=1
2444         if(self.blue and boolCol):
2445             bpy.ops.object.vertex_group_add()
2446             bpy.ops.object.vertex_group_assign()
2447             id_blue = id
2448             obj.vertex_groups[id_blue].name = col_name + '_blue'
2449             id+=1
2450         if(self.value and boolCol):
2451             bpy.ops.object.vertex_group_add()
2452             bpy.ops.object.vertex_group_assign()
2453             id_value = id
2454             obj.vertex_groups[id_value].name = col_name + '_value'
2455             id+=1
2456
2457         mult = 1
2458         if(self.invert): mult = -1
2459         bpy.ops.object.mode_set(mode='OBJECT')
2460         sub_red = 1 + self.value + self.blue + self.green
2461         sub_green = 1 + self.value + self.blue
2462         sub_blue = 1 + self.value
2463         sub_value = 1
2464
2465         id = len(obj.vertex_groups)
2466         if(id_red <= id and id_green <= id and id_blue <= id and id_value <= \
2467                 id and boolCol):
2468             v_colors = obj.data.vertex_colors.active.data
2469             i = 0
2470             for f in obj.data.polygons:
2471                 for v in f.vertices:
2472                     gr = obj.data.vertices[v].groups
2473                     if(self.red): gr[min(len(gr)-sub_red, id_red)].weight = \
2474                         self.invert + mult * v_colors[i].color[0]
2475                     if(self.green): gr[min(len(gr)-sub_green, id_green)].weight\
2476                         = self.invert + mult * v_colors[i].color[1]
2477                     if(self.blue): gr[min(len(gr)-sub_blue, id_blue)].weight = \
2478                         self.invert + mult * v_colors[i].color[2]
2479                     if(self.value):
2480                         r = v_colors[i].color[0]
2481                         g = v_colors[i].color[1]
2482                         b = v_colors[i].color[2]
2483                         gr[min(len(gr)-sub_value, id_value)].weight\
2484                         = self.invert + mult * (0.2126*r + 0.7152*g + 0.0722*b)
2485                     i+=1
2486             bpy.ops.paint.weight_paint_toggle()
2487         return {'FINISHED'}
2488
2489 class vertex_group_to_vertex_colors(Operator):
2490     bl_idname = "object.vertex_group_to_vertex_colors"
2491     bl_label = "Vertex Group"
2492     bl_options = {'REGISTER', 'UNDO'}
2493     bl_description = ("Convert the active Vertex Group into a Vertex Color")
2494
2495     channel : EnumProperty(
2496         items=[('BLUE', 'Blue Channel', 'Convert to Blue Channel'),
2497                ('GREEN', 'Green Channel', 'Convert to Green Channel'),
2498                ('RED', 'Red Channel', 'Convert to Red Channel'),
2499                ('VALUE', 'Value Channel', 'Convert to Grayscale'),
2500                ('FALSE_COLORS', 'False Colors', 'Convert to False Colors')],
2501         name="Convert to", description="Choose how to convert vertex group",
2502         default="VALUE", options={'LIBRARY_EDITABLE'})
2503
2504     invert : BoolProperty(
2505         name="invert", default=False, description="invert color channel")
2506
2507     @classmethod
2508     def poll(cls, context):
2509         return len(context.object.vertex_groups) > 0
2510
2511     def execute(self, context):
2512         obj = context.active_object
2513         me = obj.data
2514         group_id = obj.vertex_groups.active_index
2515         if (group_id == -1):
2516             return {'FINISHED'}
2517
2518         bpy.ops.object.mode_set(mode='OBJECT')
2519         group_name = obj.vertex_groups[group_id].name
2520         me.vertex_colors.new()
2521         colors_id = obj.data.vertex_colors.active_index
2522
2523         colors_name = group_name
2524         if(self.channel == 'FALSE_COLORS'): colors_name += "_false_colors"
2525         elif(self.channel == 'VALUE'):  colors_name += "_value"
2526         elif(self.channel == 'RED'):  colors_name += "_red"
2527         elif(self.channel == 'GREEN'):  colors_name += "_green"
2528         elif(self.channel == 'BLUE'):  colors_name += "_blue"
2529         context.object.data.vertex_colors[colors_id].name = colors_name
2530
2531         v_colors = obj.data.vertex_colors.active.data
2532
2533         bm = bmesh.new()
2534         bm.from_mesh(me)
2535         dvert_lay = bm.verts.layers.deform.active
2536         weight = bmesh_get_weight_numpy(group_id,dvert_lay,bm.verts)
2537         if self.invert: weight = 1-weight
2538         loops_size = get_attribute_numpy(me.polygons, attribute='loop_total', mult=1)
2539         n_colors = np.sum(loops_size)
2540         verts = np.ones(n_colors)
2541         me.polygons.foreach_get('vertices',verts)
2542         splitted_weight = weight[verts.astype(int)][:,None]
2543         r = np.zeros(splitted_weight.shape)
2544         g = np.zeros(splitted_weight.shape)
2545         b = np.zeros(splitted_weight.shape)
2546         a = np.ones(splitted_weight.shape)
2547         if(self.channel == 'FALSE_COLORS'):
2548             mult = 0.6+0.4*splitted_weight
2549             mask = splitted_weight < 0.25
2550             g[mask] = splitted_weight[mask]*4
2551             b[mask] = np.ones(splitted_weight.shape)[mask]
2552
2553             mask = np.where(np.logical_and(splitted_weight>=0.25, splitted_weight<0.5))
2554             g[mask] = np.ones(splitted_weight.shape)[mask]
2555             b[mask] = (1-(splitted_weight[mask]-0.25)*4)
2556
2557             mask = np.where(np.logical_and(splitted_weight>=0.5, splitted_weight<0.75))
2558             r[mask] = (splitted_weight[mask]-0.5)*4
2559             g[mask] = np.ones(splitted_weight.shape)[mask]
2560
2561             mask = 0.75 <= splitted_weight
2562             r[mask] = np.ones(splitted_weight.shape)[mask]
2563             g[mask] = (1-(splitted_weight[mask]-0.75)*4)
2564         elif(self.channel == 'VALUE'):
2565             r = splitted_weight
2566             g = splitted_weight
2567             b = splitted_weight
2568         elif(self.channel == 'RED'):
2569             r = splitted_weight
2570         elif(self.channel == 'GREEN'):
2571             g = splitted_weight
2572         elif(self.channel == 'BLUE'):
2573             b = splitted_weight
2574
2575         colors = np.concatenate((r,g,b,a),axis=1).flatten()
2576         v_colors.foreach_set('color',colors)
2577
2578         bpy.ops.paint.vertex_paint_toggle()
2579         context.object.data.vertex_colors[colors_id].active_render = True
2580         return {'FINISHED'}
2581
2582 class vertex_group_to_uv(Operator):
2583     bl_idname = "object.vertex_group_to_uv"
2584     bl_label = "Vertex Group"
2585     bl_options = {'REGISTER', 'UNDO'}
2586     bl_description = ("Combine two Vertex Groups as UV Map Layer")
2587
2588     vertex_group_u : StringProperty(
2589         name="U", default='',
2590         description="Vertex Group used for the U coordinate")
2591     vertex_group_v : StringProperty(
2592         name="V", default='',
2593         description="Vertex Group used for the V coordinate")
2594     normalize_weight : BoolProperty(
2595         name="Normalize Weight", default=True,
2596         description="Normalize weight values")
2597     invert_u : BoolProperty(
2598         name="Invert U", default=False, description="Invert U")
2599     invert_v : BoolProperty(
2600         name="Invert V", default=False, description="Invert V")
2601
2602     @classmethod
2603     def poll(cls, context):
2604         return len(context.object.vertex_groups) > 0
2605
2606     def invoke(self, context, event):
2607         return context.window_manager.invoke_props_dialog(self, width=250)
2608
2609     def draw(self, context):
2610         ob = context.object
2611         layout = self.layout
2612         col = layout.column(align=True)
2613         row = col.row(align=True)
2614         row.prop_search(self, 'vertex_group_u', ob, "vertex_groups", text='')
2615         row.separator()
2616         row.prop_search(self, 'vertex_group_v', ob, "vertex_groups", text='')
2617         row = col.row(align=True)
2618         row.prop(self, "invert_u")
2619         row.separator()
2620         row.prop(self, "invert_v")
2621         row = col.row(align=True)
2622         row.prop(self, "normalize_weight")
2623
2624     def execute(self, context):
2625         ob = context.active_object
2626         me = ob.data
2627         n_verts = len(me.vertices)
2628         vg_keys = ob.vertex_groups.keys()
2629         bool_u = self.vertex_group_u in vg_keys
2630         bool_v = self.vertex_group_v in vg_keys
2631         if bool_u or bool_v:
2632             bm = bmesh.new()
2633             bm.from_mesh(me)
2634             dvert_lay = bm.verts.layers.deform.active
2635             if bool_u:
2636                 u_index = ob.vertex_groups[self.vertex_group_u].index
2637                 u = bmesh_get_weight_numpy(u_index, dvert_lay, bm.verts)
2638                 if self.invert_u:
2639                     u = 1-u
2640                 if self.normalize_weight:
2641                     u = np.interp(u, (u.min(), u.max()), (0, 1))
2642             else:
2643                 u = np.zeros(n_verts)
2644             if bool_v:
2645                 v_index = ob.vertex_groups[self.vertex_group_v].index
2646                 v = bmesh_get_weight_numpy(v_index, dvert_lay, bm.verts)
2647                 if self.invert_v:
2648                     v = 1-v
2649                 if self.normalize_weight:
2650                     v = np.interp(v, (v.min(), v.max()), (0, 1))
2651             else:
2652                 v = np.zeros(n_verts)
2653         else:
2654             u = v = np.zeros(n_verts)
2655
2656         uv_layer = me.uv_layers.new(name='Weight_to_UV')
2657         loops_size = get_attribute_numpy(me.polygons, attribute='loop_total', mult=1)
2658         n_data = np.sum(loops_size)
2659         v_id = np.ones(n_data)
2660         me.polygons.foreach_get('vertices',v_id)
2661         v_id = v_id.astype(int)
2662         split_u = u[v_id,None]
2663         split_v = v[v_id,None]
2664         uv = np.concatenate((split_u,split_v),axis=1).flatten()
2665         uv_layer.data.foreach_set('uv',uv)
2666         me.uv_layers.update()
2667         return {'FINISHED'}
2668
2669 class curvature_to_vertex_groups(Operator):
2670     bl_idname = "object.curvature_to_vertex_groups"
2671     bl_label = "Curvature"
2672     bl_options = {'REGISTER', 'UNDO'}
2673     bl_description = ("Generate a Vertex Group based on the curvature of the"
2674                       "mesh. Is based on Dirty Vertex Color")
2675
2676     invert : BoolProperty(
2677         name="invert", default=False, description="invert values")
2678
2679     blur_strength : FloatProperty(
2680       name="Blur Strength", default=1, min=0.001,
2681       max=1, description="Blur strength per iteration")
2682
2683     blur_iterations : IntProperty(
2684       name="Blur Iterations", default=1, min=0,
2685       max=40, description="Number of times to blur the values")
2686
2687     min_angle : FloatProperty(
2688       name="Min Angle", default=0, min=0,
2689       max=pi/2, subtype='ANGLE', description="Minimum angle")
2690
2691     max_angle : FloatProperty(
2692       name="Max Angle", default=pi, min=pi/2,
2693       max=pi, subtype='ANGLE', description="Maximum angle")
2694
2695     invert : BoolProperty(
2696         name="Invert", default=False,
2697         description="Invert the curvature map")
2698
2699     def execute(self, context):
2700         bpy.ops.object.mode_set(mode='OBJECT')
2701         vertex_colors = context.active_object.data.vertex_colors
2702         vertex_colors.new()
2703         vertex_colors[-1].active = True
2704         vertex_colors[-1].active_render = True
2705         vertex_colors[-1].name = "Curvature"
2706         for c in vertex_colors[-1].data: c.color = (1,1,1,1)
2707         bpy.ops.object.mode_set(mode='VERTEX_PAINT')
2708         bpy.ops.paint.vertex_color_dirt(
2709             blur_strength=self.blur_strength,
2710             blur_iterations=self.blur_iterations, clean_angle=self.max_angle,
2711             dirt_angle=self.min_angle)
2712         bpy.ops.object.vertex_colors_to_vertex_groups(invert=self.invert)
2713         vertex_colors.remove(vertex_colors.active)
2714         return {'FINISHED'}
2715
2716 class face_area_to_vertex_groups(Operator):
2717     bl_idname = "object.face_area_to_vertex_groups"
2718     bl_label = "Area"
2719     bl_options = {'REGISTER', 'UNDO'}
2720     bl_description = ("Generate a Vertex Group based on the area of individual"
2721                       "faces")
2722
2723     invert : BoolProperty(
2724         name="invert", default=False, description="invert values")
2725     bounds : EnumProperty(
2726         items=(('MANUAL', "Manual Bounds", ""),
2727             ('AUTOMATIC', "Automatic Bounds", "")),
2728         default='AUTOMATIC', name="Bounds")
2729
2730     min_area : FloatProperty(
2731         name="Min", default=0.01, soft_min=0, soft_max=1,
2732         description="Faces with 0 weight")
2733
2734     max_area : FloatProperty(
2735         name="Max", default=0.1, soft_min=0, soft_max=1,
2736         description="Faces with 1 weight")
2737
2738     def draw(self, context):
2739         layout = self.layout
2740         layout.label(text="Bounds")
2741         layout.prop(self, "bounds", text="")
2742         if self.bounds == 'MANUAL':
2743             layout.prop(self, "min_area")
2744             layout.prop(self, "max_area")
2745
2746     def execute(self, context):
2747         try: ob = context.object
2748         except:
2749             self.report({'ERROR'}, "Please select an Object")
2750             return {'CANCELLED'}
2751         ob.vertex_groups.new(name="Faces Area")
2752
2753         areas = [[] for v in ob.data.vertices]
2754
2755         for p in ob.data.polygons:
2756             for v in p.vertices:
2757                 areas[v].append(p.area)
2758
2759         for i in range(len(areas)):
2760             areas[i] = mean(areas[i])
2761         if self.bounds == 'MANUAL':
2762             min_area = self.min_area
2763             max_area = self.max_area
2764         elif self.bounds == 'AUTOMATIC':
2765             min_area = min(areas)
2766             max_area = max(areas)
2767         elif self.bounds == 'COMPRESSION':
2768             min_area = 1
2769             max_area = min(areas)
2770         elif self.bounds == 'TENSION':
2771             min_area = 1
2772             max_area = max(areas)
2773         delta_area = max_area - min_area
2774         if delta_area == 0:
2775             delta_area = 0.0001
2776             if self.bounds == 'MANUAL':
2777                 delta_area = 0.0001
2778             else:
2779                 self.report({'ERROR'}, "The faces have the same areas")
2780                 #return {'CANCELLED'}
2781         for i in range(len(areas)):
2782             weight = (areas[i] - min_area)/delta_area
2783             ob.vertex_groups[-1].add([i], weight, 'REPLACE')
2784         ob.vertex_groups.update()
2785         ob.data.update()
2786         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
2787         return {'FINISHED'}
2788
2789 class random_weight(Operator):
2790     bl_idname = "object.random_weight"
2791     bl_label = "Random"
2792     bl_options = {'REGISTER', 'UNDO'}
2793     bl_description = ("Generate a random Vertex Group")
2794
2795     min_val : FloatProperty(
2796         name="Min", default=0, soft_min=0, soft_max=1,
2797         description="Minimum Value")
2798
2799     max_val : FloatProperty(
2800         name="Max", default=1, soft_min=0, soft_max=1,
2801         description="Maximum Value")
2802
2803     #def draw(self, context):
2804     #    layout = self.layout
2805     #    layout.prop(self, "min_area")
2806     #    layout.prop(self, "max_area")
2807     @classmethod
2808     def poll(cls, context):
2809         return len(context.object.vertex_groups) > 0
2810
2811     def execute(self, context):
2812         try: ob = context.object
2813         except:
2814             self.report({'ERROR'}, "Please select an Object")
2815             return {'CANCELLED'}
2816         #ob.vertex_groups.new(name="Random")
2817         n_verts = len(ob.data.vertices)
2818         weight = np.random.uniform(low=self.min_val, high=self.max_val, size=(n_verts,))
2819         np.clip(weight, 0, 1, out=weight)
2820
2821         group_id = ob.vertex_groups.active_index
2822         for i in range(n_verts):
2823             ob.vertex_groups[group_id].add([i], weight[i], 'REPLACE')
2824         ob.vertex_groups.update()
2825         ob.data.update()
2826         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
2827         return {'FINISHED'}
2828
2829
2830 class harmonic_weight(Operator):
2831     bl_idname = "object.harmonic_weight"
2832     bl_label = "Harmonic"
2833     bl_options = {'REGISTER', 'UNDO'}
2834     bl_description = ("Create an harmonic variation of the active Vertex Group")
2835
2836     freq : FloatProperty(
2837         name="Frequency", default=20, soft_min=0,
2838         soft_max=100, description="Wave frequency")
2839
2840     amp : FloatProperty(
2841         name="Amplitude", default=1, soft_min=0,
2842         soft_max=10, description="Wave amplitude")
2843
2844     midlevel : FloatProperty(
2845         name="Midlevel", default=0, min=-1,
2846         max=1, description="Midlevel")
2847
2848     add : FloatProperty(
2849         name="Add", default=0, min=-1,
2850         max=1, description="Add to the Weight")
2851
2852     mult : FloatProperty(
2853         name="Multiply", default=0, min=0,
2854         max=1, description="Multiply for he Weight")
2855
2856     @classmethod
2857     def poll(cls, context):
2858         return len(context.object.vertex_groups) > 0
2859
2860     def execute(self, context):
2861         ob = context.active_object
2862         if len(ob.vertex_groups) > 0:
2863             group_id = ob.vertex_groups.active_index
2864             ob.vertex_groups.new(name="Harmonic")
2865             for i in range(len(ob.data.vertices)):
2866                 try: val = ob.vertex_groups[group_id].weight(i)
2867                 except: val = 0
2868                 weight = self.amp*(math.sin(val*self.freq) - self.midlevel)/2 + 0.5 + self.add*val*(1-(1-val)*self.mult)
2869                 ob.vertex_groups[-1].add([i], weight, 'REPLACE')
2870             ob.data.update()
2871         else:
2872             self.report({'ERROR'}, "Active object doesn't have vertex groups")
2873             return {'CANCELLED'}
2874         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
2875         return {'FINISHED'}
2876
2877
2878 class tissue_weight_distance(Operator):
2879     bl_idname = "object.tissue_weight_distance"
2880     bl_label = "Weight Distance"
2881     bl_options = {'REGISTER', 'UNDO'}
2882     bl_description = ("Create a weight map according to the distance from the "
2883                     "selected vertices along the mesh surface")
2884
2885     mode : EnumProperty(
2886         items=(('GEOD', "Geodesic Distance", ""),
2887             ('EUCL', "Euclidean Distance", ""),
2888             ('TOPO', "Topology Distance", "")),
2889         default='GEOD', name="Distance Method")
2890
2891     normalize : BoolProperty(
2892         name="Normalize", default=True,
2893         description="Automatically remap the distance values from 0 to 1")
2894
2895     min_value : FloatProperty(
2896         name="Min", default=0, min=0,
2897         soft_max=100, description="Minimum Distance")
2898
2899     max_value : FloatProperty(
2900         name="Max", default=10, min=0,
2901         soft_max=100, description="Max Distance")
2902
2903     def invoke(self, context, event):
2904         return context.window_manager.invoke_props_dialog(self, width=250)
2905
2906     def fill_neighbors(self,verts,weight):
2907         neigh = {}
2908         for v0 in verts:
2909             for f in v0.link_faces:
2910                 for v1 in f.verts:
2911                     if self.mode == 'GEOD':
2912                         dist = weight[v0.index] + (v0.co-v1.co).length
2913                     elif self.mode == 'TOPO':
2914                         dist = weight[v0.index] + 1.0
2915                     w1 = weight[v1.index]
2916                     if w1 == None or w1 > dist:
2917                         weight[v1.index] = dist
2918                         neigh[v1] = 0
2919         if len(neigh) == 0: return weight
2920         else: return self.fill_neighbors(neigh.keys(), weight)
2921
2922     def execute(self, context):
2923         ob = context.object
2924         old_mode = ob.mode
2925         if old_mode != 'OBJECT':
2926             bpy.ops.object.mode_set(mode='OBJECT')
2927
2928         me = ob.data
2929
2930         # store weight values
2931         weight = [None]*len(me.vertices)
2932
2933         if self.mode != 'EUCL':
2934             bm = bmesh.new()
2935             bm.from_mesh(me)
2936             bm.verts.ensure_lookup_table()
2937             bm.edges.ensure_lookup_table()
2938             bm.faces.ensure_lookup_table()
2939             selected = [v for v in bm.verts if v.select]
2940             if len(selected) == 0:
2941                 bpy.ops.object.mode_set(mode=old_mode)
2942                 message = "Please, select one or more vertices"
2943                 self.report({'ERROR'}, message)
2944                 return {'CANCELLED'}
2945             for v in selected: weight[v.index] = 0
2946             weight = self.fill_neighbors(selected, weight)
2947             bm.free()
2948         else:
2949             selected = [v for v in me.vertices if v.select]
2950             kd = KDTree(len(selected))
2951             for i, v in enumerate(selected):
2952                 kd.insert(v.co, i)
2953             kd.balance()
2954             for i,v in enumerate(me.vertices):
2955                 co, index, dist = kd.find(v.co)
2956                 weight[i] = dist
2957
2958
2959         for i in range(len(weight)):
2960             if weight[i] == None: weight[i] = 0
2961         weight = np.array(weight)
2962         max_dist = np.max(weight)
2963         if self.normalize:
2964             if max_dist > 0:
2965                 weight /= max_dist
2966         else:
2967             delta_value = self.max_value - self.min_value
2968             if delta_value == 0: delta_value = 0.0000001
2969             weight = (weight-self.min_value)/delta_value
2970
2971         if self.mode == 'TOPO':
2972             vg = ob.vertex_groups.new(name='Distance: {:d}'.format(int(max_dist)))
2973         else:
2974             vg = ob.vertex_groups.new(name='Distance: {:.4f}'.format(max_dist))
2975         for i, w in enumerate(weight):
2976             vg.add([i], w, 'REPLACE')
2977         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
2978         return {'FINISHED'}
2979
2980 class TISSUE_PT_color(Panel):
2981     bl_label = "Tissue Tools"
2982     bl_category = "Tissue"
2983     bl_space_type = "VIEW_3D"
2984     bl_region_type = "UI"
2985     #bl_options = {'DEFAULT_CLOSED'}
2986     bl_context = "vertexpaint"
2987
2988     def draw(self, context):
2989         layout = self.layout
2990         col = layout.column(align=True)
2991         col.operator("object.vertex_colors_to_vertex_groups",
2992             icon="GROUP_VERTEX", text="Convert to Weight")
2993
2994 class TISSUE_PT_weight(Panel):
2995     bl_label = "Tissue Tools"
2996     bl_category = "Tissue"
2997     bl_space_type = "VIEW_3D"
2998     bl_region_type = "UI"
2999     #bl_options = {'DEFAULT_CLOSED'}
3000     bl_context = "weightpaint"
3001
3002     def draw(self, context):
3003         layout = self.layout
3004         col = layout.column(align=True)
3005         #if context.object.type == 'MESH' and context.mode == 'OBJECT':
3006             #col.label(text="Transform:")
3007             #col.separator()
3008         #elif bpy.context.mode == 'PAINT_WEIGHT':
3009         col.label(text="Weight Generate:")
3010         #col.operator(
3011         #    "object.vertex_colors_to_vertex_groups", icon="GROUP_VCOL")
3012         col.operator("object.face_area_to_vertex_groups", icon="FACESEL")
3013         col.operator("object.curvature_to_vertex_groups", icon="SMOOTHCURVE")
3014         col.operator("object.tissue_weight_distance", icon="TRACKING")
3015         row = col.row(align=True)
3016         try: row.operator("object.weight_formula", icon="CON_TRANSFORM")
3017         except: row.operator("object.weight_formula")#, icon="CON_TRANSFORM")
3018         row.operator("object.update_weight_formula", icon="FILE_REFRESH", text='')#, icon="CON_TRANSFORM")
3019         #col.label(text="Weight Processing:")
3020         col.separator()
3021
3022         # TO BE FIXED
3023         col.operator("object.weight_laplacian", icon="SMOOTHCURVE")
3024
3025         col.label(text="Weight Edit:")
3026         col.operator("object.harmonic_weight", icon="IPO_ELASTIC")
3027         col.operator("object.random_weight", icon="RNDCURVE")
3028         col.separator()
3029         col.label(text="Deformation Analysis:")
3030         col.operator("object.edges_deformation", icon="DRIVER_DISTANCE")#FULLSCREEN_ENTER")
3031         col.operator("object.edges_bending", icon="DRIVER_ROTATIONAL_DIFFERENCE")#"MOD_SIMPLEDEFORM")
3032         col.separator()
3033         col.label(text="Weight Curves:")
3034         #col.operator("object.weight_contour_curves", icon="MOD_CURVE")
3035         col.operator("object.tissue_weight_streamlines", icon="ANIM")
3036         col.operator("object.tissue_weight_contour_curves_pattern", icon="FORCE_TURBULENCE")
3037         col.separator()
3038         col.operator("object.weight_contour_displace", icon="MOD_DISPLACE")
3039         col.operator("object.weight_contour_mask", icon="MOD_MASK")
3040         col.separator()
3041         col.label(text="Simulations:")
3042         #col.operator("object.reaction_diffusion", icon="MOD_OCEAN")
3043         col.operator("object.start_reaction_diffusion",
3044                     icon="EXPERIMENTAL",
3045                     text="Reaction-Diffusion")
3046         col.separator()
3047         col.label(text="Materials:")
3048         col.operator("object.random_materials", icon='COLOR')
3049         col.operator("object.weight_to_materials", icon='GROUP_VERTEX')
3050         col.separator()
3051         col.label(text="Weight Convert:")
3052         col.operator("object.vertex_group_to_vertex_colors", icon="GROUP_VCOL",
3053             text="Convert to Colors")
3054         col.operator("object.vertex_group_to_uv", icon="UV",
3055             text="Convert to UV")
3056
3057         #col.prop(context.object, "reaction_diffusion_run", icon="PLAY", text="Run Simulation")
3058         ####col.prop(context.object, "reaction_diffusion_run")
3059         #col.separator()
3060         #col.label(text="Vertex Color from:")
3061         #col.operator("object.vertex_group_to_vertex_colors", icon="GROUP_VERTEX")
3062
3063
3064
3065
3066 class start_reaction_diffusion(Operator):
3067     bl_idname = "object.start_reaction_diffusion"
3068     bl_label = "Start Reaction Diffusion"
3069     bl_description = ("Run a Reaction-Diffusion based on existing Vertex Groups: A and B")
3070     bl_options = {'REGISTER', 'UNDO'}
3071
3072     run : BoolProperty(
3073         name="Run Reaction-Diffusion", default=True, description="Compute a new iteration on frame changes")
3074
3075     time_steps : IntProperty(
3076         name="Steps", default=10, min=0, soft_max=50,
3077         description="Number of Steps")
3078
3079     dt : FloatProperty(
3080         name="dt", default=1, min=0, soft_max=0.2,
3081         description="Time Step")
3082
3083     diff_a : FloatProperty(
3084         name="Diff A", default=0.18, min=0, soft_max=2,
3085         description="Diffusion A")
3086
3087     diff_b : FloatProperty(
3088         name="Diff B", default=0.09, min=0, soft_max=2,
3089         description="Diffusion B")
3090
3091     f : FloatProperty(
3092         name="f", default=0.055, min=0, soft_min=0.01, soft_max=0.06, max=0.1, precision=4,
3093         description="Feed Rate")
3094
3095     k : FloatProperty(
3096         name="k", default=0.062, min=0, soft_min=0.035, soft_max=0.065, max=0.1, precision=4,
3097         description="Kill Rate")
3098
3099     @classmethod
3100     def poll(cls, context):
3101         return context.object.type == 'MESH' and context.mode != 'EDIT_MESH'
3102
3103     def execute(self, context):
3104         reaction_diffusion_add_handler(self, context)
3105         set_animatable_fix_handler(self, context)
3106
3107         ob = context.object
3108
3109         ob.reaction_diffusion_settings.run = self.run
3110         ob.reaction_diffusion_settings.dt = self.dt
3111         ob.reaction_diffusion_settings.time_steps = self.time_steps
3112         ob.reaction_diffusion_settings.f = self.f
3113         ob.reaction_diffusion_settings.k = self.k
3114         ob.reaction_diffusion_settings.diff_a = self.diff_a
3115         ob.reaction_diffusion_settings.diff_b = self.diff_b
3116
3117
3118         # check vertex group A
3119         try:
3120             vg = ob.vertex_groups['A']
3121         except:
3122             ob.vertex_groups.new(name='A')
3123         # check vertex group B
3124         try:
3125             vg = ob.vertex_groups['B']
3126         except:
3127             ob.vertex_groups.new(name='B')
3128
3129         for v in ob.data.vertices:
3130             ob.vertex_groups['A'].add([v.index], 1, 'REPLACE')
3131             ob.vertex_groups['B'].add([v.index], 0, 'REPLACE')
3132
3133         ob.vertex_groups.update()
3134         ob.data.update()
3135         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
3136
3137         return {'FINISHED'}
3138
3139 class reset_reaction_diffusion_weight(Operator):
3140     bl_idname = "object.reset_reaction_diffusion_weight"
3141     bl_label = "Reset Reaction Diffusion Weight"
3142     bl_description = ("Set A and B weight to default values")
3143     bl_options = {'REGISTER', 'UNDO'}
3144
3145     @classmethod
3146     def poll(cls, context):
3147         return context.object.type == 'MESH' and context.mode != 'EDIT_MESH'
3148
3149     def execute(self, context):
3150         reaction_diffusion_add_handler(self, context)
3151         set_animatable_fix_handler(self, context)
3152
3153         ob = context.object
3154
3155         # check vertex group A
3156         try:
3157             vg = ob.vertex_groups['A']
3158         except:
3159             ob.vertex_groups.new(name='A')
3160         # check vertex group B
3161         try:
3162             vg = ob.vertex_groups['B']
3163         except:
3164             ob.vertex_groups.new(name='B')
3165
3166         for v in ob.data.vertices:
3167             ob.vertex_groups['A'].add([v.index], 1, 'REPLACE')
3168             ob.vertex_groups['B'].add([v.index], 0, 'REPLACE')
3169
3170         ob.vertex_groups.update()
3171         ob.data.update()
3172         bpy.ops.object.mode_set(mode='WEIGHT_PAINT')
3173
3174         return {'FINISHED'}
3175
3176 class bake_reaction_diffusion(Operator):
3177     bl_idname = "object.bake_reaction_diffusion"
3178     bl_label = "Bake Data"
3179     bl_description = ("Bake the Reaction-Diffusion to the cache directory")
3180     bl_options = {'REGISTER', 'UNDO'}
3181
3182     @classmethod
3183     def poll(cls, context):
3184         return context.object.type == 'MESH' and context.mode != 'EDIT_MESH'
3185
3186     def execute(self, context):
3187         ob = context.object
3188         props = ob.reaction_diffusion_settings
3189         if props.fast_bake:
3190             bool_run = props.run
3191             props.run = False
3192             context.scene.frame_current = props.cache_frame_start
3193             fast_bake_def(ob, frame_start=props.cache_frame_start, frame_end=props.cache_frame_end)
3194             #create_fast_bake_def(ob, frame_start=props.cache_frame_start, frame_end=props.cache_frame_end)
3195             context.scene.frame_current = props.cache_frame_end
3196             props.run = bool_run
3197         else:
3198             for i in range(props.cache_frame_start, props.cache_frame_end):
3199                 context.scene.frame_current = i
3200                 reaction_diffusion_def(ob, bake=True)
3201         props.bool_cache = True
3202
3203         return {'FINISHED'}
3204
3205 class reaction_diffusion_free_data(Operator):
3206     bl_idname = "object.reaction_diffusion_free_data"
3207     bl_label = "Free Data"
3208     bl_description = ("Free Reaction-Diffusion data")
3209     bl_options = {'REGISTER', 'UNDO'}
3210
3211     @classmethod
3212     def poll(cls, context):
3213         return context.object.type == 'MESH'
3214
3215     def execute(self, context):
3216         ob = context.object
3217         props = ob.reaction_diffusion_settings
3218         props.bool_cache = False
3219
3220         folder = Path(props.cache_dir)
3221         for i in range(props.cache_frame_start, props.cache_frame_end):
3222             data_a = folder / "a_{:04d}".format(i)
3223             if os.path.exists(data_a):
3224                 os.remove(data_a)
3225             data_a = folder / "b_{:04d}".format(i)
3226             if os.path.exists(data_a):
3227                 os.remove(data_a)
3228         return {'FINISHED'}
3229
3230 from bpy.app.handlers import persistent
3231
3232 def reaction_diffusion_scene(scene, bake=False):
3233     for ob in scene.objects:
3234         if ob.reaction_diffusion_settings.run:
3235             reaction_diffusion_def(ob)
3236
3237 def reaction_diffusion_def(ob, bake=False):
3238
3239     scene = bpy.context.scene
3240     start = time.time()
3241     if type(ob) == bpy.types.Scene: return None
3242     props = ob.reaction_diffusion_settings
3243
3244     if bake or props.bool_cache:
3245         if props.cache_dir == '':
3246             letters = string.ascii_letters
3247             random_name = ''.join(rnd.choice(letters) for i in range(6))
3248             if bpy.context.blend_data.filepath == '':
3249                 folder = Path(bpy.context.preferences.filepaths.temporary_directory)
3250                 folder = folder / 'reaction_diffusion_cache' / random_name
3251             else:
3252                 folder = '//' + Path(bpy.context.blend_data.filepath).stem
3253                 folder = Path(bpy.path.abspath(folder)) / 'reaction_diffusion_cache' / random_name
3254             folder.mkdir(parents=True, exist_ok=True)
3255             props.cache_dir = str(folder)
3256         else:
3257             folder = Path(props.cache_dir)
3258
3259     me = ob.data
3260     n_edges = len(me.edges)
3261     n_verts = len(me.vertices)
3262     a = np.zeros(n_verts)
3263     b = np.zeros(n_verts)
3264
3265     print("{:6d} Reaction-Diffusion: {}".format(scene.frame_current, ob.name))
3266
3267     if not props.bool_cache:
3268
3269         if props.bool_mod:
3270             # hide deforming modifiers
3271             mod_visibility = []
3272             for m in ob.modifiers:
3273                 mod_visibility.append(m.show_viewport)
3274                 if not mod_preserve_shape(m): m.show_viewport = False
3275
3276             # evaluated mesh
3277             dg = bpy.context.evaluated_depsgraph_get()
3278             ob_eval = ob.evaluated_get(dg)
3279             me = bpy.data.meshes.new_from_object(ob_eval, preserve_all_data_layers=True, depsgraph=dg)
3280
3281             # set original visibility
3282             for v, m in zip(mod_visibility, ob.modifiers):
3283                 m.show_viewport = v
3284             ob.modifiers.update()
3285
3286         bm = bmesh.new()   # create an empty BMesh
3287         bm.from_mesh(me)   # fill it in from a Mesh
3288         dvert_lay = bm.verts.layers.deform.active
3289
3290         dt = props.dt
3291         time_steps = props.time_steps
3292         f = props.f
3293         k = props.k
3294         diff_a = props.diff_a
3295         diff_b = props.diff_b
3296         scale = props.diff_mult
3297
3298         brush_mult = props.brush_mult
3299
3300         # store weight values
3301         if 'dB' in ob.vertex_groups: db = np.zeros(n_verts)
3302         if 'grad' in ob.vertex_groups: grad = np.zeros(n_verts)
3303
3304         if props.vertex_group_diff_a != '': diff_a = np.zeros(n_verts)
3305         if props.vertex_group_diff_b != '': diff_b = np.zeros(n_verts)
3306         if props.vertex_group_scale != '': scale = np.zeros(n_verts)
3307         if props.vertex_group_f != '': f = np.zeros(n_verts)
3308         if props.vertex_group_k != '': k = np.zeros(n_verts)
3309         if props.vertex_group_brush != '': brush = np.zeros(n_verts)
3310         else: brush = 0
3311
3312         group_index_a = ob.vertex_groups["A"].index
3313         group_index_b = ob.vertex_groups["B"].index
3314         a = bmesh_get_weight_numpy(group_index_a, dvert_lay, bm.verts)
3315         b = bmesh_get_weight_numpy(group_index_b, dvert_lay, bm.verts)
3316
3317         if props.vertex_group_diff_a != '':
3318             group_index = ob.vertex_groups[props.vertex_group_diff_a].index
3319             diff_a = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3320             if props.invert_vertex_group_diff_a:
3321                 vg_bounds = (props.min_diff_a, props.max_diff_a)
3322             else:
3323                 vg_bounds = (props.max_diff_a, props.min_diff_a)
3324             diff_a = np.interp(diff_a, (0,1), vg_bounds)
3325
3326         if props.vertex_group_diff_b != '':
3327             group_index = ob.vertex_groups[props.vertex_group_diff_b].index
3328             diff_b = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3329             if props.invert_vertex_group_diff_b:
3330                 vg_bounds = (props.max_diff_b, props.min_diff_b)
3331             else:
3332                 vg_bounds = (props.min_diff_b, props.max_diff_b)
3333             diff_b = np.interp(diff_b, (0,1), vg_bounds)
3334
3335         if props.vertex_group_scale != '':
3336             group_index = ob.vertex_groups[props.vertex_group_scale].index
3337             scale = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3338             if props.invert_vertex_group_scale:
3339                 vg_bounds = (props.max_scale, props.min_scale)
3340             else:
3341                 vg_bounds = (props.min_scale, props.max_scale)
3342             scale = np.interp(scale, (0,1), vg_bounds)
3343
3344         if props.vertex_group_f != '':
3345             group_index = ob.vertex_groups[props.vertex_group_f].index
3346             f = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3347             if props.invert_vertex_group_f:
3348                 vg_bounds = (props.max_f, props.min_f)
3349             else:
3350                 vg_bounds = (props.min_f, props.max_f)
3351             f = np.interp(f, (0,1), vg_bounds, )
3352
3353         if props.vertex_group_k != '':
3354             group_index = ob.vertex_groups[props.vertex_group_k].index
3355             k = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3356             if props.invert_vertex_group_k:
3357                 vg_bounds = (props.max_k, props.min_k)
3358             else:
3359                 vg_bounds = (props.min_k, props.max_k)
3360             k = np.interp(k, (0,1), vg_bounds)
3361
3362         if props.vertex_group_brush != '':
3363             group_index = ob.vertex_groups[props.vertex_group_brush].index
3364             brush = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3365             brush *= brush_mult
3366
3367
3368         #timeElapsed = time.time() - start
3369         #print('RD - Read Vertex Groups:',timeElapsed)
3370         #start = time.time()
3371
3372         diff_a *= scale
3373         diff_b *= scale
3374
3375         edge_verts = [0]*n_edges*2
3376         me.edges.foreach_get("vertices", edge_verts)
3377         edge_verts = np.array(edge_verts)
3378
3379         if 'gradient' in ob.vertex_groups.keys() and False:
3380             group_index = ob.vertex_groups['gradient'].index
3381             gradient = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3382
3383             arr = (np.arange(n_edges)*2).astype(int)
3384             id0 = edge_verts[arr]
3385             id1 = edge_verts[arr+1]
3386
3387             #gradient = np.abs(gradient[id0] - gradient[id1])
3388             gradient = gradient[id1] - gradient[id0]
3389             gradient /= np.max(gradient)
3390             sign = np.sign(gradient)
3391             sign[sign==0] = 1
3392             gradient = (0.05*abs(gradient) + 0.95)*sign
3393             #gradient *= (1-abs(gradient)
3394             #gradient = 0.2*(1-gradient) + 0.95
3395
3396         #gradient = get_uv_edge_vectors(me)
3397         #uv_dir = Vector((0.5,0.5,0)).normalized()
3398         #gradient = np.array([abs(g.dot(uv_dir.normalized())) for g in gradient])
3399         #gradient = (gradient + 0.5)/2
3400         #gradient = np.array([max(0,g.dot(uv_dir.normalized())) for g in gradient])
3401
3402         timeElapsed = time.time() - start
3403         print('       Preparation Time:',timeElapsed)
3404         start = time.time()
3405
3406         try:
3407             _f = f if type(f) is np.ndarray else np.array((f,))
3408             _k = k if type(k) is np.ndarray else np.array((k,))
3409             _diff_a = diff_a if type(diff_a) is np.ndarray else np.array((diff_a,))
3410             _diff_b = diff_b if type(diff_b) is np.ndarray else np.array((diff_b,))
3411             _brush = brush if type(brush) is np.ndarray else np.array((brush,))
3412
3413             #a, b = numba_reaction_diffusion_anisotropic(n_verts, n_edges, edge_verts, a, b, _brush, _diff_a, _diff_b, _f, _k, dt, time_steps, gradient)
3414             a, b = numba_reaction_diffusion(n_verts, n_edges, edge_verts, a, b, _brush, _diff_a, _diff_b, _f, _k, dt, time_steps)
3415         except:
3416             print('Not using Numba! The simulation could be slow.')
3417             arr = np.arange(n_edges)*2
3418             id0 = edge_verts[arr]     # first vertex indices for each edge
3419             id1 = edge_verts[arr+1]   # second vertex indices for each edge
3420             for i in range(time_steps):
3421                 b += brush
3422                 lap_a = np.zeros(n_verts)
3423                 lap_b = np.zeros(n_verts)
3424                 lap_a0 =  a[id1] -  a[id0]   # laplacian increment for first vertex of each edge
3425                 lap_b0 =  b[id1] -  b[id0]   # laplacian increment for first vertex of each edge
3426
3427                 np.add.at(lap_a, id0, lap_a0)
3428                 np.add.at(lap_b, id0, lap_b0)
3429                 np.add.at(lap_a, id1, -lap_a0)
3430                 np.add.at(lap_b, id1, -lap_b0)
3431
3432                 ab2 = a*b**2
3433                 a += eval("(diff_a*lap_a - ab2 + f*(1-a))*dt")
3434                 b += eval("(diff_b*lap_b + ab2 - (k+f)*b)*dt")
3435                 #a += (diff_a*lap_a - ab2 + f*(1-a))*dt
3436                 #b += (diff_b*lap_b + ab2 - (k+f)*b)*dt
3437
3438                 a = nan_to_num(a)
3439                 b = nan_to_num(b)
3440
3441         timeElapsed = time.time() - start
3442         print('       Simulation Time:',timeElapsed)
3443
3444     if bake:
3445         if not(os.path.exists(folder)):
3446             os.mkdir(folder)
3447         file_name = folder / "a_{:04d}".format(scene.frame_current)
3448         a.tofile(file_name)
3449         file_name = folder / "b_{:04d}".format(scene.frame_current)
3450         b.tofile(file_name)
3451     elif props.bool_cache:
3452         try:
3453             file_name = folder / "a_{:04d}".format(scene.frame_current)
3454             a = np.fromfile(file_name)
3455             file_name = folder / "b_{:04d}".format(scene.frame_current)
3456             b = np.fromfile(file_name)
3457         except:
3458             print('       Cannot read cache.')
3459             return
3460
3461     if props.update_weight_a or props.update_weight_b:
3462         start = time.time()
3463         if props.update_weight_a:
3464             if 'A' in ob.vertex_groups.keys():
3465                 vg_a = ob.vertex_groups['A']
3466             else:
3467                 vg_a = ob.vertex_groups.new(name='A')
3468         else:
3469             vg_a = None
3470         if props.update_weight_b:
3471             if 'B' in ob.vertex_groups.keys():
3472                 vg_b = ob.vertex_groups['B']
3473             else:
3474                 vg_b = ob.vertex_groups.new(name='B')
3475         else:
3476             vg_b = None
3477         if vg_a == vg_b == None:
3478             pass
3479         else:
3480             if ob.mode == 'WEIGHT_PAINT':# or props.bool_cache:
3481                 # slower, but prevent crashes
3482                 for i in range(n_verts):
3483                     if vg_a: vg_a.add([i], a[i], 'REPLACE')
3484                     if vg_b: vg_b.add([i], b[i], 'REPLACE')
3485             else:
3486                 if props.bool_mod or props.bool_cache:
3487                     #bm.free()               # release old bmesh
3488                     bm = bmesh.new()        # create an empty BMesh
3489                     bm.from_mesh(ob.data)   # fill it in from a Mesh
3490                     dvert_lay = bm.verts.layers.deform.active
3491                 # faster, but can cause crashes while painting weight
3492                 if vg_a: index_a = vg_a.index
3493                 if vg_b: index_b = vg_b.index
3494                 for i, v in enumerate(bm.verts):
3495                     dvert = v[dvert_lay]
3496                     if vg_a: dvert[index_a] = a[i]
3497                     if vg_b: dvert[index_b] = b[i]
3498                 bm.to_mesh(ob.data)
3499                 bm.free()
3500         print('       Writing Vertex Groups Time:',time.time() - start)
3501     if props.normalize:
3502         min_a = np.min(a)
3503         max_a = np.max(a)
3504         min_b = np.min(b)
3505         max_b = np.max(b)
3506         a = (a - min_a)/(max_a - min_a)
3507         b = (b - min_b)/(max_b - min_b)
3508     split_a = None
3509     split_b = None
3510     splitted = False
3511     if props.update_colors:#_a or props.update_colors_b:
3512         start = time.time()
3513         loops_size = get_attribute_numpy(me.polygons, attribute='loop_total', mult=1)
3514         n_colors = np.sum(loops_size)
3515         v_id = np.ones(n_colors)
3516         me.polygons.foreach_get('vertices',v_id)
3517         v_id = v_id.astype(int)
3518         #v_id = np.array([v for p in ob.data.polygons for v in p.vertices])
3519         '''
3520         if props.update_colors_b:
3521             if 'B' in ob.data.vertex_colors.keys():
3522                 vc = ob.data.vertex_colors['B']
3523             else:
3524                 vc = ob.data.vertex_colors.new(name='B')
3525             c_val = b[v_id]
3526             c_val = np.repeat(c_val, 4, axis=0)
3527             vc.data.foreach_set('color',c_val)
3528
3529         if props.update_colors_a:
3530             if 'A' in ob.data.vertex_colors.keys():
3531                 vc = ob.data.vertex_colors['A']
3532             else:
3533                 vc = ob.data.vertex_colors.new(name='A')
3534             c_val = a[v_id]
3535             c_val = np.repeat(c_val, 4, axis=0)
3536             vc.data.foreach_set('color',c_val)
3537         '''
3538         split_a = a[v_id,None]
3539         split_b = b[v_id,None]
3540         splitted = True
3541         ones = np.ones((n_colors,1))
3542         #rgba = np.concatenate((split_a,split_b,-split_b+split_a,ones),axis=1).flatten()
3543         rgba = np.concatenate((split_a,split_b,ones,ones),axis=1).flatten()
3544         if 'AB' in ob.data.vertex_colors.keys():
3545             vc = ob.data.vertex_colors['AB']
3546         else:
3547             vc = ob.data.vertex_colors.new(name='AB')
3548         vc.data.foreach_set('color',rgba)
3549         ob.data.vertex_colors.update()
3550
3551         print('       Writing Vertex Colors Time:',time.time() - start)
3552     if props.update_uv:
3553         start = time.time()
3554         if 'AB' in me.uv_layers.keys():
3555             uv_layer = me.uv_layers['AB']
3556         else:
3557             uv_layer = me.uv_layers.new(name='AB')
3558         if not splitted:
3559             loops_size = get_attribute_numpy(me.polygons, attribute='loop_total', mult=1)
3560             n_data = np.sum(loops_size)
3561             v_id = np.ones(n_data)
3562             me.polygons.foreach_get('vertices',v_id)
3563             v_id = v_id.astype(int)
3564             split_a = a[v_id,None]
3565             split_b = b[v_id,None]
3566         uv = np.concatenate((split_a,split_b),axis=1).flatten()
3567         uv_layer.data.foreach_set('uv',uv)
3568         me.uv_layers.update()
3569         print('       Writing UV Map Time:',time.time() - start)
3570
3571     for ps in ob.particle_systems:
3572         if ps.vertex_group_density == 'B' or ps.vertex_group_density == 'A':
3573             ps.invert_vertex_group_density = not ps.invert_vertex_group_density
3574             ps.invert_vertex_group_density = not ps.invert_vertex_group_density
3575
3576     if props.bool_mod and not props.bool_cache: bpy.data.meshes.remove(me)
3577
3578 def fast_bake_def(ob, frame_start=1, frame_end=250):
3579     scene = bpy.context.scene
3580     start = time.time()
3581     if type(ob) == bpy.types.Scene: return None
3582     props = ob.reaction_diffusion_settings
3583
3584     # Define cache folder
3585     if props.cache_dir == '':
3586         letters = string.ascii_letters
3587         random_name = ''.join(rnd.choice(letters) for i in range(6))
3588         if bpy.context.blend_data.filepath == '':
3589             folder = Path(bpy.context.preferences.filepaths.temporary_directory)
3590             folder = folder / 'reaction_diffusion_cache' / random_name
3591         else:
3592             folder = '//' + Path(bpy.context.blend_data.filepath).stem
3593             folder = Path(bpy.path.abspath(folder)) / 'reaction_diffusion_cache' / random_name
3594         folder.mkdir(parents=True, exist_ok=True)
3595         props.cache_dir = str(folder)
3596     else:
3597         folder = Path(props.cache_dir)
3598
3599     if props.bool_mod:
3600         # hide deforming modifiers
3601         mod_visibility = []
3602         for m in ob.modifiers:
3603             mod_visibility.append(m.show_viewport)
3604             if not mod_preserve_shape(m): m.show_viewport = False
3605
3606         # evaluated mesh
3607         dg = bpy.context.evaluated_depsgraph_get()
3608         ob_eval = ob.evaluated_get(dg)
3609         me = bpy.data.meshes.new_from_object(ob_eval, preserve_all_data_layers=True, depsgraph=dg)
3610
3611         # set original visibility
3612         for v, m in zip(mod_visibility, ob.modifiers):
3613             m.show_viewport = v
3614         ob.modifiers.update()
3615     else:
3616         me = ob.data
3617
3618     bm = bmesh.new()   # create an empty BMesh
3619     bm.from_mesh(me)   # fill it in from a Mesh
3620     dvert_lay = bm.verts.layers.deform.active
3621     n_edges = len(me.edges)
3622     n_verts = len(me.vertices)
3623     a = np.zeros(n_verts)
3624     b = np.zeros(n_verts)
3625     group_index_a = ob.vertex_groups["A"].index
3626     group_index_b = ob.vertex_groups["B"].index
3627
3628     dt = props.dt
3629     time_steps = props.time_steps
3630     f = props.f
3631     k = props.k
3632     diff_a = props.diff_a
3633     diff_b = props.diff_b
3634     scale = props.diff_mult
3635
3636     brush_mult = props.brush_mult
3637
3638     # store weight values
3639     if 'dB' in ob.vertex_groups: db = np.zeros(n_verts)
3640     if 'grad' in ob.vertex_groups: grad = np.zeros(n_verts)
3641
3642     if props.vertex_group_diff_a != '': diff_a = np.zeros(n_verts)
3643     if props.vertex_group_diff_b != '': diff_b = np.zeros(n_verts)
3644     if props.vertex_group_scale != '': scale = np.zeros(n_verts)
3645     if props.vertex_group_f != '': f = np.zeros(n_verts)
3646     if props.vertex_group_k != '': k = np.zeros(n_verts)
3647     if props.vertex_group_brush != '': brush = np.zeros(n_verts)
3648     else: brush = 0
3649
3650     a = bmesh_get_weight_numpy(group_index_a, dvert_lay, bm.verts)
3651     b = bmesh_get_weight_numpy(group_index_b, dvert_lay, bm.verts)
3652
3653     if props.vertex_group_diff_a != '':
3654         group_index = ob.vertex_groups[props.vertex_group_diff_a].index
3655         diff_a = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3656         if props.invert_vertex_group_diff_a:
3657             vg_bounds = (props.min_diff_a, props.max_diff_a)
3658         else:
3659             vg_bounds = (props.max_diff_a, props.min_diff_a)
3660         diff_a = np.interp(diff_a, (0,1), vg_bounds)
3661
3662     if props.vertex_group_diff_b != '':
3663         group_index = ob.vertex_groups[props.vertex_group_diff_b].index
3664         diff_b = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3665         if props.invert_vertex_group_diff_b:
3666             vg_bounds = (props.max_diff_b, props.min_diff_b)
3667         else:
3668             vg_bounds = (props.min_diff_b, props.max_diff_b)
3669         diff_b = np.interp(diff_b, (0,1), vg_bounds)
3670
3671     if props.vertex_group_scale != '':
3672         group_index = ob.vertex_groups[props.vertex_group_scale].index
3673         scale = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3674         if props.invert_vertex_group_scale:
3675             vg_bounds = (props.max_scale, props.min_scale)
3676         else:
3677             vg_bounds = (props.min_scale, props.max_scale)
3678         scale = np.interp(scale, (0,1), vg_bounds)
3679
3680     if props.vertex_group_f != '':
3681         group_index = ob.vertex_groups[props.vertex_group_f].index
3682         f = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3683         if props.invert_vertex_group_f:
3684             vg_bounds = (props.max_f, props.min_f)
3685         else:
3686             vg_bounds = (props.min_f, props.max_f)
3687         f = np.interp(f, (0,1), vg_bounds, )
3688
3689     if props.vertex_group_k != '':
3690         group_index = ob.vertex_groups[props.vertex_group_k].index
3691         k = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3692         if props.invert_vertex_group_k:
3693             vg_bounds = (props.max_k, props.min_k)
3694         else:
3695             vg_bounds = (props.min_k, props.max_k)
3696         k = np.interp(k, (0,1), vg_bounds)
3697
3698     if props.vertex_group_brush != '':
3699         group_index = ob.vertex_groups[props.vertex_group_brush].index
3700         brush = bmesh_get_weight_numpy(group_index, dvert_lay, bm.verts)
3701         brush *= brush_mult
3702
3703     diff_a *= scale
3704     diff_b *= scale
3705
3706     edge_verts = [0]*n_edges*2
3707     me.edges.foreach_get("vertices", edge_verts)
3708
3709     gradient = get_uv_edge_vectors(me)
3710     uv_dir = Vector((0.5,0.5,0))
3711     #gradient = [abs(g.dot(uv_dir)) for g in gradient]
3712     gradient = [max(0,g.dot(uv_dir)) for g in gradient]
3713
3714     timeElapsed = time.time() - start
3715     print('       Preparation Time:',timeElapsed)
3716     start = time.time()
3717
3718     try:
3719         edge_verts = np.array(edge_verts)
3720         _f = f if type(f) is np.ndarray else np.array((f,))
3721         _k = k if type(k) is np.ndarray else np.array((k,))
3722         _diff_a = diff_a if type(diff_a) is np.ndarray else np.array((diff_a,))
3723         _diff_b = diff_b if type(diff_b) is np.ndarray else np.array((diff_b,))
3724         _brush = brush if type(brush) is np.ndarray else np.array((brush,))
3725
3726         run_rd = False
3727         for j in range(props.cache_frame_start, props.cache_frame_end+1):
3728             start2 = time.time()
3729             print("{:6d} Reaction-Diffusion: {}".format(j, ob.name))
3730             if run_rd:
3731                 b += _brush
3732                 a, b = numba_reaction_diffusion(n_verts, n_edges, edge_verts, a, b, _brush, _diff_a, _diff_b, _f, _k, dt, time_steps)
3733             else:
3734                 run_rd = True
3735
3736             if not(os.path.exists(folder)):
3737                 os.mkdir(folder)
3738             file_name = folder / "a_{:04d}".format(j)
3739             a.tofile(file_name)
3740             file_name = folder / "b_{:04d}".format(j)
3741             b.tofile(file_name)
3742
3743             timeElapsed = time.time() - start2
3744             print('       Simulation Time:',timeElapsed)
3745
3746     except:
3747         print('Not using Numba! The simulation could be slow.')
3748         edge_verts = np.array(edge_verts)
3749         arr = np.arange(n_edges)*2
3750         id0 = edge_verts[arr]     # first vertex indices for each edge
3751         id1 = edge_verts[arr+1]   # second vertex indices for each edge
3752         for j in range(props.cache_frame_start, props.cache_frame_end):
3753             for i in range(time_steps):
3754                 b += brush
3755                 lap_a = np.zeros(n_verts)
3756                 lap_b = np.zeros(n_verts)
3757                 lap_a0 =  a[id1] -  a[id0]   # laplacian increment for first vertex of each edge
3758                 lap_b0 =  b[id1] -  b[id0]   # laplacian increment for first vertex of each edge
3759
3760                 np.add.at(lap_a, id0, lap_a0)
3761                 np.add.at(lap_b, id0, lap_b0)
3762                 np.add.at(lap_a, id1, -lap_a0)
3763                 np.add.at(lap_b, id1, -lap_b0)
3764
3765                 ab2 = a*b**2
3766                 a += eval("(diff_a*lap_a - ab2 + f*(1-a))*dt")
3767                 b += eval("(diff_b*lap_b + ab2 - (k+f)*b)*dt")
3768
3769                 a = nan_to_num(a)
3770                 b = nan_to_num(b)
3771
3772             if not(os.path.exists(folder)):
3773                 os.mkdir(folder)
3774             file_name = folder / "a_{:04d}".format(j)
3775             a.tofile(file_name)
3776             file_name = folder / "b_{:04d}".format(j)
3777             b.tofile(file_name)
3778
3779     if ob.mode == 'WEIGHT_PAINT':
3780         # slower, but prevent crashes
3781         vg_a = ob.vertex_groups['A']
3782         vg_b = ob.vertex_groups['B']
3783         for i in range(n_verts):
3784             vg_a.add([i], a[i], 'REPLACE')
3785             vg_b.add([i], b[i], 'REPLACE')
3786     else:
3787         if props.bool_mod:
3788             bm.free()               # release old bmesh
3789             bm = bmesh.new()        # create an empty BMesh
3790             bm.from_mesh(ob.data)   # fill it in from a Mesh
3791             dvert_lay = bm.verts.layers.deform.active
3792         # faster, but can cause crashes while painting weight
3793         for i, v in enumerate(bm.verts):
3794             dvert = v[dvert_lay]
3795             dvert[group_index_a] = a[i]
3796             dvert[group_index_b] = b[i]
3797         bm.to_mesh(ob.data)
3798
3799     # Update Vertex Colors
3800     if 'A' in ob.data.vertex_colors or 'B' in ob.data.vertex_colors:
3801         v_id = np.array([v for p in ob.data.polygons for v in p.vertices])
3802
3803         if 'B' in ob.data.vertex_colors:
3804             c_val = b[v_id]
3805             c_val = np.repeat(c_val, 4, axis=0)
3806             vc = ob.data.vertex_colors['B']
3807             vc.data.foreach_set('color',c_val.tolist())
3808
3809         if 'A' in ob.data.vertex_colors:
3810             c_val = a[v_id]
3811             c_val = np.repeat(c_val, 4, axis=0)
3812             vc = ob.data.vertex_colors['A']
3813             vc.data.foreach_set('color',c_val.tolist())
3814
3815     for ps in ob.particle_systems:
3816         if ps.vertex_group_density == 'B' or ps.vertex_group_density == 'A':
3817             ps.invert_vertex_group_density = not ps.invert_vertex_group_density
3818             ps.invert_vertex_group_density = not ps.invert_vertex_group_density
3819
3820     if props.bool_mod: bpy.data.meshes.remove(me)
3821     bm.free()
3822     timeElapsed = time.time() - start
3823     print('       Closing Time:',timeElapsed)
3824
3825 def create_fast_bake_def(ob, frame_start=1, frame_end=250):
3826     scene = bpy.context.scene
3827     start = time.time()
3828     if type(ob) == bpy.types.Scene: return None
3829     props = ob.reaction_diffusion_settings
3830
3831     dt = props.dt
3832     time_steps = props.time_steps
3833     scale = props.diff_mult
3834
3835     if props.cache_dir == '':
3836         letters = string.ascii_letters
3837         random_name = ''.join(rnd.choice(letters) for i in range(6))
3838         if bpy.context.blend_data.filepath == '':
3839             folder = Path(bpy.context.preferences.filepaths.temporary_directory)
3840             folder = folder / 'reaction_diffusion_cache' / random_name
3841         else:
3842             folder = '//' + Path(bpy.context.blend_data.filepath).stem
3843             folder = Path(bpy.path.abspath(folder)) / 'reaction_diffusion_cache' / random_name
3844         folder.mkdir(parents=True, exist_ok=True)
3845         props.cache_dir = str(folder)
3846     else:
3847         folder = Path(props.cache_dir)
3848
3849     if props.bool_mod:
3850         # hide deforming modifiers
3851         mod_visibility = []
3852         for m in ob.modifiers:
3853             mod_visibility.append(m.show_viewport)
3854             if not mod_preserve_shape(m): m.show_viewport = False
3855
3856         # evaluated mesh
3857         dg = bpy.context.evaluated_depsgraph_get()
3858         ob_eval = ob.evaluated_get(dg)
3859         me = bpy.data.meshes.new_from_object(ob_eval, preserve_all_data_layers=True, depsgraph=dg)
3860
3861         # set original visibility
3862         for v, m in zip(mod_visibility, ob.modifiers):
3863             m.show_viewport = v
3864         ob.modifiers.update()
3865     else:
3866         me = ob.data
3867
3868     bm = bmesh.new()   # create an empty BMesh
3869     bm.from_mesh(me)   # fill it in from a Mesh
3870     verts = get_vertices_numpy(me)
3871     dvert_lay = bm.verts.layers.deform.active
3872     n_edges = len(me.edges)
3873     n_verts = len(me.vertices)
3874     group_index_x = ob.vertex_groups["x"].index
3875     group_index_y = ob.vertex_groups["y"].index
3876     group_index_module = ob.vertex_groups["module"].index
3877     group_index_values = ob.vertex_groups["values"].index
3878
3879     if not props.bool_cache:
3880         time_steps = props.time_steps
3881
3882         # store weight values
3883         if 'dB' in ob.vertex_groups: db = np.zeros(n_verts)
3884         if 'grad' in ob.vertex_groups: grad = np.zeros(n_verts)
3885         vec_x = np.zeros(n_verts)
3886         vec_y = np.zeros(n_verts)
3887         vec_module = np.zeros(n_verts)
3888         values = np.zeros(n_verts)
3889
3890         vec_x = bmesh_get_weight_numpy(group_index_x, dvert_lay, bm.verts)
3891         vec_y = bmesh_get_weight_numpy(group_index_y, dvert_lay, bm.verts)
3892         vec_module = bmesh_get_weight_numpy(group_index_module, dvert_lay, bm.verts)
3893         values = bmesh_get_weight_numpy(group_index_values, dvert_lay, bm.verts)
3894         field = np.concatenate((vec_x[:,None],vec_y[:,None],vec_y[:,None]*0),axis=1)
3895         field = field*2-1
3896         field[:,2] = 0
3897         edge_verts = get_edges_numpy(me)
3898
3899         id0 = edge_verts[:,0]
3900         id1 = edge_verts[:,1]
3901         vert0 = verts[id0]
3902         vert1 = verts[id1]
3903         vec = vert1-vert0
3904         edge_field = (field[id0] + field[id1])/2    # average vector associated to the edge
3905         print(vert0.shape)
3906         print(field.shape)
3907         print(edge_field.shape)
3908         # normalize vectors
3909         vec /= np.linalg.norm(vec,axis=1)[:,None]
3910         edge_field /= np.linalg.norm(edge_field,axis=1)[:,None]
3911         edge_flow = np.einsum('...j,...j', vec, edge_field)
3912         #sign = (edge_flow>0).astype(int)
3913         #values[edge_verts[sign]] += values[edge_verts[1-sign]]*
3914         #values[verts0] += values[verts1]*edge_flow
3915
3916         timeElapsed = time.time() - start
3917         print('       Preparation Time:',timeElapsed)
3918         start = time.time()
3919
3920         # Preserve energy
3921         mult = np.zeros(values.shape)
3922         #mult[id0] -= edge_flow
3923         #mult[id1] += edge_flow
3924         np.add.at(mult,id0,-edge_flow)
3925         np.add.at(mult,id1,edge_flow)
3926         print("mult")
3927         mult = scale/mult
3928         print(mult)
3929         print(np.sum(mult))
3930
3931
3932         #try:
3933         print(vec)
3934         print(edge_flow)
3935         print(edge_flow)
3936
3937         bool_run = False
3938         for j in range(props.cache_frame_start, props.cache_frame_end+1):
3939             start2 = time.time()
3940             print("{:6d} Reaction-Diffusion: {}".format(j, ob.name))
3941             if bool_run:
3942                 print(values)
3943                 #for i in range(1):
3944                 values = integrate_field(n_edges,id0,id1,values,edge_flow,mult,time_steps)
3945                 #values0 = values
3946                 #np.add.at(values, id0, values0[id1]*edge_flow*mult[id1])
3947                 #np.add.at(values, id1, -values0[id0]*edge_flow*mult[id0])
3948                 #np.add.at(values, id0, values0[id1]*edge_flow*mult)
3949                 #np.add.at(values, id1, -values0[id0]*edge_flow*mult)
3950                 #values[id1] += values0[id0]*edge_flow/mult[id1]*dt
3951                 #values[id0] -= values0[id1]*edge_flow/mult[id0]*dt
3952                 #values[id1] = edge_flow
3953                 #values[id1] += edge_flow
3954                 #a, b = numba_reaction_diffusion(n_verts, n_edges, edge_verts, a, b, _brush, _diff_a, _diff_b, _f, _k, dt, time_steps)
3955
3956                 '''
3957                 lap_a = np.zeros(n_verts)
3958                 lap_b = np.zeros(n_verts)
3959                 lap_a0 =  a[id1] -  a[id0]   # laplacian increment for first vertex of each edge
3960                 lap_b0 =  b[id1] -  b[id0]   # laplacian increment for first vertex of each edge
3961
3962                 np.add.at(lap_a, id0, lap_a0)
3963                 np.add.at(lap_b, id0, lap_b0)
3964                 np.add.at(lap_a, id1, -lap_a0)
3965                 np.add.at(lap_b, id1, -lap_b0)
3966                 '''
3967             else:
3968                 bool_run = True
3969
3970             if not(os.path.exists(folder)):
3971                 os.mkdir(folder)
3972             file_name = folder / "a_{:04d}".format(j)
3973             values.tofile(file_name)
3974             file_name = folder / "b_{:04d}".format(j)
3975             values.tofile(file_name)
3976
3977
3978             timeElapsed = time.time() - start2
3979             print('       Simulation Time:',timeElapsed)
3980
3981     if props.bool_mod: bpy.data.meshes.remove(me)
3982     bm.free()
3983     timeElapsed = time.time() - start
3984     print('       Closing Time:',timeElapsed)
3985
3986
3987
3988
3989
3990 class TISSUE_PT_reaction_diffusion(Panel):
3991     bl_space_type = 'PROPERTIES'
3992     bl_region_type = 'WINDOW'
3993     bl_context = "data"
3994     bl_label = "Tissue Reaction-Diffusion"
3995     bl_options = {'DEFAULT_CLOSED'}
3996
3997     @classmethod
3998     def poll(cls, context):
3999         return 'A' and 'B' in context.object.vertex_groups
4000
4001     def draw(self, context):
4002         reaction_diffusion_add_handler(self, context)
4003
4004         ob = context.object
4005         props = ob.reaction_diffusion_settings
4006         layout = self.layout
4007         col = layout.column(align=True)
4008         row = col.row(align=True)
4009         if not ("A" and "B" in ob.vertex_groups):
4010             row.operator("object.start_reaction_diffusion",
4011                         icon="EXPERIMENTAL",
4012                         text="Reaction-Diffusion")
4013         else:
4014             row.operator("object.start_reaction_diffusion",
4015                         icon="EXPERIMENTAL",
4016                         text="Reset Reaction-Diffusion")
4017             row = col.row(align=True)
4018             row.prop(props, "run", text="Run Reaction-Diffusion")
4019             col = layout.column(align=True)
4020             row = col.row(align=True)
4021             row.prop(props, "time_steps")
4022             row.prop(props, "dt")
4023             row.enabled = not props.bool_cache
4024             col.separator()
4025             row = col.row(align=True)
4026             col1 = row.column(align=True)
4027             col1.prop(props, "diff_a")
4028             col1.enabled = props.vertex_group_diff_a == '' and not props.bool_cache
4029             col1 = row.column(align=True)
4030             col1.prop(props, "diff_b")
4031             col1.enabled = props.vertex_group_diff_b == '' and not props.bool_cache
4032             row = col.row(align=True)
4033             row.prop(props, "diff_mult")
4034             row.enabled = props.vertex_group_scale == '' and not props.bool_cache
4035             #col.separator()
4036             row = col.row(align=True)
4037             col1 = row.column(align=True)
4038             col1.prop(props, "f")
4039             col1.enabled = props.vertex_group_f == '' and not props.bool_cache
4040             col1 = row.column(align=True)
4041             col1.prop(props, "k")
4042             col1.enabled = props.vertex_group_k == '' and not props.bool_cache
4043             col.separator()
4044             col.label(text='Cache:')
4045             #col.prop(props, "bool_cache")
4046             col.prop(props, "cache_dir", text='')
4047             col.separator()
4048             row = col.row(align=True)
4049             row.prop(props, "cache_frame_start")
4050             row.prop(props, "cache_frame_end")
4051             col.separator()
4052             if props.bool_cache:
4053                 col.operator("object.reaction_diffusion_free_data")
4054             else:
4055                 row = col.row(align=True)
4056                 row.operator("object.bake_reaction_diffusion")
4057                 file = bpy.context.blend_data.filepath
4058                 temp = bpy.context.preferences.filepaths.temporary_directory
4059                 if file == temp == props.cache_dir == '':
4060                     row.enabled = False
4061                     col.label(text="Cannot use cache", icon='ERROR')
4062                     col.label(text='please save the Blender or set a Cache directory')
4063                 col.prop(props, "fast_bake")
4064
4065             col.separator()
4066             col.label(text='Output attributes:')
4067             row = col.row(align=True)
4068             col2 = row.column(align=True)
4069             row2 = col2.row(align=True)
4070             row2.prop(props, "update_weight_a", icon='GROUP_VERTEX', text='A')
4071             row2.prop(props, "update_weight_b", icon='GROUP_VERTEX', text='B')
4072             col2.enabled = props.bool_cache
4073             row.separator()
4074             #row.prop(props, "update_colors_a", icon='GROUP_VCOL', text='A')
4075             #row.prop(props, "update_colors_b", icon='GROUP_VCOL', text='B')
4076             row.prop(props, "update_colors", icon='GROUP_VCOL', text='AB')
4077             row.separator()
4078             row.prop(props, "update_uv", icon='GROUP_UVS', text='AB')
4079             col.prop(props,'normalize')
4080
4081             #col.prop_search(props, 'vertex_group_diff_a', ob, "vertex_groups", text='Diff A')
4082             #col.prop_search(props, 'vertex_group_diff_b', ob, "vertex_groups", text='Diff B')
4083             #col.prop_search(props, 'vertex_group_scale', ob, "vertex_groups", text='Scale')
4084             #col.prop_search(props, 'vertex_group_f', ob, "vertex_groups", text='f')
4085             #col.prop_search(props, 'vertex_group_k', ob, "vertex_groups", text='k')
4086
4087
4088 class TISSUE_PT_reaction_diffusion_weight(Panel):
4089     bl_space_type = 'PROPERTIES'
4090     bl_region_type = 'WINDOW'
4091     bl_context = "data"
4092     bl_parent_id = "TISSUE_PT_reaction_diffusion"
4093     bl_label = "Vertex Groups"
4094     bl_options = {'DEFAULT_CLOSED'}
4095
4096     @classmethod
4097     def poll(cls, context):
4098         return 'A' and 'B' in context.object.vertex_groups
4099
4100     def draw(self, context):
4101         ob = context.object
4102         props = ob.reaction_diffusion_settings
4103         layout = self.layout
4104         #layout.use_property_split = True
4105         col = layout.column(align=True)
4106         col.prop(props, "bool_mod")
4107         if props.bool_mod and props.fast_bake:
4108             col.label(text="When Fast Bake is on, the modifiers", icon='ERROR')
4109             col.label(text="          are used only for the first frame")
4110         col.separator()
4111         insert_weight_parameter(col, ob, 'brush', text='Brush:')
4112         insert_weight_parameter(col, ob, 'diff_a', text='Diff A:')
4113         insert_weight_parameter(col, ob, 'diff_b', text='Diff B:')
4114         insert_weight_parameter(col, ob, 'scale', text='Scale:')
4115         insert_weight_parameter(col, ob, 'f', text='f:')
4116         insert_weight_parameter(col, ob, 'k', text='k:')
4117         col.enabled = not props.bool_cache
4118
4119 def insert_weight_parameter(col, ob, name, text=''):
4120     props = ob.reaction_diffusion_settings
4121     split = col.split(factor=0.25, align=True)
4122     col2 = split.column(align=True)
4123     col2.label(text=text)
4124     col2 = split.column(align=True)
4125     row2 = col2.row(align=True)
4126     row2.prop_search(props, 'vertex_group_' + name, ob, "vertex_groups", text='')
4127     if name != 'brush':
4128         row2.prop(props, "invert_vertex_group_" + name, text="", toggle=True, icon='ARROW_LEFTRIGHT')
4129     if 'vertex_group_' + name in props:
4130         if props['vertex_group_' + name] != '':
4131             if name == 'brush':
4132                 col2.prop(props, "brush_mult")
4133             else:
4134                 row2 = col2.row(align=True)
4135                 row2.prop(props, "min_" + name, text="Min")
4136                 row2 = col2.row(align=True)
4137                 row2.prop(props, "max_" + name, text="Max")
4138     col.separator()
4139
4140 def contour_edges_pattern(operator, c, verts_count, iso_val, vertices, normals, filtered_edges, weight, pattern_weight, bevel_weight):
4141     # vertices indexes
4142     id0 = filtered_edges[:,0]
4143     id1 = filtered_edges[:,1]
4144     # vertices weight
4145     w0 = weight[id0]
4146     w1 = weight[id1]
4147     # weight condition
4148     bool_w0 = w0 < iso_val
4149     bool_w1 = w1 < iso_val
4150
4151     # mask all edges that have one weight value below the iso value
4152     mask_new_verts = np.logical_xor(bool_w0, bool_w1)
4153     if not mask_new_verts.any():
4154         return np.array([[None]]), {}, np.array([[None]]), np.array([[None]])
4155
4156     id0 = id0[mask_new_verts]
4157     id1 = id1[mask_new_verts]
4158     # filter arrays
4159     v0 = vertices[id0]
4160     v1 = vertices[id1]
4161     n0 = normals[id0]
4162     n1 = normals[id1]
4163     w0 = w0[mask_new_verts]
4164     w1 = w1[mask_new_verts]
4165     pattern0 = pattern_weight[id0]
4166     pattern1 = pattern_weight[id1]
4167     try:
4168         bevel0 = bevel_weight[id0]
4169         bevel1 = bevel_weight[id1]
4170     except: pass
4171
4172     ### Spiral
4173     #edge_nor = (n0+n1)/2
4174     #shift = np.arctan2(edge_nor[:,0], edge_nor[:,1])/2/pi*delta_iso
4175
4176     #param = (iso_val + shift - w0)/(w1-w0)
4177     param = (iso_val - w0)/(w1-w0)
4178     # pattern displace
4179     #mult = 1 if c%2 == 0 else -1
4180     if c%(operator.in_steps + operator.out_steps) < operator.in_steps:
4181         mult = operator.in_displace
4182     else:
4183         mult = operator.out_displace
4184     pattern_value = pattern0 + (pattern1-pattern0)*param
4185     try:
4186         bevel_value = bevel0 + (bevel1-bevel0)*param
4187         bevel_value = np.expand_dims(bevel_value,axis=1)
4188     except: bevel_value = None
4189     disp = pattern_value * mult
4190
4191     param = np.expand_dims(param,axis=1)
4192     disp = np.expand_dims(disp,axis=1)
4193     verts = v0 + (v1-v0)*param
4194     norm = n0 + (n1-n0)*param
4195     if operator.limit_z: disp *= 1-abs(np.expand_dims(norm[:,2], axis=1))
4196     verts = verts + norm*disp
4197     #verts = verts[np.flip(np.argsort(shift))]
4198     #verts = verts[np.argsort(shift)]
4199
4200     # indexes of edges with new vertices
4201     edges_index = filtered_edges[mask_new_verts][:,2]
4202
4203     # remove all edges completely below the iso value
4204     #mask_edges = np.logical_not(np.logical_and(bool_w0, bool_w1))
4205     #filtered_edges = filtered_edges[mask_edges]
4206     return filtered_edges, edges_index, verts, bevel_value
4207
4208 def contour_bmesh(me, bm, weight, iso_val):
4209     bm.verts.ensure_lookup_table()
4210     bm.edges.ensure_lookup_table()
4211     bm.faces.ensure_lookup_table()
4212
4213     # store weight values
4214
4215     vertices = get_vertices_numpy(me)
4216     faces_mask = np.array(bm.faces)
4217     filtered_edges = get_edges_id_numpy(me)
4218     n_verts = len(bm.verts)
4219
4220     #############################
4221
4222     # vertices indexes
4223     id0 = filtered_edges[:,0]
4224     id1 = filtered_edges[:,1]
4225     # vertices weight
4226     w0 = weight[id0]
4227     w1 = weight[id1]
4228     # weight condition
4229     bool_w0 = w0 < iso_val
4230     bool_w1 = w1 < iso_val
4231
4232     # mask all edges that have one weight value below the iso value
4233     mask_new_verts = np.logical_xor(bool_w0, bool_w1)
4234     if not mask_new_verts.any(): return np.array([[None]]), {}, np.array([[None]])
4235
4236     id0 = id0[mask_new_verts]
4237     id1 = id1[mask_new_verts]
4238     # filter arrays
4239     v0 = vertices[id0]
4240     v1 = vertices[id1]
4241     w0 = w0[mask_new_verts]
4242     w1 = w1[mask_new_verts]
4243     param = (iso_val-w0)/(w1-w0)
4244     param = np.expand_dims(param,axis=1)
4245     verts = v0 + (v1-v0)*param
4246
4247     # indexes of edges with new vertices
4248     #edges_index = filtered_edges[mask_new_verts][:,2]
4249
4250     edges_id = {}
4251     for i, e in enumerate(filtered_edges):
4252         #edges_id[id] = i + n_verts
4253         edges_id['{}_{}'.format(e[0],e[1])] = i + n_verts
4254         edges_id['{}_{}'.format(e[1],e[0])] = i + n_verts
4255
4256
4257     '''
4258     for e in filtered_edges:
4259         id0 = e.verts[0].index
4260         id1 = e.verts[1].index
4261         w0 = weight[id0]
4262         w1 = weight[id1]
4263
4264         if w0 == w1: continue
4265         elif w0 > iso_val and w1 > iso_val:
4266             continue
4267         elif w0 < iso_val and w1 < iso_val: continue
4268         elif w0 == iso_val or w1 == iso_val: continue
4269         else:
4270             v0 = me0.vertices[id0].co
4271             v1 = me0.vertices[id1].co
4272             v = v0.lerp(v1, (iso_val-w0)/(w1-w0))
4273             delete_edges.append(e)
4274             verts.append(v)
4275             edges_id[str(id0)+"_"+str(id1)] = count
4276             edges_id[str(id1)+"_"+str(id0)] = count
4277             count += 1
4278     '''
4279
4280     splitted_faces = []
4281
4282     switch = False
4283     # splitting faces
4284     for f in faces_mask:
4285         # create sub-faces slots. Once a new vertex is reached it will
4286         # change slot, storing the next vertices for a new face.
4287         build_faces = [[],[]]
4288         #switch = False
4289         verts0 = list(me.polygons[f.index].vertices)
4290         verts1 = list(verts0)
4291         verts1.append(verts1.pop(0)) # shift list
4292         for id0, id1 in zip(verts0, verts1):
4293
4294             # add first vertex to active slot
4295             build_faces[switch].append(id0)
4296
4297             # try to split edge
4298             try:
4299                 # check if the edge must be splitted
4300                 new_vert = edges_id['{}_{}'.format(id0,id1)]
4301                 # add new vertex
4302                 build_faces[switch].append(new_vert)
4303                 # if there is an open face on the other slot
4304                 if len(build_faces[not switch]) > 0:
4305                     # store actual face
4306                     splitted_faces.append(build_faces[switch])
4307                     # reset actual faces and switch
4308                     build_faces[switch] = []
4309                     # change face slot
4310                 switch = not switch
4311                 # continue previous face
4312                 build_faces[switch].append(new_vert)
4313             except: pass
4314         if len(build_faces[not switch]) == 2:
4315             build_faces[not switch].append(id0)
4316         if len(build_faces[not switch]) > 2:
4317             splitted_faces.append(build_faces[not switch])
4318         # add last face
4319         splitted_faces.append(build_faces[switch])
4320
4321     # adding new vertices use fast local method access
4322     _new_vert = bm.verts.new
4323     for v in verts: _new_vert(v)
4324     bm.verts.ensure_lookup_table()
4325
4326     # deleting old edges/faces
4327     bm.edges.ensure_lookup_table()
4328     remove_edges = [bm.edges[i] for i in filtered_edges[:,2]]
4329     #for e in remove_edges: bm.edges.remove(e)
4330     #for e in delete_edges: bm.edges.remove(e)
4331
4332     bm.verts.ensure_lookup_table()
4333     # adding new faces use fast local method access
4334     _new_face = bm.faces.new
4335     missed_faces = []
4336     for f in splitted_faces:
4337         try:
4338             face_verts = [bm.verts[i] for i in f]
4339             _new_face(face_verts)
4340         except:
4341             missed_faces.append(f)
4342
4343     #me = bpy.data.meshes.new('_tissue_tmp_')
4344     bm.to_mesh(me)
4345     weight = np.concatenate((weight, np.ones(len(verts))*iso_val))
4346
4347     return me, bm, weight
4348
4349
4350
4351
4352 class tissue_weight_streamlines(Operator):
4353     bl_idname = "object.tissue_weight_streamlines"
4354     bl_label = "Streamlines Curves"
4355     bl_description = ("")
4356     bl_options = {'REGISTER', 'UNDO'}
4357
4358     mode : EnumProperty(
4359         items=(
4360             ('VERTS', "Verts", "Follow vertices"),
4361             ('EDGES', "Edges", "Follow Edges")
4362             ),
4363         default='VERTS',
4364         name="Streamlines path mode"
4365         )
4366
4367     interpolation : EnumProperty(
4368         items=(
4369             ('POLY', "Poly", "Generate Polylines"),
4370             ('NURBS', "NURBS", "Generate Nurbs curves")
4371             ),
4372         default='POLY',
4373         name="Interpolation mode"
4374         )
4375
4376     use_modifiers : BoolProperty(
4377         name="Use Modifiers", default=True,
4378         description="Apply all the modifiers")
4379
4380     use_selected : BoolProperty(
4381         name="Use Selected Vertices", default=False,
4382         description="Use selected vertices as Seed")
4383
4384     same_weight : BoolProperty(
4385         name="Same Weight", default=True,
4386         description="Continue the streamlines when the weight is the same")
4387
4388     min_iso : FloatProperty(
4389         name="Min Value", default=0., soft_min=0, soft_max=1,
4390         description="Minimum weight value")
4391     max_iso : FloatProperty(
4392         name="Max Value", default=1, soft_min=0, soft_max=1,
4393         description="Maximum weight value")
4394
4395     rand_seed : IntProperty(
4396         name="Seed", default=0, min=0, soft_max=10,
4397         description="Random Seed")
4398     n_curves : IntProperty(
4399         name="Curves", default=50, soft_min=1, soft_max=100000,
4400         description="Number of Curves")
4401     min_rad = 1
4402     max_rad = 1
4403
4404     pos_steps : IntProperty(
4405         name="High Steps", default=50, min=0, soft_max=100,
4406         description="Number of steps in the direction of high weight")
4407     neg_steps : IntProperty(
4408         name="Low Steps", default=50, min=0, soft_max=100,
4409         description="Number of steps in the direction of low weight")
4410
4411     bevel_depth : FloatProperty(
4412         name="Bevel Depth", default=0, min=0, soft_max=1,
4413         description="")
4414     min_bevel_depth : FloatProperty(
4415         name="Min Bevel Depth", default=0.1, min=0, soft_max=1,
4416         description="")
4417     max_bevel_depth : FloatProperty(
4418         name="Max Bevel Depth", default=1, min=0, soft_max=1,
4419         description="")
4420
4421     rand_dir : FloatProperty(
4422         name="Randomize", default=0, min=0, max=1,
4423         description="Randomize streamlines directions (Slower)")
4424
4425     vertex_group_seeds : StringProperty(
4426         name="Displace", default='',
4427         description="Vertex Group used for pattern displace")
4428
4429     vertex_group_bevel : StringProperty(
4430         name="Bevel", default='',
4431         description="Variable Bevel depth")
4432
4433     object_name : StringProperty(
4434         name="Active Object", default='',
4435         description="")
4436
4437     try: vg_name = bpy.context.object.vertex_groups.active.name
4438     except: vg_name = ''
4439
4440     vertex_group_streamlines : StringProperty(
4441         name="Flow", default=vg_name,
4442         description="Vertex Group used for streamlines")
4443
4444     @classmethod
4445     def poll(cls, context):
4446         ob = context.object
4447         return ob and len(ob.vertex_groups) > 0 or ob.type == 'CURVE'
4448
4449     def invoke(self, context, event):
4450         return context.window_manager.invoke_props_dialog(self, width=250)
4451
4452     def draw(self, context):
4453         if not context.object.type == 'CURVE':
4454             self.object_name = context.object.name
4455         ob = bpy.data.objects[self.object_name]
4456         if self.vertex_group_streamlines not in [vg.name for vg in ob.vertex_groups]:
4457             self.vertex_group_streamlines = ob.vertex_groups.active.name
4458         layout = self.layout
4459         col = layout.column(align=True)
4460         row = col.row(align=True)
4461         row.prop(self, 'mode', expand=True,
4462             slider=True, toggle=False, icon_only=False, event=False,
4463             full_event=False, emboss=True, index=-1)
4464         col.prop(self, "use_modifiers")
4465         col.label(text="Streamlines Curves:")
4466         row = col.row(align=True)
4467         row.prop(self, 'interpolation', expand=True,
4468             slider=True, toggle=False, icon_only=False, event=False,
4469             full_event=False, emboss=True, index=-1)
4470         col.separator()
4471         col.prop_search(self, 'vertex_group_streamlines', ob, "vertex_groups", text='')
4472         if not (self.use_selected or context.mode == 'EDIT_MESH'):
4473             row = col.row(align=True)
4474             row.prop(self,'n_curves')
4475             #row.enabled = context.mode != 'EDIT_MESH'
4476             row = col.row(align=True)
4477             row.prop(self,'rand_seed')
4478             #row.enabled = context.mode != 'EDIT_MESH'
4479         row = col.row(align=True)
4480         row.prop(self,'neg_steps')
4481         row.prop(self,'pos_steps')
4482         #row = col.row(align=True)
4483         #row.prop(self,'min_iso')
4484         #row.prop(self,'max_iso')
4485         col.prop(self, "same_weight")
4486         col.separator()
4487         col.label(text='Curves Bevel:')
4488         col.prop_search(self, 'vertex_group_bevel', ob, "vertex_groups", text='')
4489         if self.vertex_group_bevel != '':
4490             row = col.row(align=True)
4491             row.prop(self,'min_bevel_depth')
4492             row.prop(self,'max_bevel_depth')
4493         else:
4494             col.prop(self,'bevel_depth')
4495         col.separator()
4496         col.prop(self, "rand_dir")
4497
4498     def execute(self, context):
4499         start_time = timeit.default_timer()
4500         try:
4501             check = context.object.vertex_groups[0]
4502         except:
4503             self.report({'ERROR'}, "The object doesn't have Vertex Groups")
4504             return {'CANCELLED'}
4505         ob = bpy.data.objects[self.object_name]
4506         ob.select_set(False)
4507
4508
4509         seeds = []
4510
4511         if bpy.context.mode == 'EDIT_MESH':
4512             self.use_selected = True
4513             bpy.ops.object.mode_set(mode='OBJECT')
4514             #ob = bpy.context.object
4515             #me = simple_to_mesh(ob)
4516         ob = convert_object_to_mesh(ob, apply_modifiers=self.use_modifiers)
4517         #dg = context.evaluated_depsgraph_get()
4518         #ob = ob.evaluated_get(dg)
4519         me = ob.data
4520
4521         if self.use_selected:
4522             # generate new bmesh
4523             bm = bmesh.new()
4524             bm.from_mesh(me)
4525             print(len(me.vertices))
4526             #for v in me.vertices:
4527             #    if v.select: seeds.append(v.index)
4528             for v in bm.verts:
4529                 if v.select: seeds.append(v.index)
4530             bm.free()
4531         n_verts = len(me.vertices)
4532         n_edges = len(me.edges)
4533         n_faces = len(me.polygons)
4534
4535         # store weight values
4536         try:
4537             weight = get_weight_numpy(ob.vertex_groups[self.vertex_group_streamlines], n_verts)
4538         except:
4539             bpy.data.objects.remove(ob)
4540             self.report({'ERROR'}, "Please select a Vertex Group for streamlines")
4541             return {'CANCELLED'}
4542
4543         variable_bevel = False
4544         bevel_weight = None
4545         bevel_depth = self.bevel_depth
4546         try:
4547             if self.min_bevel_depth == self.max_bevel_depth:
4548                 #bevel_weight = np.ones((n_verts))
4549                 bevel_depth = self.min_bevel_depth
4550             else:
4551                 b0 = min(self.min_bevel_depth, self.max_bevel_depth)
4552                 b1 = max(self.min_bevel_depth, self.max_bevel_depth)
4553                 bevel_weight = get_weight_numpy(ob.vertex_groups[self.vertex_group_bevel], n_verts)
4554                 if self.min_bevel_depth > self.max_bevel_depth:
4555                     bevel_weight = 1-bevel_weight
4556                 bevel_weight = b0/b1 + bevel_weight*((b1-b0)/b1)
4557                 bevel_depth = b1
4558             variable_bevel = True
4559         except:
4560             pass#bevel_weight = np.ones((n_verts))
4561
4562
4563         if not seeds:
4564             np.random.seed(self.rand_seed)
4565             seeds = np.random.randint(n_verts, size=self.n_curves)
4566
4567         #weight = np.array(get_weight(ob.vertex_groups.active, n_verts))
4568
4569         curves_pts = []
4570         curves_weight = []
4571
4572         neigh = [[] for i in range(n_verts)]
4573         if self.mode == 'EDGES':
4574             # store neighbors
4575             for e in me.edges:
4576                 ev = e.vertices
4577                 neigh[ev[0]].append(ev[1])
4578                 neigh[ev[1]].append(ev[0])
4579
4580         elif self.mode == 'VERTS':
4581             # store neighbors
4582             for p in me.polygons:
4583                 face_verts = [v for v in p.vertices]
4584                 n_face_verts = len(face_verts)
4585                 for i in range(n_face_verts):
4586                     fv = face_verts.copy()
4587                     neigh[fv.pop(i)] += fv
4588
4589         neigh_weight = [weight[n].tolist() for n in neigh]
4590
4591         # evaluate direction
4592         next_vert = [-1]*n_verts
4593
4594         if self.rand_dir > 0:
4595             for i in range(n_verts):
4596                 n = neigh[i]
4597                 nw = neigh_weight[i]
4598                 sorted_nw = neigh_weight[i].copy()
4599                 sorted_nw.sort()
4600                 for w in sorted_nw:
4601                     neigh[i] = [n[nw.index(w)] for w in sorted_nw]
4602         else:
4603             if self.pos_steps > 0:
4604                 for i in range(n_verts):
4605                     n = neigh[i]
4606                     nw = neigh_weight[i]
4607                     max_w = max(nw)
4608                     if self.same_weight:
4609                         if max_w >= weight[i]:
4610                             next_vert[i] = n[nw.index(max(nw))]
4611                     else:
4612                         if max_w > weight[i]:
4613                             next_vert[i] = n[nw.index(max(nw))]
4614
4615             if self.neg_steps > 0:
4616                 prev_vert = [-1]*n_verts
4617                 for i in range(n_verts):
4618                     n = neigh[i]
4619                     nw = neigh_weight[i]
4620                     min_w = min(nw)
4621                     if self.same_weight:
4622                         if min_w <= weight[i]:
4623                             prev_vert[i] = n[nw.index(min(nw))]
4624                     else:
4625                         if min_w < weight[i]:
4626                             prev_vert[i] = n[nw.index(min(nw))]
4627
4628         co = [0]*3*n_verts
4629         me.vertices.foreach_get('co', co)
4630         co = np.array(co).reshape((-1,3))
4631
4632         # create streamlines
4633         curves = []
4634         for i in seeds:
4635             next_pts = [i]
4636             for j in range(self.pos_steps):
4637                 if self.rand_dir > 0:
4638                     n = neigh[next_pts[-1]]
4639                     next = n[int(len(n) * (1-random.random() * self.rand_dir))]
4640                 else:
4641                     next = next_vert[next_pts[-1]]
4642                 if next > 0:
4643                     if next not in next_pts: next_pts.append(next)
4644                 else: break
4645
4646             prev_pts = [i]
4647             for j in range(self.neg_steps):
4648                 if self.rand_dir > 0:
4649                     n = neigh[prev_pts[-1]]
4650                     prev = n[int(len(n) * random.random() * self.rand_dir)]
4651                 else:
4652                     prev = prev_vert[prev_pts[-1]]
4653                 if prev > 0:
4654                     if prev not in prev_pts:
4655                         prev_pts.append(prev)
4656                 else: break
4657
4658             next_pts = np.array(next_pts).astype('int')
4659             prev_pts = np.flip(prev_pts[1:]).astype('int')
4660             all_pts = np.concatenate((prev_pts, next_pts))
4661             if len(all_pts) > 1:
4662                 curves.append(all_pts)
4663         crv = nurbs_from_vertices(curves, co, bevel_weight, ob.name + '_Streamlines', True, self.interpolation)
4664         crv.data.bevel_depth = bevel_depth
4665         crv.matrix_world = ob.matrix_world
4666         bpy.data.objects.remove(ob)
4667
4668         print("Streamlines Curves, total time: " + str(timeit.default_timer() - start_time) + " sec")
4669         return {'FINISHED'}