Cleanup: style, use braces for blenkernel
[blender.git] / source / blender / blenkernel / intern / particle_child.c
1 /*
2  * This program is free software; you can redistribute it and/or
3  * modify it under the terms of the GNU General Public License
4  * as published by the Free Software Foundation; either version 2
5  * of the License, or (at your option) any later version.
6  *
7  * This program is distributed in the hope that it will be useful,
8  * but WITHOUT ANY WARRANTY; without even the implied warranty of
9  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
10  * GNU General Public License for more details.
11  *
12  * You should have received a copy of the GNU General Public License
13  * along with this program; if not, write to the Free Software Foundation,
14  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
15  *
16  * The Original Code is Copyright (C) Blender Foundation
17  * All rights reserved.
18  */
19
20 /** \file
21  * \ingroup bke
22  */
23
24 #include "BLI_math.h"
25 #include "BLI_noise.h"
26
27 #include "DNA_material_types.h"
28
29 #include "BKE_colortools.h"
30 #include "BKE_particle.h"
31
32 #include "particle_private.h"
33
34 /* ------------------------------------------------------------------------- */
35
36 typedef struct ParticlePathIterator {
37   ParticleCacheKey *key;
38   int index;
39   float time;
40
41   ParticleCacheKey *parent_key;
42   float parent_rotation[4];
43 } ParticlePathIterator;
44
45 static void psys_path_iter_get(ParticlePathIterator *iter,
46                                ParticleCacheKey *keys,
47                                int totkeys,
48                                ParticleCacheKey *parent,
49                                int index)
50 {
51   BLI_assert(index >= 0 && index < totkeys);
52
53   iter->key = keys + index;
54   iter->index = index;
55   iter->time = (float)index / (float)(totkeys - 1);
56
57   if (parent) {
58     iter->parent_key = parent + index;
59     if (index > 0) {
60       mul_qt_qtqt(iter->parent_rotation, iter->parent_key->rot, parent->rot);
61     }
62     else {
63       copy_qt_qt(iter->parent_rotation, parent->rot);
64     }
65   }
66   else {
67     iter->parent_key = NULL;
68     unit_qt(iter->parent_rotation);
69   }
70 }
71
72 typedef struct ParticlePathModifier {
73   struct ParticlePathModifier *next, *prev;
74
75   void (*apply)(ParticleCacheKey *keys, int totkeys, ParticleCacheKey *parent_keys);
76 } ParticlePathModifier;
77
78 /* ------------------------------------------------------------------------- */
79
80 static void do_kink_spiral_deform(ParticleKey *state,
81                                   const float dir[3],
82                                   const float kink[3],
83                                   float time,
84                                   float freq,
85                                   float shape,
86                                   float amplitude,
87                                   const float spiral_start[3])
88 {
89   float result[3];
90
91   CLAMP(time, 0.f, 1.f);
92
93   copy_v3_v3(result, state->co);
94
95   {
96     /* Creates a logarithmic spiral:
97      *   r(theta) = a * exp(b * theta)
98      *
99      * The "density" parameter b is defined by the shape parameter
100      * and goes up to the Golden Spiral for 1.0
101      * https://en.wikipedia.org/wiki/Golden_spiral
102      */
103     const float b = shape * (1.0f + sqrtf(5.0f)) / (float)M_PI * 0.25f;
104     /* angle of the spiral against the curve (rotated opposite to make a smooth transition) */
105     const float start_angle = ((b != 0.0f) ? atanf(1.0f / b) : (float)-M_PI_2) +
106                               (b > 0.0f ? -(float)M_PI_2 : (float)M_PI_2);
107
108     float spiral_axis[3], rot[3][3];
109     float vec[3];
110
111     float theta = freq * time * 2.0f * (float)M_PI;
112     float radius = amplitude * expf(b * theta);
113
114     /* a bit more intuitive than using negative frequency for this */
115     if (amplitude < 0.0f) {
116       theta = -theta;
117     }
118
119     cross_v3_v3v3(spiral_axis, dir, kink);
120     normalize_v3(spiral_axis);
121
122     mul_v3_v3fl(vec, kink, -radius);
123
124     axis_angle_normalized_to_mat3(rot, spiral_axis, theta);
125     mul_m3_v3(rot, vec);
126
127     madd_v3_v3fl(vec, kink, amplitude);
128
129     axis_angle_normalized_to_mat3(rot, spiral_axis, -start_angle);
130     mul_m3_v3(rot, vec);
131
132     add_v3_v3v3(result, spiral_start, vec);
133   }
134
135   copy_v3_v3(state->co, result);
136 }
137
138 static void do_kink_spiral(ParticleThreadContext *ctx,
139                            ParticleTexture *ptex,
140                            const float parent_orco[3],
141                            ChildParticle *cpa,
142                            const float orco[3],
143                            float hairmat[4][4],
144                            ParticleCacheKey *keys,
145                            ParticleCacheKey *parent_keys,
146                            int *r_totkeys,
147                            float *r_max_length)
148 {
149   struct ParticleSettings *part = ctx->sim.psys->part;
150   const int seed = ctx->sim.psys->child_seed + (int)(cpa - ctx->sim.psys->child);
151   const int totkeys = ctx->segments + 1;
152   const int extrakeys = ctx->extra_segments;
153
154   float kink_amp_random = part->kink_amp_random;
155   float kink_amp = part->kink_amp *
156                    (1.0f - kink_amp_random * psys_frand(ctx->sim.psys, 93541 + seed));
157   float kink_freq = part->kink_freq;
158   float kink_shape = part->kink_shape;
159   float kink_axis_random = part->kink_axis_random;
160   float rough1 = part->rough1;
161   float rough2 = part->rough2;
162   float rough_end = part->rough_end;
163
164   ParticlePathIterator iter;
165   ParticleCacheKey *key;
166   int k;
167
168   float dir[3];
169   float spiral_start[3] = {0.0f, 0.0f, 0.0f};
170   float spiral_start_time = 0.0f;
171   float spiral_par_co[3] = {0.0f, 0.0f, 0.0f};
172   float spiral_par_vel[3] = {0.0f, 0.0f, 0.0f};
173   float spiral_par_rot[4] = {1.0f, 0.0f, 0.0f, 0.0f};
174   float totlen;
175   float cut_time;
176   int start_index = 0, end_index = 0;
177   float kink_base[3];
178
179   if (ptex) {
180     kink_amp *= ptex->kink_amp;
181     kink_freq *= ptex->kink_freq;
182     rough1 *= ptex->rough1;
183     rough2 *= ptex->rough2;
184     rough_end *= ptex->roughe;
185   }
186
187   cut_time = (totkeys - 1) * ptex->length;
188   zero_v3(spiral_start);
189
190   for (k = 0, key = keys; k < totkeys - 1; k++, key++) {
191     if ((float)(k + 1) >= cut_time) {
192       float fac = cut_time - (float)k;
193       ParticleCacheKey *par = parent_keys + k;
194
195       start_index = k + 1;
196       end_index = start_index + extrakeys;
197
198       spiral_start_time = ((float)k + fac) / (float)(totkeys - 1);
199       interp_v3_v3v3(spiral_start, key->co, (key + 1)->co, fac);
200
201       interp_v3_v3v3(spiral_par_co, par->co, (par + 1)->co, fac);
202       interp_v3_v3v3(spiral_par_vel, par->vel, (par + 1)->vel, fac);
203       interp_qt_qtqt(spiral_par_rot, par->rot, (par + 1)->rot, fac);
204
205       break;
206     }
207   }
208
209   zero_v3(dir);
210
211   zero_v3(kink_base);
212   kink_base[part->kink_axis] = 1.0f;
213   mul_mat3_m4_v3(ctx->sim.ob->obmat, kink_base);
214
215   /* Fill in invariant part of modifier context. */
216   ParticleChildModifierContext modifier_ctx = {NULL};
217   modifier_ctx.thread_ctx = ctx;
218   modifier_ctx.sim = &ctx->sim;
219   modifier_ctx.ptex = ptex;
220   modifier_ctx.cpa = cpa;
221   modifier_ctx.orco = orco;
222   modifier_ctx.parent_keys = parent_keys;
223
224   for (k = 0, key = keys; k < end_index; k++, key++) {
225     float par_time;
226     float *par_co, *par_vel, *par_rot;
227
228     psys_path_iter_get(&iter, keys, end_index, NULL, k);
229     if (k < start_index) {
230       sub_v3_v3v3(dir, (key + 1)->co, key->co);
231       normalize_v3(dir);
232
233       par_time = (float)k / (float)(totkeys - 1);
234       par_co = parent_keys[k].co;
235       par_vel = parent_keys[k].vel;
236       par_rot = parent_keys[k].rot;
237     }
238     else {
239       float spiral_time = (float)(k - start_index) / (float)(extrakeys - 1);
240       float kink[3], tmp[3];
241
242       /* use same time value for every point on the spiral */
243       par_time = spiral_start_time;
244       par_co = spiral_par_co;
245       par_vel = spiral_par_vel;
246       par_rot = spiral_par_rot;
247
248       project_v3_v3v3(tmp, kink_base, dir);
249       sub_v3_v3v3(kink, kink_base, tmp);
250       normalize_v3(kink);
251
252       if (kink_axis_random > 0.0f) {
253         float a = kink_axis_random * (psys_frand(ctx->sim.psys, 7112 + seed) * 2.0f - 1.0f) *
254                   (float)M_PI;
255         float rot[3][3];
256
257         axis_angle_normalized_to_mat3(rot, dir, a);
258         mul_m3_v3(rot, kink);
259       }
260
261       do_kink_spiral_deform((ParticleKey *)key,
262                             dir,
263                             kink,
264                             spiral_time,
265                             kink_freq,
266                             kink_shape,
267                             kink_amp,
268                             spiral_start);
269     }
270
271     /* Fill in variant part of modifier context. */
272     modifier_ctx.par_co = par_co;
273     modifier_ctx.par_vel = par_vel;
274     modifier_ctx.par_rot = par_rot;
275     modifier_ctx.par_orco = parent_orco;
276
277     /* Apply different deformations to the child path/ */
278     do_child_modifiers(&modifier_ctx, hairmat, (ParticleKey *)key, par_time);
279   }
280
281   totlen = 0.0f;
282   for (k = 0, key = keys; k < end_index - 1; k++, key++) {
283     totlen += len_v3v3((key + 1)->co, key->co);
284   }
285
286   *r_totkeys = end_index;
287   *r_max_length = totlen;
288 }
289
290 /* ------------------------------------------------------------------------- */
291
292 static bool check_path_length(int k,
293                               ParticleCacheKey *keys,
294                               ParticleCacheKey *key,
295                               float max_length,
296                               float step_length,
297                               float *cur_length,
298                               float dvec[3])
299 {
300   if (*cur_length + step_length > max_length) {
301     sub_v3_v3v3(dvec, key->co, (key - 1)->co);
302     mul_v3_fl(dvec, (max_length - *cur_length) / step_length);
303     add_v3_v3v3(key->co, (key - 1)->co, dvec);
304     keys->segments = k;
305     /* something over the maximum step value */
306     return false;
307   }
308   else {
309     *cur_length += step_length;
310     return true;
311   }
312 }
313
314 void psys_apply_child_modifiers(ParticleThreadContext *ctx,
315                                 struct ListBase *modifiers,
316                                 ChildParticle *cpa,
317                                 ParticleTexture *ptex,
318                                 const float orco[3],
319                                 float hairmat[4][4],
320                                 ParticleCacheKey *keys,
321                                 ParticleCacheKey *parent_keys,
322                                 const float parent_orco[3])
323 {
324   struct ParticleSettings *part = ctx->sim.psys->part;
325   struct Material *ma = ctx->ma;
326   const bool draw_col_ma = (part->draw_col == PART_DRAW_COL_MAT);
327   const bool use_length_check = !ELEM(part->kink, PART_KINK_SPIRAL);
328
329   ParticlePathModifier *mod;
330   ParticleCacheKey *key;
331   int totkeys, k;
332   float max_length;
333
334   /* TODO for the future: use true particle modifiers that work on the whole curve */
335
336   (void)modifiers;
337   (void)mod;
338
339   if (part->kink == PART_KINK_SPIRAL) {
340     do_kink_spiral(
341         ctx, ptex, parent_orco, cpa, orco, hairmat, keys, parent_keys, &totkeys, &max_length);
342     keys->segments = totkeys - 1;
343   }
344   else {
345     /* Fill in invariant part of modifier context. */
346     ParticleChildModifierContext modifier_ctx = {NULL};
347     modifier_ctx.thread_ctx = ctx;
348     modifier_ctx.sim = &ctx->sim;
349     modifier_ctx.ptex = ptex;
350     modifier_ctx.cpa = cpa;
351     modifier_ctx.orco = orco;
352     modifier_ctx.parent_keys = parent_keys;
353
354     totkeys = ctx->segments + 1;
355     max_length = ptex->length;
356
357     for (k = 0, key = keys; k < totkeys; k++, key++) {
358       ParticlePathIterator iter;
359       psys_path_iter_get(&iter, keys, totkeys, parent_keys, k);
360
361       ParticleKey *par = (ParticleKey *)iter.parent_key;
362
363       /* Fill in variant part of modifier context. */
364       modifier_ctx.par_co = par->co;
365       modifier_ctx.par_vel = par->vel;
366       modifier_ctx.par_rot = iter.parent_rotation;
367       modifier_ctx.par_orco = parent_orco;
368
369       /* Apply different deformations to the child path. */
370       do_child_modifiers(&modifier_ctx, hairmat, (ParticleKey *)key, iter.time);
371     }
372   }
373
374   {
375     const float step_length = 1.0f / (float)(totkeys - 1);
376     float cur_length = 0.0f;
377
378     if (max_length <= 0.0f) {
379       keys->segments = -1;
380       totkeys = 0;
381     }
382
383     /* we have to correct velocity because of kink & clump */
384     for (k = 0, key = keys; k < totkeys; ++k, ++key) {
385       if (k >= 2) {
386         sub_v3_v3v3((key - 1)->vel, key->co, (key - 2)->co);
387         mul_v3_fl((key - 1)->vel, 0.5);
388       }
389
390       if (use_length_check && k > 0) {
391         float dvec[3];
392         /* check if path needs to be cut before actual end of data points */
393         if (!check_path_length(k, keys, key, max_length, step_length, &cur_length, dvec)) {
394           /* last key */
395           sub_v3_v3v3(key->vel, key->co, (key - 1)->co);
396           if (ma && draw_col_ma) {
397             copy_v3_v3(key->col, &ma->r);
398           }
399           break;
400         }
401       }
402       if (k == totkeys - 1) {
403         /* last key */
404         sub_v3_v3v3(key->vel, key->co, (key - 1)->co);
405       }
406
407       if (ma && draw_col_ma) {
408         copy_v3_v3(key->col, &ma->r);
409       }
410     }
411   }
412 }
413
414 /* ------------------------------------------------------------------------- */
415
416 void do_kink(ParticleKey *state,
417              const float par_co[3],
418              const float par_vel[3],
419              const float par_rot[4],
420              float time,
421              float freq,
422              float shape,
423              float amplitude,
424              float flat,
425              short type,
426              short axis,
427              float obmat[4][4],
428              int smooth_start)
429 {
430   float kink[3] = {1.f, 0.f, 0.f}, par_vec[3], q1[4] = {1.f, 0.f, 0.f, 0.f};
431   float t, dt = 1.f, result[3];
432
433   if (ELEM(type, PART_KINK_NO, PART_KINK_SPIRAL)) {
434     return;
435   }
436
437   CLAMP(time, 0.f, 1.f);
438
439   if (shape != 0.0f && !ELEM(type, PART_KINK_BRAID)) {
440     if (shape < 0.0f) {
441       time = (float)pow(time, 1.f + shape);
442     }
443     else {
444       time = (float)pow(time, 1.f / (1.f - shape));
445     }
446   }
447
448   t = time * freq * (float)M_PI;
449
450   if (smooth_start) {
451     dt = fabsf(t);
452     /* smooth the beginning of kink */
453     CLAMP(dt, 0.f, (float)M_PI);
454     dt = sinf(dt / 2.f);
455   }
456
457   if (!ELEM(type, PART_KINK_RADIAL)) {
458     float temp[3];
459
460     kink[axis] = 1.f;
461
462     if (obmat) {
463       mul_mat3_m4_v3(obmat, kink);
464     }
465
466     mul_qt_v3(par_rot, kink);
467
468     /* make sure kink is normal to strand */
469     project_v3_v3v3(temp, kink, par_vel);
470     sub_v3_v3(kink, temp);
471     normalize_v3(kink);
472   }
473
474   copy_v3_v3(result, state->co);
475   sub_v3_v3v3(par_vec, par_co, state->co);
476
477   switch (type) {
478     case PART_KINK_CURL: {
479       float curl_offset[3];
480
481       /* rotate kink vector around strand tangent */
482       mul_v3_v3fl(curl_offset, kink, amplitude);
483       axis_angle_to_quat(q1, par_vel, t);
484       mul_qt_v3(q1, curl_offset);
485
486       interp_v3_v3v3(par_vec, state->co, par_co, flat);
487       add_v3_v3v3(result, par_vec, curl_offset);
488       break;
489     }
490     case PART_KINK_RADIAL: {
491       if (flat > 0.f) {
492         float proj[3];
493         /* flatten along strand */
494         project_v3_v3v3(proj, par_vec, par_vel);
495         madd_v3_v3fl(result, proj, flat);
496       }
497
498       madd_v3_v3fl(result, par_vec, -amplitude * sinf(t));
499       break;
500     }
501     case PART_KINK_WAVE: {
502       madd_v3_v3fl(result, kink, amplitude * sinf(t));
503
504       if (flat > 0.f) {
505         float proj[3];
506         /* flatten along wave */
507         project_v3_v3v3(proj, par_vec, kink);
508         madd_v3_v3fl(result, proj, flat);
509
510         /* flatten along strand */
511         project_v3_v3v3(proj, par_vec, par_vel);
512         madd_v3_v3fl(result, proj, flat);
513       }
514       break;
515     }
516     case PART_KINK_BRAID: {
517       float y_vec[3] = {0.f, 1.f, 0.f};
518       float z_vec[3] = {0.f, 0.f, 1.f};
519       float vec_one[3], state_co[3];
520       float inp_y, inp_z, length;
521
522       if (par_rot) {
523         mul_qt_v3(par_rot, y_vec);
524         mul_qt_v3(par_rot, z_vec);
525       }
526
527       negate_v3(par_vec);
528       normalize_v3_v3(vec_one, par_vec);
529
530       inp_y = dot_v3v3(y_vec, vec_one);
531       inp_z = dot_v3v3(z_vec, vec_one);
532
533       if (inp_y > 0.5f) {
534         copy_v3_v3(state_co, y_vec);
535
536         mul_v3_fl(y_vec, amplitude * cosf(t));
537         mul_v3_fl(z_vec, amplitude / 2.f * sinf(2.f * t));
538       }
539       else if (inp_z > 0.0f) {
540         mul_v3_v3fl(state_co, z_vec, sinf((float)M_PI / 3.f));
541         madd_v3_v3fl(state_co, y_vec, -0.5f);
542
543         mul_v3_fl(y_vec, -amplitude * cosf(t + (float)M_PI / 3.f));
544         mul_v3_fl(z_vec, amplitude / 2.f * cosf(2.f * t + (float)M_PI / 6.f));
545       }
546       else {
547         mul_v3_v3fl(state_co, z_vec, -sinf((float)M_PI / 3.f));
548         madd_v3_v3fl(state_co, y_vec, -0.5f);
549
550         mul_v3_fl(y_vec, amplitude * -sinf(t + (float)M_PI / 6.f));
551         mul_v3_fl(z_vec, amplitude / 2.f * -sinf(2.f * t + (float)M_PI / 3.f));
552       }
553
554       mul_v3_fl(state_co, amplitude);
555       add_v3_v3(state_co, par_co);
556       sub_v3_v3v3(par_vec, state->co, state_co);
557
558       length = normalize_v3(par_vec);
559       mul_v3_fl(par_vec, MIN2(length, amplitude / 2.f));
560
561       add_v3_v3v3(state_co, par_co, y_vec);
562       add_v3_v3(state_co, z_vec);
563       add_v3_v3(state_co, par_vec);
564
565       shape = 2.f * (float)M_PI * (1.f + shape);
566
567       if (t < shape) {
568         shape = t / shape;
569         shape = (float)sqrt((double)shape);
570         interp_v3_v3v3(result, result, state_co, shape);
571       }
572       else {
573         copy_v3_v3(result, state_co);
574       }
575       break;
576     }
577   }
578
579   /* blend the start of the kink */
580   if (dt < 1.f) {
581     interp_v3_v3v3(state->co, state->co, result, dt);
582   }
583   else {
584     copy_v3_v3(state->co, result);
585   }
586 }
587
588 static float do_clump_level(float result[3],
589                             const float co[3],
590                             const float par_co[3],
591                             float time,
592                             float clumpfac,
593                             float clumppow,
594                             float pa_clump,
595                             CurveMapping *clumpcurve)
596 {
597   float clump = 0.0f;
598
599   if (clumpcurve) {
600     clump = pa_clump * (1.0f - clamp_f(curvemapping_evaluateF(clumpcurve, 0, time), 0.0f, 1.0f));
601
602     interp_v3_v3v3(result, co, par_co, clump);
603   }
604   else if (clumpfac != 0.0f) {
605     float cpow;
606
607     if (clumppow < 0.0f) {
608       cpow = 1.0f + clumppow;
609     }
610     else {
611       cpow = 1.0f + 9.0f * clumppow;
612     }
613
614     if (clumpfac < 0.0f) { /* clump roots instead of tips */
615       clump = -clumpfac * pa_clump * (float)pow(1.0 - (double)time, (double)cpow);
616     }
617     else {
618       clump = clumpfac * pa_clump * (float)pow((double)time, (double)cpow);
619     }
620
621     interp_v3_v3v3(result, co, par_co, clump);
622   }
623
624   return clump;
625 }
626
627 float do_clump(ParticleKey *state,
628                const float par_co[3],
629                float time,
630                const float orco_offset[3],
631                float clumpfac,
632                float clumppow,
633                float pa_clump,
634                bool use_clump_noise,
635                float clump_noise_size,
636                CurveMapping *clumpcurve)
637 {
638   float clump;
639
640   if (use_clump_noise && clump_noise_size != 0.0f) {
641     float center[3], noisevec[3];
642     float da[4], pa[12];
643
644     mul_v3_v3fl(noisevec, orco_offset, 1.0f / clump_noise_size);
645     voronoi(noisevec[0], noisevec[1], noisevec[2], da, pa, 1.0f, 0);
646     mul_v3_fl(&pa[0], clump_noise_size);
647     add_v3_v3v3(center, par_co, &pa[0]);
648
649     do_clump_level(state->co, state->co, center, time, clumpfac, clumppow, pa_clump, clumpcurve);
650   }
651
652   clump = do_clump_level(
653       state->co, state->co, par_co, time, clumpfac, clumppow, pa_clump, clumpcurve);
654
655   return clump;
656 }
657
658 static void do_rough(const float loc[3],
659                      float mat[4][4],
660                      float t,
661                      float fac,
662                      float size,
663                      float thres,
664                      ParticleKey *state)
665 {
666   float rough[3];
667   float rco[3];
668
669   if (thres != 0.0f) {
670     if (fabsf((float)(-1.5f + loc[0] + loc[1] + loc[2])) < 1.5f * thres) {
671       return;
672     }
673   }
674
675   copy_v3_v3(rco, loc);
676   mul_v3_fl(rco, t);
677   rough[0] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[0], rco[1], rco[2], 2, 0, 2);
678   rough[1] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[1], rco[2], rco[0], 2, 0, 2);
679   rough[2] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[2], rco[0], rco[1], 2, 0, 2);
680
681   madd_v3_v3fl(state->co, mat[0], fac * rough[0]);
682   madd_v3_v3fl(state->co, mat[1], fac * rough[1]);
683   madd_v3_v3fl(state->co, mat[2], fac * rough[2]);
684 }
685
686 static void do_rough_end(
687     const float loc[3], float mat[4][4], float t, float fac, float shape, ParticleKey *state)
688 {
689   float rough[2];
690   float roughfac;
691
692   roughfac = fac * (float)pow((double)t, shape);
693   copy_v2_v2(rough, loc);
694   rough[0] = -1.0f + 2.0f * rough[0];
695   rough[1] = -1.0f + 2.0f * rough[1];
696   mul_v2_fl(rough, roughfac);
697
698   madd_v3_v3fl(state->co, mat[0], rough[0]);
699   madd_v3_v3fl(state->co, mat[1], rough[1]);
700 }
701
702 static void do_rough_curve(const float loc[3],
703                            float mat[4][4],
704                            float time,
705                            float fac,
706                            float size,
707                            CurveMapping *roughcurve,
708                            ParticleKey *state)
709 {
710   float rough[3];
711   float rco[3];
712
713   if (!roughcurve) {
714     return;
715   }
716
717   fac *= clamp_f(curvemapping_evaluateF(roughcurve, 0, time), 0.0f, 1.0f);
718
719   copy_v3_v3(rco, loc);
720   mul_v3_fl(rco, time);
721   rough[0] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[0], rco[1], rco[2], 2, 0, 2);
722   rough[1] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[1], rco[2], rco[0], 2, 0, 2);
723   rough[2] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[2], rco[0], rco[1], 2, 0, 2);
724
725   madd_v3_v3fl(state->co, mat[0], fac * rough[0]);
726   madd_v3_v3fl(state->co, mat[1], fac * rough[1]);
727   madd_v3_v3fl(state->co, mat[2], fac * rough[2]);
728 }
729
730 static int twist_num_segments(const ParticleChildModifierContext *modifier_ctx)
731 {
732   ParticleThreadContext *thread_ctx = modifier_ctx->thread_ctx;
733   return (thread_ctx != NULL) ? thread_ctx->segments : modifier_ctx->sim->psys->part->draw_step;
734 }
735
736 static void twist_get_axis(const ParticleChildModifierContext *modifier_ctx,
737                            const float time,
738                            float r_axis[3])
739 {
740   const int num_segments = twist_num_segments(modifier_ctx);
741   const int index = clamp_i(time * num_segments, 0, num_segments);
742   if (index > 0) {
743     sub_v3_v3v3(
744         r_axis, modifier_ctx->parent_keys[index].co, modifier_ctx->parent_keys[index - 1].co);
745   }
746   else {
747     sub_v3_v3v3(
748         r_axis, modifier_ctx->parent_keys[index + 1].co, modifier_ctx->parent_keys[index].co);
749   }
750 }
751
752 static float curvemapping_integrate_clamped(CurveMapping *curve,
753                                             float start,
754                                             float end,
755                                             float step)
756 {
757   float integral = 0.0f;
758   float x = start;
759   while (x < end) {
760     float y = curvemapping_evaluateF(curve, 0, x);
761     y = clamp_f(y, 0.0f, 1.0f);
762     /* TODO(sergey): Clamp last step to end. */
763     integral += y * step;
764     x += step;
765   }
766   return integral;
767 }
768
769 static void do_twist(const ParticleChildModifierContext *modifier_ctx,
770                      ParticleKey *state,
771                      const float time)
772 {
773   ParticleThreadContext *thread_ctx = modifier_ctx->thread_ctx;
774   ParticleSimulationData *sim = modifier_ctx->sim;
775   ParticleTexture *ptex = modifier_ctx->ptex;
776   ParticleSettings *part = sim->psys->part;
777   /* Early output checks. */
778   if (modifier_ctx->parent_keys == NULL) {
779     /* Cannot get axis of rotation... */
780     return;
781   }
782   if (part->childtype != PART_CHILD_PARTICLES) {
783     /* Interpolated children behave weird with twist. */
784     return;
785   }
786   if (part->twist == 0.0f) {
787     /* No twist along the strand.  */
788     return;
789   }
790   /* Dependent on whether it's threaded update or not, curve comes
791    * from different places.
792    */
793   CurveMapping *twist_curve = NULL;
794   if (part->child_flag & PART_CHILD_USE_TWIST_CURVE) {
795     twist_curve = (thread_ctx != NULL) ? thread_ctx->twistcurve : part->twistcurve;
796   }
797   /* Axis of rotation. */
798   float axis[3];
799   twist_get_axis(modifier_ctx, time, axis);
800   /* Angle of rotation. */
801   float angle = part->twist;
802   if (ptex != NULL) {
803     angle *= (ptex->twist - 0.5f) * 2.0f;
804   }
805   if (twist_curve != NULL) {
806     const int num_segments = twist_num_segments(modifier_ctx);
807     angle *= curvemapping_integrate_clamped(twist_curve, 0.0f, time, 1.0f / num_segments);
808   }
809   else {
810     angle *= time;
811   }
812   /* Perform rotation around parent curve. */
813   float vec[3];
814   sub_v3_v3v3(vec, state->co, modifier_ctx->par_co);
815   rotate_v3_v3v3fl(state->co, vec, axis, angle * 2.0f * M_PI);
816   add_v3_v3(state->co, modifier_ctx->par_co);
817 }
818
819 void do_child_modifiers(const ParticleChildModifierContext *modifier_ctx,
820                         float mat[4][4],
821                         ParticleKey *state,
822                         float t)
823 {
824   ParticleThreadContext *ctx = modifier_ctx->thread_ctx;
825   ParticleSimulationData *sim = modifier_ctx->sim;
826   ParticleTexture *ptex = modifier_ctx->ptex;
827   ChildParticle *cpa = modifier_ctx->cpa;
828   ParticleSettings *part = sim->psys->part;
829   CurveMapping *clumpcurve = NULL, *roughcurve = NULL;
830   int i = cpa - sim->psys->child;
831   int guided = 0;
832
833   if (part->child_flag & PART_CHILD_USE_CLUMP_CURVE) {
834     clumpcurve = (ctx != NULL) ? ctx->clumpcurve : part->clumpcurve;
835   }
836   if (part->child_flag & PART_CHILD_USE_ROUGH_CURVE) {
837     roughcurve = (ctx != NULL) ? ctx->roughcurve : part->roughcurve;
838   }
839
840   float kink_amp = part->kink_amp;
841   float kink_amp_clump = part->kink_amp_clump;
842   float kink_freq = part->kink_freq;
843   float rough1 = part->rough1;
844   float rough2 = part->rough2;
845   float rough_end = part->rough_end;
846   const bool smooth_start = (sim->psys->part->childtype == PART_CHILD_FACES);
847
848   if (ptex) {
849     kink_amp *= ptex->kink_amp;
850     kink_freq *= ptex->kink_freq;
851     rough1 *= ptex->rough1;
852     rough2 *= ptex->rough2;
853     rough_end *= ptex->roughe;
854   }
855
856   do_twist(modifier_ctx, state, t);
857
858   if (part->flag & PART_CHILD_EFFECT) {
859     /* state is safe to cast, since only co and vel are used */
860     guided = do_guides(sim->depsgraph,
861                        sim->psys->part,
862                        sim->psys->effectors,
863                        (ParticleKey *)state,
864                        cpa->parent,
865                        t);
866   }
867
868   if (guided == 0) {
869     float orco_offset[3];
870     float clump;
871
872     sub_v3_v3v3(orco_offset, modifier_ctx->orco, modifier_ctx->par_orco);
873     clump = do_clump(state,
874                      modifier_ctx->par_co,
875                      t,
876                      orco_offset,
877                      part->clumpfac,
878                      part->clumppow,
879                      ptex ? ptex->clump : 1.0f,
880                      part->child_flag & PART_CHILD_USE_CLUMP_NOISE,
881                      part->clump_noise_size,
882                      clumpcurve);
883
884     if (kink_freq != 0.f) {
885       kink_amp *= (1.f - kink_amp_clump * clump);
886
887       do_kink(state,
888               modifier_ctx->par_co,
889               modifier_ctx->par_vel,
890               modifier_ctx->par_rot,
891               t,
892               kink_freq,
893               part->kink_shape,
894               kink_amp,
895               part->kink_flat,
896               part->kink,
897               part->kink_axis,
898               sim->ob->obmat,
899               smooth_start);
900     }
901   }
902
903   if (roughcurve) {
904     do_rough_curve(modifier_ctx->orco, mat, t, rough1, part->rough1_size, roughcurve, state);
905   }
906   else {
907     if (rough1 > 0.f) {
908       do_rough(modifier_ctx->orco, mat, t, rough1, part->rough1_size, 0.0, state);
909     }
910
911     if (rough2 > 0.f) {
912       float vec[3];
913       psys_frand_vec(sim->psys, i + 27, vec);
914       do_rough(vec, mat, t, rough2, part->rough2_size, part->rough2_thres, state);
915     }
916
917     if (rough_end > 0.f) {
918       float vec[3];
919       psys_frand_vec(sim->psys, i + 27, vec);
920       do_rough_end(vec, mat, t, rough_end, part->rough_end_shape, state);
921     }
922   }
923 }