← Back to branch summary

~ubuntu-branches/ubuntu/trusty/blender/trusty

« back to all changes in this revision

Viewing changes to source/blender/blenkernel/intern/ocean.c

  • Committer: Package Import Robot
  • Author(s): Jeremy Bicha
  • Date: 2013-03-06 12:08:47 UTC
  • mfrom: (1.5.1) (14.1.8 experimental)
  • Revision ID: package-import@ubuntu.com-20130306120847-frjfaryb2zrotwcg
Tags: 2.66a-1ubuntu1
https://launchpad.net/bugs/1076930
https://launchpad.net/bugs/1089256
https://launchpad.net/bugs/1052743
https://launchpad.net/bugs/999024
* Resynchronize with Debian (LP: #1076930, #1089256, #1052743, #999024,
  #1122888, #1147084)
* debian/control:
  - Lower build-depends on libavcodec-dev since we're not
    doing the libav9 transition in Ubuntu yet

Show diffs side-by-side

added added

removed removed

Lines of Context:
source/blender/blenkernel/intern/ocean.c
35
35
 
36
36
#include "DNA_scene_types.h"
37
37
 
38
 
#include "BKE_image.h"
39
 
#include "BKE_ocean.h"
40
 
#include "BKE_utildefines.h"
41
 
 
42
 
#include "BKE_global.h" // XXX TESTING
43
 
 
44
 
#include "BLI_math_base.h"
45
 
#include "BLI_math_inline.h"
 
38
#include "BLI_math.h"
 
39
#include "BLI_path_util.h"
46
40
#include "BLI_rand.h"
47
41
#include "BLI_string.h"
48
42
#include "BLI_threads.h"
49
 
#include "BLI_path_util.h"
50
43
#include "BLI_utildefines.h"
51
44
 
 
45
#include "BKE_image.h"
 
46
#include "BKE_ocean.h"
 
47
 
52
48
#include "IMB_imbuf.h"
53
49
#include "IMB_imbuf_types.h"
54
50
 
56
52
 
57
53
#ifdef WITH_OCEANSIM
58
54
 
59
 
// Ocean code
 
55
/* Ocean code */
60
56
#include "fftw3.h"
61
57
 
62
58
#define GRAVITY  9.81f
84
80
        float _Lx;
85
81
        float _Lz;
86
82
 
87
 
        float normalize_factor;                                 // init w
 
83
        float normalize_factor;                 /* init w */
88
84
        float time;
89
85
 
90
86
        short _do_disp_y;
98
94
        /* ********* sim data arrays ********* */
99
95
 
100
96
        /* two dimensional arrays of complex */
101
 
        fftw_complex *_fft_in;                  // init w       sim w
102
 
        fftw_complex *_fft_in_x;                        // init w       sim w
103
 
        fftw_complex *_fft_in_z;                        // init w       sim w
104
 
        fftw_complex *_fft_in_jxx;                      // init w       sim w
105
 
        fftw_complex *_fft_in_jzz;                      // init w       sim w
106
 
        fftw_complex *_fft_in_jxz;                      // init w       sim w
107
 
        fftw_complex *_fft_in_nx;                       // init w       sim w
108
 
        fftw_complex *_fft_in_nz;                       // init w       sim w
109
 
        fftw_complex *_htilda;                  // init w       sim w (only once)
 
97
        fftw_complex *_fft_in;          /* init w       sim w */
 
98
        fftw_complex *_fft_in_x;        /* init w       sim w */
 
99
        fftw_complex *_fft_in_z;        /* init w       sim w */
 
100
        fftw_complex *_fft_in_jxx;      /* init w       sim w */
 
101
        fftw_complex *_fft_in_jzz;      /* init w       sim w */
 
102
        fftw_complex *_fft_in_jxz;      /* init w       sim w */
 
103
        fftw_complex *_fft_in_nx;       /* init w       sim w */
 
104
        fftw_complex *_fft_in_nz;       /* init w       sim w */
 
105
        fftw_complex *_htilda;          /* init w       sim w (only once) */
110
106
 
111
107
        /* fftw "plans" */
112
 
        fftw_plan _disp_y_plan;                 // init w       sim r
113
 
        fftw_plan _disp_x_plan;                 // init w       sim r
114
 
        fftw_plan _disp_z_plan;                 // init w       sim r
115
 
        fftw_plan _N_x_plan;                    // init w       sim r
116
 
        fftw_plan _N_z_plan;                    // init w       sim r
117
 
        fftw_plan _Jxx_plan;                    // init w       sim r
118
 
        fftw_plan _Jxz_plan;                    // init w       sim r
119
 
        fftw_plan _Jzz_plan;                    // init w       sim r
 
108
        fftw_plan _disp_y_plan;         /* init w       sim r */
 
109
        fftw_plan _disp_x_plan;         /* init w       sim r */
 
110
        fftw_plan _disp_z_plan;         /* init w       sim r */
 
111
        fftw_plan _N_x_plan;            /* init w       sim r */
 
112
        fftw_plan _N_z_plan;            /* init w       sim r */
 
113
        fftw_plan _Jxx_plan;            /* init w       sim r */
 
114
        fftw_plan _Jxz_plan;            /* init w       sim r */
 
115
        fftw_plan _Jzz_plan;            /* init w       sim r */
120
116
 
121
117
        /* two dimensional arrays of float */
122
 
        double *  _disp_y;                              // init w       sim w via plan?
123
 
        double * _N_x;                                  // init w       sim w via plan?
124
 
        /*float * _N_y; all member of this array has same values, so convert this array to a float to reduce memory usage (MEM01)*/
125
 
        double _N_y;                                    //                      sim w ********* can be rearranged?
126
 
        double * _N_z;                                  // init w       sim w via plan?
127
 
        double * _disp_x;                               // init w       sim w via plan?
128
 
        double * _disp_z;                               // init w       sim w via plan?
 
118
        double *_disp_y;                /* init w       sim w via plan? */
 
119
        double *_N_x;                   /* init w       sim w via plan? */
 
120
        /* all member of this array has same values, so convert this array to a float to reduce memory usage (MEM01)*/
 
121
        /*float * _N_y; */
 
122
        double _N_y;                    /*                      sim w ********* can be rearranged? */
 
123
        double *_N_z;                   /* init w       sim w via plan? */
 
124
        double *_disp_x;                /* init w       sim w via plan? */
 
125
        double *_disp_z;                /* init w       sim w via plan? */
129
126
 
130
127
        /* two dimensional arrays of float */
131
128
        /* Jacobian and minimum eigenvalue */
132
 
        double * _Jxx;                                  // init w       sim w
133
 
        double * _Jzz;                                  // init w       sim w
134
 
        double * _Jxz;                                  // init w       sim w
 
129
        double *_Jxx;                   /* init w       sim w */
 
130
        double *_Jzz;                   /* init w       sim w */
 
131
        double *_Jxz;                   /* init w       sim w */
135
132
 
136
133
        /* one dimensional float array */
137
 
        float * _kx;                                    // init w       sim r
138
 
        float * _kz;                                    // init w       sim r
 
134
        float *_kx;                     /* init w       sim r */
 
135
        float *_kz;                     /* init w       sim r */
139
136
 
140
137
        /* two dimensional complex array */
141
 
        fftw_complex * _h0;                             // init w       sim r
142
 
        fftw_complex * _h0_minus;               // init w       sim r
 
138
        fftw_complex *_h0;              /* init w       sim r */
 
139
        fftw_complex *_h0_minus;        /* init w       sim r */
143
140
 
144
141
        /* two dimensional float array */
145
 
        float * _k;                                             // init w       sim r
 
142
        float *_k;                      /* init w       sim r */
146
143
} Ocean;
147
144
 
148
145
 
149
146
 
150
147
static float nextfr(float min, float max)
151
148
{
152
 
        return BLI_frand()*(min-max)+max;
 
149
        return BLI_frand() * (min - max) + max;
153
150
}
154
151
 
155
 
static float gaussRand (void)
 
152
static float gaussRand(void)
156
153
{
157
 
        float x;                // Note: to avoid numerical problems with very small
158
 
        float y;                // numbers, we make these variables singe-precision
159
 
        float length2;  // floats, but later we call the double-precision log()
160
 
        // and sqrt() functions instead of logf() and sqrtf().
161
 
        do
162
 
        {
163
 
                x = (float) (nextfr (-1, 1));
164
 
                y = (float)(nextfr (-1, 1));
 
154
        /* Note: to avoid numerical problems with very small numbers, we make these variables singe-precision floats,
 
155
         * but later we call the double-precision log() and sqrt() functions instead of logf() and sqrtf().
 
156
         */ 
 
157
        float x;
 
158
        float y;
 
159
        float length2;
 
160
 
 
161
        do {
 
162
                x = (float) (nextfr(-1, 1));
 
163
                y = (float)(nextfr(-1, 1));
165
164
                length2 = x * x + y * y;
166
 
        }
167
 
        while (length2 >= 1 || length2 == 0);
 
165
        } while (length2 >= 1 || length2 == 0);
168
166
 
169
167
        return x * sqrtf(-2.0f * logf(length2) / length2);
170
168
}
171
169
 
172
170
/**
173
171
 * Some useful functions
174
 
 * */
175
 
MINLINE float lerp(float a,float b,float f)
176
 
{
177
 
        return a + (b-a)*f;
178
 
}
179
 
 
180
 
MINLINE float catrom(float p0,float p1,float p2,float p3,float f)
181
 
{
182
 
        return 0.5f *((2.0f * p1) +
183
 
                      (-p0 + p2) * f +
184
 
                      (2.0f*p0 - 5.0f*p1 + 4.0f*p2 - p3) * f*f +
185
 
                      (-p0 + 3.0f*p1- 3.0f*p2 + p3) * f*f*f);
 
172
 */
 
173
MINLINE float catrom(float p0, float p1, float p2, float p3, float f)
 
174
{
 
175
        return 0.5f * ((2.0f * p1) +
 
176
                       (-p0 + p2) * f +
 
177
                       (2.0f * p0 - 5.0f * p1 + 4.0f * p2 - p3) * f * f +
 
178
                       (-p0 + 3.0f * p1 - 3.0f * p2 + p3) * f * f * f);
186
179
}
187
180
 
188
181
MINLINE float omega(float k, float depth)
189
182
{
190
 
        return sqrt(GRAVITY*k * tanh(k*depth));
 
183
        return sqrtf(GRAVITY * k * tanhf(k * depth));
191
184
}
192
185
 
193
 
// modified Phillips spectrum
194
 
static float Ph(struct Ocean* o, float kx,float kz )
 
186
/* modified Phillips spectrum */
 
187
static float Ph(struct Ocean *o, float kx, float kz)
195
188
{
196
189
        float tmp;
197
 
        float k2 = kx*kx + kz*kz;
 
190
        float k2 = kx * kx + kz * kz;
198
191
 
199
 
        if (k2 == 0.0f)
200
 
        {
201
 
                return 0.0f; // no DC component
 
192
        if (k2 == 0.0f) {
 
193
                return 0.0f; /* no DC component */
202
194
        }
203
195
 
204
 
        // damp out the waves going in the direction opposite the wind
205
 
        tmp = (o->_wx * kx  + o->_wz * kz)/sqrtf(k2);
206
 
        if (tmp < 0)
207
 
        {
 
196
        /* damp out the waves going in the direction opposite the wind */
 
197
        tmp = (o->_wx * kx + o->_wz * kz) / sqrtf(k2);
 
198
        if (tmp < 0) {
208
199
                tmp *= o->_damp_reflections;
209
200
        }
210
201
 
211
 
        return o->_A * expf( -1.0f / (k2*(o->_L*o->_L))) * expf(-k2 * (o->_l*o->_l)) * powf(fabsf(tmp),o->_wind_alignment) / (k2*k2);
 
202
        return o->_A * expf(-1.0f / (k2 * (o->_L * o->_L))) * expf(-k2 * (o->_l * o->_l)) *
 
203
               powf(fabsf(tmp), o->_wind_alignment) / (k2 * k2);
212
204
}
213
205
 
214
 
static void compute_eigenstuff(struct OceanResult *ocr, float jxx,float jzz,float jxz)
 
206
static void compute_eigenstuff(struct OceanResult *ocr, float jxx, float jzz, float jxz)
215
207
{
216
 
        float a,b,qplus,qminus;
 
208
        float a, b, qplus, qminus;
217
209
        a = jxx + jzz;
218
 
        b = sqrt((jxx - jzz)*(jxx - jzz) + 4 * jxz * jxz);
219
 
 
220
 
        ocr->Jminus = 0.5f*(a-b);
221
 
        ocr->Jplus  = 0.5f*(a+b);
222
 
 
223
 
        qplus  = (ocr->Jplus  - jxx)/jxz;
224
 
        qminus = (ocr->Jminus - jxx)/jxz;
225
 
 
226
 
        a = sqrt(1 + qplus*qplus);
227
 
        b = sqrt(1 + qminus*qminus);
228
 
 
229
 
        ocr->Eplus[0] = 1.0f/ a;
 
210
        b = sqrt((jxx - jzz) * (jxx - jzz) + 4 * jxz * jxz);
 
211
 
 
212
        ocr->Jminus = 0.5f * (a - b);
 
213
        ocr->Jplus  = 0.5f * (a + b);
 
214
 
 
215
        qplus  = (ocr->Jplus  - jxx) / jxz;
 
216
        qminus = (ocr->Jminus - jxx) / jxz;
 
217
 
 
218
        a = sqrt(1 + qplus * qplus);
 
219
        b = sqrt(1 + qminus * qminus);
 
220
 
 
221
        ocr->Eplus[0] = 1.0f / a;
230
222
        ocr->Eplus[1] = 0.0f;
231
 
        ocr->Eplus[2] = qplus/a;
 
223
        ocr->Eplus[2] = qplus / a;
232
224
 
233
 
        ocr->Eminus[0] = 1.0f/b;
 
225
        ocr->Eminus[0] = 1.0f / b;
234
226
        ocr->Eminus[1] = 0.0f;
235
 
        ocr->Eminus[2] = qminus/b;
 
227
        ocr->Eminus[2] = qminus / b;
236
228
}
237
229
 
238
230
/*
246
238
        cmpl[1] = image;
247
239
}
248
240
 
249
 
#if 0   // unused
 
241
#if 0   /* unused */
250
242
static void add_complex_f(fftw_complex res, fftw_complex cmpl, float f)
251
243
{
252
244
        res[0] = cmpl[0] + f;
262
254
 
263
255
static void mul_complex_f(fftw_complex res, fftw_complex cmpl, float f)
264
256
{
265
 
        res[0] = cmpl[0]*f;
266
 
        res[1] = cmpl[1]*f;
 
257
        res[0] = cmpl[0] * (double)f;
 
258
        res[1] = cmpl[1] * (double)f;
267
259
}
268
260
 
269
261
static void mul_complex_c(fftw_complex res, fftw_complex cmpl1, fftw_complex cmpl2)
270
262
{
271
263
        fftwf_complex temp;
272
 
        temp[0] = cmpl1[0]*cmpl2[0]-cmpl1[1]*cmpl2[1];
273
 
        temp[1] = cmpl1[0]*cmpl2[1]+cmpl1[1]*cmpl2[0];
 
264
        temp[0] = cmpl1[0] * cmpl2[0] - cmpl1[1] * cmpl2[1];
 
265
        temp[1] = cmpl1[0] * cmpl2[1] + cmpl1[1] * cmpl2[0];
274
266
        res[0] = temp[0];
275
267
        res[1] = temp[1];
276
268
}
295
287
{
296
288
        float r = expf(cmpl[0]);
297
289
 
298
 
        res[0] = cos(cmpl[1])*r;
299
 
        res[1] = sin(cmpl[1])*r;
 
290
        res[0] = cosf(cmpl[1]) * r;
 
291
        res[1] = sinf(cmpl[1]) * r;
300
292
}
301
293
 
302
294
float BKE_ocean_jminus_to_foam(float jminus, float coverage)
303
295
{
304
296
        float foam = jminus * -0.005f + coverage;
305
297
        CLAMP(foam, 0.0f, 1.0f);
306
 
        return foam*foam;
 
298
        return foam * foam;
307
299
}
308
300
 
309
 
void BKE_ocean_eval_uv(struct Ocean *oc, struct OceanResult *ocr, float u,float v)
 
301
void BKE_ocean_eval_uv(struct Ocean *oc, struct OceanResult *ocr, float u, float v)
310
302
{
311
 
        int i0,i1,j0,j1;
312
 
        float frac_x,frac_z;
313
 
        float uu,vv;
 
303
        int i0, i1, j0, j1;
 
304
        float frac_x, frac_z;
 
305
        float uu, vv;
314
306
 
315
 
        // first wrap the texture so 0 <= (u,v) < 1
316
 
        u = fmodf(u,1.0f);
317
 
        v = fmodf(v,1.0f);
 
307
        /* first wrap the texture so 0 <= (u, v) < 1 */
 
308
        u = fmodf(u, 1.0f);
 
309
        v = fmodf(v, 1.0f);
318
310
 
319
311
        if (u < 0) u += 1.0f;
320
312
        if (v < 0) v += 1.0f;
340
332
        j1 = j1 % oc->_N;
341
333
 
342
334
 
343
 
#define BILERP(m) (lerp(lerp(m[i0*oc->_N+j0],m[i1*oc->_N+j0],frac_x),lerp(m[i0*oc->_N+j1],m[i1*oc->_N+j1],frac_x),frac_z))
 
335
#define BILERP(m) (interpf(interpf(m[i1 * oc->_N + j1], m[i0 * oc->_N + j1], frac_x), \
 
336
                           interpf(m[i1 * oc->_N + j0], m[i0 * oc->_N + j0], frac_x), \
 
337
                           frac_z))
344
338
        {
345
339
                if (oc->_do_disp_y) {
346
340
                        ocr->disp[1] = BILERP(oc->_disp_y);
348
342
 
349
343
                if (oc->_do_normals) {
350
344
                        ocr->normal[0] = BILERP(oc->_N_x);
351
 
                        ocr->normal[1] = oc->_N_y/*BILERP(oc->_N_y) (MEM01)*/;
 
345
                        ocr->normal[1] = oc->_N_y /*BILERP(oc->_N_y) (MEM01)*/;
352
346
                        ocr->normal[2] = BILERP(oc->_N_z);
353
347
                }
354
348
 
362
356
                }
363
357
 
364
358
                if (oc->_do_jacobian) {
365
 
                        compute_eigenstuff(ocr, BILERP(oc->_Jxx),BILERP(oc->_Jzz),BILERP(oc->_Jxz));
 
359
                        compute_eigenstuff(ocr, BILERP(oc->_Jxx), BILERP(oc->_Jzz), BILERP(oc->_Jxz));
366
360
                }
367
361
        }
368
362
#undef BILERP
370
364
        BLI_rw_mutex_unlock(&oc->oceanmutex);
371
365
}
372
366
 
373
 
// use catmullrom interpolation rather than linear
374
 
void BKE_ocean_eval_uv_catrom(struct Ocean *oc, struct OceanResult *ocr, float u,float v)
 
367
/* use catmullrom interpolation rather than linear */
 
368
void BKE_ocean_eval_uv_catrom(struct Ocean *oc, struct OceanResult *ocr, float u, float v)
375
369
{
376
 
        int i0,i1,i2,i3,j0,j1,j2,j3;
377
 
        float frac_x,frac_z;
378
 
        float uu,vv;
 
370
        int i0, i1, i2, i3, j0, j1, j2, j3;
 
371
        float frac_x, frac_z;
 
372
        float uu, vv;
379
373
 
380
 
        // first wrap the texture so 0 <= (u,v) < 1
381
 
        u = fmod(u,1.0f);
382
 
        v = fmod(v,1.0f);
 
374
        /* first wrap the texture so 0 <= (u, v) < 1 */
 
375
        u = fmod(u, 1.0f);
 
376
        v = fmod(v, 1.0f);
383
377
 
384
378
        if (u < 0) u += 1.0f;
385
379
        if (v < 0) v += 1.0f;
404
398
        i2 = i2 % oc->_M;
405
399
        j2 = j2 % oc->_N;
406
400
 
407
 
        i0 = (i1-1);
408
 
        i3 = (i2+1);
 
401
        i0 = (i1 - 1);
 
402
        i3 = (i2 + 1);
409
403
        i0 = i0 <   0 ? i0 + oc->_M : i0;
410
404
        i3 = i3 >= oc->_M ? i3 - oc->_M : i3;
411
405
 
412
 
        j0 = (j1-1);
413
 
        j3 = (j2+1);
 
406
        j0 = (j1 - 1);
 
407
        j3 = (j2 + 1);
414
408
        j0 = j0 <   0 ? j0 + oc->_N : j0;
415
409
        j3 = j3 >= oc->_N ? j3 - oc->_N : j3;
416
410
 
417
 
#define INTERP(m) catrom(catrom(m[i0*oc->_N+j0],m[i1*oc->_N+j0],m[i2*oc->_N+j0],m[i3*oc->_N+j0],frac_x),\
418
 
        catrom(m[i0*oc->_N+j1],m[i1*oc->_N+j1],m[i2*oc->_N+j1],m[i3*oc->_N+j1],frac_x),\
419
 
        catrom(m[i0*oc->_N+j2],m[i1*oc->_N+j2],m[i2*oc->_N+j2],m[i3*oc->_N+j2],frac_x),\
420
 
        catrom(m[i0*oc->_N+j3],m[i1*oc->_N+j3],m[i2*oc->_N+j3],m[i3*oc->_N+j3],frac_x),\
421
 
        frac_z)
 
411
#define INTERP(m) catrom(catrom(m[i0 * oc->_N + j0], m[i1 * oc->_N + j0], \
 
412
                                m[i2 * oc->_N + j0], m[i3 * oc->_N + j0], frac_x), \
 
413
                         catrom(m[i0 * oc->_N + j1], m[i1 * oc->_N + j1], \
 
414
                                m[i2 * oc->_N + j1], m[i3 * oc->_N + j1], frac_x), \
 
415
                         catrom(m[i0 * oc->_N + j2], m[i1 * oc->_N + j2], \
 
416
                                m[i2 * oc->_N + j2], m[i3 * oc->_N + j2], frac_x), \
 
417
                         catrom(m[i0 * oc->_N + j3], m[i1 * oc->_N + j3], \
 
418
                                m[i2 * oc->_N + j3], m[i3 * oc->_N + j3], frac_x), \
 
419
                         frac_z)
422
420
 
423
421
        {
424
 
                if (oc->_do_disp_y)
425
 
                {
 
422
                if (oc->_do_disp_y) {
426
423
                        ocr->disp[1] = INTERP(oc->_disp_y);
427
424
                }
428
 
                if (oc->_do_normals)
429
 
                {
 
425
                if (oc->_do_normals) {
430
426
                        ocr->normal[0] = INTERP(oc->_N_x);
431
 
                        ocr->normal[1] = oc->_N_y/*INTERP(oc->_N_y) (MEM01)*/;
 
427
                        ocr->normal[1] = oc->_N_y /*INTERP(oc->_N_y) (MEM01)*/;
432
428
                        ocr->normal[2] = INTERP(oc->_N_z);
433
429
                }
434
 
                if (oc->_do_chop)
435
 
                {
 
430
                if (oc->_do_chop) {
436
431
                        ocr->disp[0] = INTERP(oc->_disp_x);
437
432
                        ocr->disp[2] = INTERP(oc->_disp_z);
438
433
                }
441
436
                        ocr->disp[2] = 0.0;
442
437
                }
443
438
 
444
 
                if (oc->_do_jacobian)
445
 
                {
446
 
                        compute_eigenstuff(ocr, INTERP(oc->_Jxx),INTERP(oc->_Jzz),INTERP(oc->_Jxz));
 
439
                if (oc->_do_jacobian) {
 
440
                        compute_eigenstuff(ocr, INTERP(oc->_Jxx), INTERP(oc->_Jzz), INTERP(oc->_Jxz));
447
441
                }
448
442
        }
449
443
#undef INTERP
452
446
 
453
447
}
454
448
 
455
 
void BKE_ocean_eval_xz(struct Ocean *oc, struct OceanResult *ocr, float x,float z)
456
 
{
457
 
        BKE_ocean_eval_uv(oc, ocr, x/oc->_Lx,z/oc->_Lz);
458
 
}
459
 
 
460
 
void BKE_ocean_eval_xz_catrom(struct Ocean *oc, struct OceanResult *ocr, float x,float z)
461
 
{
462
 
        BKE_ocean_eval_uv_catrom(oc, ocr, x/oc->_Lx,z/oc->_Lz);
463
 
}
464
 
 
465
 
// note that this doesn't wrap properly for i,j < 0, but its
466
 
// not really meant for that being just a way to get the raw data out
467
 
// to save in some image format.
468
 
void BKE_ocean_eval_ij(struct Ocean *oc, struct OceanResult *ocr, int i,int j)
 
449
void BKE_ocean_eval_xz(struct Ocean *oc, struct OceanResult *ocr, float x, float z)
 
450
{
 
451
        BKE_ocean_eval_uv(oc, ocr, x / oc->_Lx, z / oc->_Lz);
 
452
}
 
453
 
 
454
void BKE_ocean_eval_xz_catrom(struct Ocean *oc, struct OceanResult *ocr, float x, float z)
 
455
{
 
456
        BKE_ocean_eval_uv_catrom(oc, ocr, x / oc->_Lx, z / oc->_Lz);
 
457
}
 
458
 
 
459
/* note that this doesn't wrap properly for i, j < 0, but its not really meant for that being just a way to get
 
460
 * the raw data out to save in some image format.
 
461
 */
 
462
void BKE_ocean_eval_ij(struct Ocean *oc, struct OceanResult *ocr, int i, int j)
469
463
{
470
464
        BLI_rw_mutex_lock(&oc->oceanmutex, THREAD_LOCK_READ);
471
465
 
472
466
        i = abs(i) % oc->_M;
473
467
        j = abs(j) % oc->_N;
474
468
 
475
 
        ocr->disp[1] = oc->_do_disp_y ? oc->_disp_y[i*oc->_N+j] : 0.0f;
 
469
        ocr->disp[1] = oc->_do_disp_y ? (float)oc->_disp_y[i * oc->_N + j] : 0.0f;
476
470
 
477
 
        if (oc->_do_chop)
478
 
        {
479
 
                ocr->disp[0] = oc->_disp_x[i*oc->_N+j];
480
 
                ocr->disp[2] = oc->_disp_z[i*oc->_N+j];
 
471
        if (oc->_do_chop) {
 
472
                ocr->disp[0] = oc->_disp_x[i * oc->_N + j];
 
473
                ocr->disp[2] = oc->_disp_z[i * oc->_N + j];
481
474
        }
482
475
        else {
483
476
                ocr->disp[0] = 0.0f;
484
477
                ocr->disp[2] = 0.0f;
485
478
        }
486
479
 
487
 
        if (oc->_do_normals)
488
 
        {
489
 
                ocr->normal[0] = oc->_N_x[i*oc->_N+j];
490
 
                ocr->normal[1] = oc->_N_y/*oc->_N_y[i*oc->_N+j] (MEM01)*/;
491
 
                ocr->normal[2] = oc->_N_z[i*oc->_N+j];
 
480
        if (oc->_do_normals) {
 
481
                ocr->normal[0] = oc->_N_x[i * oc->_N + j];
 
482
                ocr->normal[1] = oc->_N_y  /* oc->_N_y[i * oc->_N + j] (MEM01) */;
 
483
                ocr->normal[2] = oc->_N_z[i * oc->_N + j];
492
484
 
493
485
                normalize_v3(ocr->normal);
494
486
        }
495
487
 
496
 
        if (oc->_do_jacobian)
497
 
        {
498
 
                compute_eigenstuff(ocr, oc->_Jxx[i*oc->_N+j],oc->_Jzz[i*oc->_N+j],oc->_Jxz[i*oc->_N+j]);
 
488
        if (oc->_do_jacobian) {
 
489
                compute_eigenstuff(ocr, oc->_Jxx[i * oc->_N + j], oc->_Jzz[i * oc->_N + j], oc->_Jxz[i * oc->_N + j]);
499
490
        }
500
491
 
501
492
        BLI_rw_mutex_unlock(&oc->oceanmutex);
509
500
 
510
501
        BLI_rw_mutex_lock(&o->oceanmutex, THREAD_LOCK_WRITE);
511
502
 
512
 
        // compute a new htilda
 
503
        /* compute a new htilda */
513
504
#pragma omp parallel for private(i, j)
514
 
        for (i = 0 ; i  < o->_M ; ++i)
515
 
        {
516
 
                // note the <= _N/2 here, see the fftw doco about
517
 
                // the mechanics of the complex->real fft storage
518
 
                for ( j  = 0 ; j  <= o->_N / 2 ; ++j)
519
 
                {
 
505
        for (i = 0; i < o->_M; ++i) {
 
506
                /* note the <= _N/2 here, see the fftw doco about the mechanics of the complex->real fft storage */
 
507
                for (j = 0; j <= o->_N / 2; ++j) {
520
508
                        fftw_complex exp_param1;
521
509
                        fftw_complex exp_param2;
522
510
                        fftw_complex conj_param;
523
511
 
524
512
 
525
 
                        init_complex(exp_param1, 0.0, omega(o->_k[i*(1+o->_N/2)+j],o->_depth)*t);
526
 
                        init_complex(exp_param2, 0.0, -omega(o->_k[i*(1+o->_N/2)+j],o->_depth)*t);
 
513
                        init_complex(exp_param1, 0.0, omega(o->_k[i * (1 + o->_N / 2) + j], o->_depth) * t);
 
514
                        init_complex(exp_param2, 0.0, -omega(o->_k[i * (1 + o->_N / 2) + j], o->_depth) * t);
527
515
                        exp_complex(exp_param1, exp_param1);
528
516
                        exp_complex(exp_param2, exp_param2);
529
 
                        conj_complex(conj_param, o->_h0_minus[i*o->_N+j]);
 
517
                        conj_complex(conj_param, o->_h0_minus[i * o->_N + j]);
530
518
 
531
 
                        mul_complex_c(exp_param1, o->_h0[i*o->_N+j], exp_param1);
 
519
                        mul_complex_c(exp_param1, o->_h0[i * o->_N + j], exp_param1);
532
520
                        mul_complex_c(exp_param2, conj_param, exp_param2);
533
521
 
534
 
                        add_comlex_c(o->_htilda[i*(1+o->_N/2)+j], exp_param1, exp_param2);
535
 
                        mul_complex_f(o->_fft_in[i*(1+o->_N/2)+j], o->_htilda[i*(1+o->_N/2)+j], scale);
 
522
                        add_comlex_c(o->_htilda[i * (1 + o->_N / 2) + j], exp_param1, exp_param2);
 
523
                        mul_complex_f(o->_fft_in[i * (1 + o->_N / 2) + j], o->_htilda[i * (1 + o->_N / 2) + j], scale);
536
524
                }
537
525
        }
538
526
 
541
529
 
542
530
#pragma omp section
543
531
                {
544
 
                        if (o->_do_disp_y)
545
 
                        {
546
 
                                // y displacement
 
532
                        if (o->_do_disp_y) {
 
533
                                /* y displacement */
547
534
                                fftw_execute(o->_disp_y_plan);
548
535
                        }
549
 
                } // section 1
 
536
                } /* section 1 */
550
537
 
551
538
#pragma omp section
552
539
                {
553
 
                        if (o->_do_chop)
554
 
                        {
555
 
                                // x displacement
556
 
                                for ( i = 0 ; i  < o->_M ; ++i)
557
 
                                {
558
 
                                        for ( j  = 0 ; j  <= o->_N / 2 ; ++j)
559
 
                                        {
 
540
                        if (o->_do_chop) {
 
541
                                /* x displacement */
 
542
                                for (i = 0; i < o->_M; ++i) {
 
543
                                        for (j = 0; j <= o->_N / 2; ++j) {
560
544
                                                fftw_complex mul_param;
561
545
                                                fftw_complex minus_i;
562
546
 
564
548
                                                init_complex(mul_param, -scale, 0);
565
549
                                                mul_complex_f(mul_param, mul_param, chop_amount);
566
550
                                                mul_complex_c(mul_param, mul_param, minus_i);
567
 
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i*(1+o->_N/2)+j]);
568
 
                                                mul_complex_f(mul_param, mul_param, (o->_k[i*(1+o->_N/2)+j] == 0.0 ? 0.0 : o->_kx[i] / o->_k[i*(1+o->_N/2)+j]));
569
 
                                                init_complex(o->_fft_in_x[i*(1+o->_N/2)+j], real_c(mul_param), image_c(mul_param));
 
551
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i * (1 + o->_N / 2) + j]);
 
552
                                                mul_complex_f(mul_param, mul_param,
 
553
                                                              ((o->_k[i * (1 + o->_N / 2) + j] == 0.0f) ?
 
554
                                                               0.0f :
 
555
                                                               o->_kx[i] / o->_k[i * (1 + o->_N / 2) + j]));
 
556
                                                init_complex(o->_fft_in_x[i * (1 + o->_N / 2) + j], real_c(mul_param), image_c(mul_param));
570
557
                                        }
571
558
                                }
572
559
                                fftw_execute(o->_disp_x_plan);
573
560
                        }
574
 
                } //section 2
 
561
                } /* section 2 */
575
562
 
576
563
#pragma omp section
577
564
                {
578
 
                        if (o->_do_chop)
579
 
                        {
580
 
                                // z displacement
581
 
                                for ( i = 0 ; i  < o->_M ; ++i)
582
 
                                {
583
 
                                        for ( j  = 0 ; j  <= o->_N / 2 ; ++j)
584
 
                                        {
 
565
                        if (o->_do_chop) {
 
566
                                /* z displacement */
 
567
                                for (i = 0; i < o->_M; ++i) {
 
568
                                        for (j = 0; j <= o->_N / 2; ++j) {
585
569
                                                fftw_complex mul_param;
586
570
                                                fftw_complex minus_i;
587
571
 
589
573
                                                init_complex(mul_param, -scale, 0);
590
574
                                                mul_complex_f(mul_param, mul_param, chop_amount);
591
575
                                                mul_complex_c(mul_param, mul_param, minus_i);
592
 
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i*(1+o->_N/2)+j]);
593
 
                                                mul_complex_f(mul_param, mul_param, (o->_k[i*(1+o->_N/2)+j] == 0.0 ? 0.0 : o->_kz[j] / o->_k[i*(1+o->_N/2)+j]));
594
 
                                                init_complex(o->_fft_in_z[i*(1+o->_N/2)+j], real_c(mul_param), image_c(mul_param));
 
576
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i * (1 + o->_N / 2) + j]);
 
577
                                                mul_complex_f(mul_param, mul_param,
 
578
                                                              ((o->_k[i * (1 + o->_N / 2) + j] == 0.0f) ?
 
579
                                                               0.0f :
 
580
                                                               o->_kz[j] / o->_k[i * (1 + o->_N / 2) + j]));
 
581
                                                init_complex(o->_fft_in_z[i * (1 + o->_N / 2) + j], real_c(mul_param), image_c(mul_param));
595
582
                                        }
596
583
                                }
597
584
                                fftw_execute(o->_disp_z_plan);
598
585
                        }
599
 
                } // section 3
 
586
                } /* section 3 */
600
587
 
601
588
#pragma omp section
602
589
                {
603
 
                        if (o->_do_jacobian)
604
 
                        {
605
 
                                // Jxx
606
 
                                for ( i = 0 ; i  < o->_M ; ++i)
607
 
                                {
608
 
                                        for ( j  = 0 ; j  <= o->_N / 2 ; ++j)
609
 
                                        {
 
590
                        if (o->_do_jacobian) {
 
591
                                /* Jxx */
 
592
                                for (i = 0; i < o->_M; ++i) {
 
593
                                        for (j = 0; j <= o->_N / 2; ++j) {
610
594
                                                fftw_complex mul_param;
611
595
 
612
 
                                                //init_complex(mul_param, -scale, 0);
 
596
                                                /* init_complex(mul_param, -scale, 0); */
613
597
                                                init_complex(mul_param, -1, 0);
614
598
 
615
599
                                                mul_complex_f(mul_param, mul_param, chop_amount);
616
 
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i*(1+o->_N/2)+j]);
617
 
                                                mul_complex_f(mul_param, mul_param, (o->_k[i*(1+o->_N/2)+j] == 0.0 ? 0.0 : o->_kx[i]*o->_kx[i] / o->_k[i*(1+o->_N/2)+j]));
618
 
                                                init_complex(o->_fft_in_jxx[i*(1+o->_N/2)+j], real_c(mul_param), image_c(mul_param));
 
600
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i * (1 + o->_N / 2) + j]);
 
601
                                                mul_complex_f(mul_param, mul_param,
 
602
                                                              ((o->_k[i * (1 + o->_N / 2) + j] == 0.0f) ?
 
603
                                                               0.0f :
 
604
                                                               o->_kx[i] * o->_kx[i] / o->_k[i * (1 + o->_N / 2) + j]));
 
605
                                                init_complex(o->_fft_in_jxx[i * (1 + o->_N / 2) + j], real_c(mul_param), image_c(mul_param));
619
606
                                        }
620
607
                                }
621
608
                                fftw_execute(o->_Jxx_plan);
622
609
 
623
 
                                for ( i = 0 ; i  < o->_M ; ++i)
624
 
                                {
625
 
                                        for ( j  = 0 ; j  < o->_N ; ++j)
626
 
                                        {
627
 
                                                o->_Jxx[i*o->_N+j] += 1.0;
 
610
                                for (i = 0; i < o->_M; ++i) {
 
611
                                        for (j = 0; j < o->_N; ++j) {
 
612
                                                o->_Jxx[i * o->_N + j] += 1.0;
628
613
                                        }
629
614
                                }
630
615
                        }
631
 
                } // section 4
 
616
                } /* section 4 */
632
617
 
633
618
#pragma omp section
634
619
                {
635
 
                        if (o->_do_jacobian)
636
 
                        {
637
 
                                // Jzz
638
 
                                for ( i = 0 ; i  < o->_M ; ++i)
639
 
                                {
640
 
                                        for ( j  = 0 ; j  <= o->_N / 2 ; ++j)
641
 
                                        {
 
620
                        if (o->_do_jacobian) {
 
621
                                /* Jzz */
 
622
                                for (i = 0; i < o->_M; ++i) {
 
623
                                        for (j = 0; j <= o->_N / 2; ++j) {
642
624
                                                fftw_complex mul_param;
643
625
 
644
 
                                                //init_complex(mul_param, -scale, 0);
 
626
                                                /* init_complex(mul_param, -scale, 0); */
645
627
                                                init_complex(mul_param, -1, 0);
646
628
 
647
629
                                                mul_complex_f(mul_param, mul_param, chop_amount);
648
 
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i*(1+o->_N/2)+j]);
649
 
                                                mul_complex_f(mul_param, mul_param, (o->_k[i*(1+o->_N/2)+j] == 0.0 ? 0.0 : o->_kz[j]*o->_kz[j] / o->_k[i*(1+o->_N/2)+j]));
650
 
                                                init_complex(o->_fft_in_jzz[i*(1+o->_N/2)+j], real_c(mul_param), image_c(mul_param));
 
630
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i * (1 + o->_N / 2) + j]);
 
631
                                                mul_complex_f(mul_param, mul_param,
 
632
                                                              ((o->_k[i * (1 + o->_N / 2) + j] == 0.0f) ?
 
633
                                                               0.0f :
 
634
                                                               o->_kz[j] * o->_kz[j] / o->_k[i * (1 + o->_N / 2) + j]));
 
635
                                                init_complex(o->_fft_in_jzz[i * (1 + o->_N / 2) + j], real_c(mul_param), image_c(mul_param));
651
636
                                        }
652
637
                                }
653
638
                                fftw_execute(o->_Jzz_plan);
654
 
                                for ( i = 0 ; i  < o->_M ; ++i)
655
 
                                {
656
 
                                        for ( j  = 0 ; j  < o->_N ; ++j)
657
 
                                        {
658
 
                                                o->_Jzz[i*o->_N+j] += 1.0;
 
639
                                for (i = 0; i < o->_M; ++i) {
 
640
                                        for (j = 0; j < o->_N; ++j) {
 
641
                                                o->_Jzz[i * o->_N + j] += 1.0;
659
642
                                        }
660
643
                                }
661
644
                        }
662
 
                } // section 5
 
645
                } /* section 5 */
663
646
 
664
647
#pragma omp section
665
648
                {
666
 
                        if (o->_do_jacobian)
667
 
                        {
668
 
                                // Jxz
669
 
                                for ( i = 0 ; i  < o->_M ; ++i)
670
 
                                {
671
 
                                        for ( j  = 0 ; j  <= o->_N / 2 ; ++j)
672
 
                                        {
 
649
                        if (o->_do_jacobian) {
 
650
                                /* Jxz */
 
651
                                for (i = 0; i < o->_M; ++i) {
 
652
                                        for (j = 0; j <= o->_N / 2; ++j) {
673
653
                                                fftw_complex mul_param;
674
654
 
675
 
                                                //init_complex(mul_param, -scale, 0);
 
655
                                                /* init_complex(mul_param, -scale, 0); */
676
656
                                                init_complex(mul_param, -1, 0);
677
657
 
678
658
                                                mul_complex_f(mul_param, mul_param, chop_amount);
679
 
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i*(1+o->_N/2)+j]);
680
 
                                                mul_complex_f(mul_param, mul_param, (o->_k[i*(1+o->_N/2)+j] == 0.0f ? 0.0f : o->_kx[i]*o->_kz[j] / o->_k[i*(1+o->_N/2)+j]));
681
 
                                                init_complex(o->_fft_in_jxz[i*(1+o->_N/2)+j], real_c(mul_param), image_c(mul_param));
 
659
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i * (1 + o->_N / 2) + j]);
 
660
                                                mul_complex_f(mul_param, mul_param,
 
661
                                                              ((o->_k[i * (1 + o->_N / 2) + j] == 0.0f) ?
 
662
                                                               0.0f :
 
663
                                                               o->_kx[i] * o->_kz[j] / o->_k[i * (1 + o->_N / 2) + j]));
 
664
                                                init_complex(o->_fft_in_jxz[i * (1 + o->_N / 2) + j], real_c(mul_param), image_c(mul_param));
682
665
                                        }
683
666
                                }
684
667
                                fftw_execute(o->_Jxz_plan);
685
668
                        }
686
 
                } // section 6
 
669
                } /* section 6 */
687
670
 
688
671
#pragma omp section
689
672
                {
690
 
                        // fft normals
691
 
                        if (o->_do_normals)
692
 
                        {
693
 
                                for ( i = 0 ; i  < o->_M ; ++i)
694
 
                                {
695
 
                                        for ( j  = 0 ; j  <= o->_N / 2 ; ++j)
696
 
                                        {
 
673
                        /* fft normals */
 
674
                        if (o->_do_normals) {
 
675
                                for (i = 0; i < o->_M; ++i) {
 
676
                                        for (j = 0; j <= o->_N / 2; ++j) {
697
677
                                                fftw_complex mul_param;
698
678
 
699
679
                                                init_complex(mul_param, 0.0, -1.0);
700
 
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i*(1+o->_N/2)+j]);
 
680
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i * (1 + o->_N / 2) + j]);
701
681
                                                mul_complex_f(mul_param, mul_param, o->_kx[i]);
702
 
                                                init_complex(o->_fft_in_nx[i*(1+o->_N/2)+j], real_c(mul_param), image_c(mul_param));
 
682
                                                init_complex(o->_fft_in_nx[i * (1 + o->_N / 2) + j], real_c(mul_param), image_c(mul_param));
703
683
                                        }
704
684
                                }
705
685
                                fftw_execute(o->_N_x_plan);
706
686
 
707
687
                        }
708
 
                } // section 7
 
688
                } /* section 7 */
709
689
 
710
690
#pragma omp section
711
691
                {
712
 
                        if (o->_do_normals)
713
 
                        {
714
 
                                for ( i = 0 ; i  < o->_M ; ++i)
715
 
                                {
716
 
                                        for ( j  = 0 ; j  <= o->_N / 2 ; ++j)
717
 
                                        {
 
692
                        if (o->_do_normals) {
 
693
                                for (i = 0; i < o->_M; ++i) {
 
694
                                        for (j = 0; j <= o->_N / 2; ++j) {
718
695
                                                fftw_complex mul_param;
719
696
 
720
697
                                                init_complex(mul_param, 0.0, -1.0);
721
 
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i*(1+o->_N/2)+j]);
 
698
                                                mul_complex_c(mul_param, mul_param, o->_htilda[i * (1 + o->_N / 2) + j]);
722
699
                                                mul_complex_f(mul_param, mul_param, o->_kz[i]);
723
 
                                                init_complex(o->_fft_in_nz[i*(1+o->_N/2)+j], real_c(mul_param), image_c(mul_param));
 
700
                                                init_complex(o->_fft_in_nz[i * (1 + o->_N / 2) + j], real_c(mul_param), image_c(mul_param));
724
701
                                        }
725
702
                                }
726
703
                                fftw_execute(o->_N_z_plan);
727
704
 
728
705
#if 0
729
 
                                for ( i = 0 ; i  < o->_M ; ++i) {
730
 
                                        for ( j  = 0 ; j  < o->_N ; ++j) {
731
 
                                                o->_N_y[i*o->_N+j] = 1.0f/scale;
 
706
                                for (i = 0; i < o->_M; ++i) {
 
707
                                        for (j = 0; j < o->_N; ++j) {
 
708
                                                o->_N_y[i * o->_N + j] = 1.0f / scale;
732
709
                                        }
733
710
                                }
734
711
                                (MEM01)
735
712
#endif
736
 
                        o->_N_y = 1.0f/scale;
 
713
                                o->_N_y = 1.0f / scale;
737
714
                        }
738
 
                } // section 8
 
715
                } /* section 8 */
739
716
 
740
 
        } // omp sections
 
717
        } /* omp sections */
741
718
 
742
719
        BLI_rw_mutex_unlock(&o->oceanmutex);
743
720
}
747
724
        float res = 1.0;
748
725
        float max_h = 0.0;
749
726
 
750
 
        int i,j;
 
727
        int i, j;
751
728
 
752
729
        if (!oc->_do_disp_y) return;
753
730
 
757
734
 
758
735
        BLI_rw_mutex_lock(&oc->oceanmutex, THREAD_LOCK_READ);
759
736
 
760
 
        for (i = 0; i < oc->_M; ++i)
761
 
        {
762
 
                for (j = 0; j < oc->_N; ++j)
763
 
                {
764
 
                        if ( max_h < fabsf(oc->_disp_y[i*oc->_N+j]))
765
 
                        {
766
 
                                max_h = fabsf(oc->_disp_y[i*oc->_N+j]);
 
737
        for (i = 0; i < oc->_M; ++i) {
 
738
                for (j = 0; j < oc->_N; ++j) {
 
739
                        if (max_h < fabsf(oc->_disp_y[i * oc->_N + j])) {
 
740
                                max_h = fabsf(oc->_disp_y[i * oc->_N + j]);
767
741
                        }
768
742
                }
769
743
        }
770
744
 
771
745
        BLI_rw_mutex_unlock(&oc->oceanmutex);
772
746
 
773
 
        if (max_h == 0.0f) max_h = 0.00001f; // just in case ...
 
747
        if (max_h == 0.0f)
 
748
                max_h = 0.00001f;  /* just in case ... */
774
749
 
775
750
        res = 1.0f / (max_h);
776
751
 
786
761
        return oc;
787
762
}
788
763
 
789
 
void BKE_init_ocean(struct Ocean* o, int M,int N, float Lx, float Lz, float V, float l, float A, float w, float damp,
790
 
                                        float alignment, float depth, float time, short do_height_field, short do_chop, short do_normals, short do_jacobian, int seed)
 
764
void BKE_init_ocean(struct Ocean *o, int M, int N, float Lx, float Lz, float V, float l, float A, float w, float damp,
 
765
                    float alignment, float depth, float time, short do_height_field, short do_chop, short do_normals,
 
766
                    short do_jacobian, int seed)
791
767
{
792
 
        int i,j,ii;
 
768
        int i, j, ii;
793
769
 
794
770
        BLI_rw_mutex_lock(&o->oceanmutex, THREAD_LOCK_WRITE);
795
771
 
805
781
        o->_Lx = Lx;
806
782
        o->_Lz = Lz;
807
783
        o->_wx = cos(w);
808
 
        o->_wz = -sin(w); // wave direction
809
 
        o->_L = V*V / GRAVITY;  // largest wave for a given velocity V
 
784
        o->_wz = -sin(w); /* wave direction */
 
785
        o->_L = V * V / GRAVITY;  /* largest wave for a given velocity V */
810
786
        o->time = time;
811
787
 
812
788
        o->_do_disp_y = do_height_field;
814
790
        o->_do_chop = do_chop;
815
791
        o->_do_jacobian = do_jacobian;
816
792
 
817
 
        o->_k = (float*) MEM_mallocN(M * (1+N/2) * sizeof(float), "ocean_k");
818
 
        o->_h0 = (fftw_complex*) MEM_mallocN(M * N * sizeof(fftw_complex), "ocean_h0");
819
 
        o->_h0_minus = (fftw_complex*) MEM_mallocN(M * N * sizeof(fftw_complex), "ocean_h0_minus");
820
 
        o->_kx = (float*) MEM_mallocN(o->_M * sizeof(float), "ocean_kx");
821
 
        o->_kz = (float*) MEM_mallocN(o->_N * sizeof(float), "ocean_kz");
 
793
        o->_k = (float *) MEM_mallocN(M * (1 + N / 2) * sizeof(float), "ocean_k");
 
794
        o->_h0 = (fftw_complex *) MEM_mallocN(M * N * sizeof(fftw_complex), "ocean_h0");
 
795
        o->_h0_minus = (fftw_complex *) MEM_mallocN(M * N * sizeof(fftw_complex), "ocean_h0_minus");
 
796
        o->_kx = (float *) MEM_mallocN(o->_M * sizeof(float), "ocean_kx");
 
797
        o->_kz = (float *) MEM_mallocN(o->_N * sizeof(float), "ocean_kz");
822
798
 
823
 
        // make this robust in the face of erroneous usage
 
799
        /* make this robust in the face of erroneous usage */
824
800
        if (o->_Lx == 0.0f)
825
801
                o->_Lx = 0.001f;
826
802
 
827
803
        if (o->_Lz == 0.0f)
828
804
                o->_Lz = 0.001f;
829
805
 
830
 
        // the +ve components and DC
831
 
        for (i = 0 ; i <= o->_M/2 ; ++i)
 
806
        /* the +ve components and DC */
 
807
        for (i = 0; i <= o->_M / 2; ++i)
832
808
                o->_kx[i] = 2.0f * (float)M_PI * i / o->_Lx;
833
809
 
834
 
        // the -ve components
835
 
        for (i = o->_M-1,ii=0 ; i > o->_M/2 ; --i,++ii)
 
810
        /* the -ve components */
 
811
        for (i = o->_M - 1, ii = 0; i > o->_M / 2; --i, ++ii)
836
812
                o->_kx[i] = -2.0f * (float)M_PI * ii / o->_Lx;
837
813
 
838
 
        // the +ve components and DC
839
 
        for (i = 0 ; i <= o->_N/2 ; ++i)
 
814
        /* the +ve components and DC */
 
815
        for (i = 0; i <= o->_N / 2; ++i)
840
816
                o->_kz[i] = 2.0f * (float)M_PI * i / o->_Lz;
841
817
 
842
 
        // the -ve components
843
 
        for (i = o->_N-1,ii=0 ; i > o->_N/2 ; --i,++ii)
 
818
        /* the -ve components */
 
819
        for (i = o->_N - 1, ii = 0; i > o->_N / 2; --i, ++ii)
844
820
                o->_kz[i] = -2.0f * (float)M_PI * ii / o->_Lz;
845
821
 
846
 
        // pre-calculate the k matrix
847
 
        for (i = 0 ; i  < o->_M ; ++i)
848
 
                for (j  = 0 ; j  <= o->_N / 2 ; ++j)
849
 
                        o->_k[i*(1+o->_N/2)+j] = sqrt(o->_kx[i]*o->_kx[i] + o->_kz[j]*o->_kz[j] );
 
822
        /* pre-calculate the k matrix */
 
823
        for (i = 0; i < o->_M; ++i)
 
824
                for (j = 0; j <= o->_N / 2; ++j)
 
825
                        o->_k[i * (1 + o->_N / 2) + j] = sqrt(o->_kx[i] * o->_kx[i] + o->_kz[j] * o->_kz[j]);
850
826
 
851
827
        /*srand(seed);*/
852
828
        BLI_srand(seed);
853
829
 
854
 
        for (i = 0 ; i  < o->_M ; ++i)
855
 
        {
856
 
                for (j = 0 ; j  < o->_N ; ++j)
857
 
                {
 
830
        for (i = 0; i < o->_M; ++i) {
 
831
                for (j = 0; j < o->_N; ++j) {
858
832
                        float r1 = gaussRand();
859
833
                        float r2 = gaussRand();
860
834
 
861
835
                        fftw_complex r1r2;
862
836
                        init_complex(r1r2, r1, r2);
863
 
                        mul_complex_f(o->_h0[i*o->_N+j], r1r2, (float)(sqrt(Ph(o,  o->_kx[i], o->_kz[j]) / 2.0f)));
864
 
                        mul_complex_f(o->_h0_minus[i*o->_N+j], r1r2, (float)(sqrt(Ph(o, -o->_kx[i],-o->_kz[j]) / 2.0f)));
 
837
                        mul_complex_f(o->_h0[i * o->_N + j], r1r2, (float)(sqrt(Ph(o, o->_kx[i], o->_kz[j]) / 2.0f)));
 
838
                        mul_complex_f(o->_h0_minus[i * o->_N + j], r1r2, (float)(sqrt(Ph(o, -o->_kx[i], -o->_kz[j]) / 2.0f)));
865
839
                }
866
840
        }
867
841
 
868
 
        o->_fft_in = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_fft_in");
869
 
        o->_htilda = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_htilda");
 
842
        o->_fft_in = (fftw_complex *)MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex), "ocean_fft_in");
 
843
        o->_htilda = (fftw_complex *)MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex), "ocean_htilda");
870
844
 
871
845
        if (o->_do_disp_y) {
872
 
                o->_disp_y = (double*) MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_disp_y");
873
 
                o->_disp_y_plan = fftw_plan_dft_c2r_2d(o->_M,o->_N, o->_fft_in, o->_disp_y, FFTW_ESTIMATE);
 
846
                o->_disp_y = (double *)MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_disp_y");
 
847
                o->_disp_y_plan = fftw_plan_dft_c2r_2d(o->_M, o->_N, o->_fft_in, o->_disp_y, FFTW_ESTIMATE);
874
848
        }
875
849
 
876
850
        if (o->_do_normals) {
877
 
                o->_fft_in_nx = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_fft_in_nx");
878
 
                o->_fft_in_nz = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_fft_in_nz");
879
 
 
880
 
                o->_N_x = (double*) MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_N_x");
881
 
                /*o->_N_y = (float*) fftwf_malloc(o->_M * o->_N * sizeof(float)); (MEM01)*/
882
 
                o->_N_z = (double*) MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_N_z");
883
 
 
884
 
                o->_N_x_plan = fftw_plan_dft_c2r_2d(o->_M,o->_N, o->_fft_in_nx, o->_N_x, FFTW_ESTIMATE);
885
 
                o->_N_z_plan = fftw_plan_dft_c2r_2d(o->_M,o->_N, o->_fft_in_nz, o->_N_z, FFTW_ESTIMATE);
 
851
                o->_fft_in_nx = (fftw_complex *) MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex), "ocean_fft_in_nx");
 
852
                o->_fft_in_nz = (fftw_complex *) MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex), "ocean_fft_in_nz");
 
853
 
 
854
                o->_N_x = (double *)MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_N_x");
 
855
                /* o->_N_y = (float *) fftwf_malloc(o->_M * o->_N * sizeof(float)); (MEM01) */
 
856
                o->_N_z = (double *)MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_N_z");
 
857
 
 
858
                o->_N_x_plan = fftw_plan_dft_c2r_2d(o->_M, o->_N, o->_fft_in_nx, o->_N_x, FFTW_ESTIMATE);
 
859
                o->_N_z_plan = fftw_plan_dft_c2r_2d(o->_M, o->_N, o->_fft_in_nz, o->_N_z, FFTW_ESTIMATE);
886
860
        }
887
861
 
888
862
        if (o->_do_chop) {
889
 
                o->_fft_in_x = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_fft_in_x");
890
 
                o->_fft_in_z = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_fft_in_z");
891
 
 
892
 
                o->_disp_x = (double*) MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_disp_x");
893
 
                o->_disp_z = (double*) MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_disp_z");
894
 
 
895
 
                o->_disp_x_plan = fftw_plan_dft_c2r_2d(o->_M,o->_N, o->_fft_in_x, o->_disp_x, FFTW_ESTIMATE);
896
 
                o->_disp_z_plan = fftw_plan_dft_c2r_2d(o->_M,o->_N, o->_fft_in_z, o->_disp_z, FFTW_ESTIMATE);
 
863
                o->_fft_in_x = (fftw_complex *)MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex), "ocean_fft_in_x");
 
864
                o->_fft_in_z = (fftw_complex *)MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex), "ocean_fft_in_z");
 
865
 
 
866
                o->_disp_x = (double *)MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_disp_x");
 
867
                o->_disp_z = (double *)MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_disp_z");
 
868
 
 
869
                o->_disp_x_plan = fftw_plan_dft_c2r_2d(o->_M, o->_N, o->_fft_in_x, o->_disp_x, FFTW_ESTIMATE);
 
870
                o->_disp_z_plan = fftw_plan_dft_c2r_2d(o->_M, o->_N, o->_fft_in_z, o->_disp_z, FFTW_ESTIMATE);
897
871
        }
898
872
        if (o->_do_jacobian) {
899
 
                o->_fft_in_jxx = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_fft_in_jxx");
900
 
                o->_fft_in_jzz = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_fft_in_jzz");
901
 
                o->_fft_in_jxz = (fftw_complex*) MEM_mallocN(o->_M * (1+o->_N/2) * sizeof(fftw_complex), "ocean_fft_in_jxz");
902
 
 
903
 
                o->_Jxx = (double*) MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_Jxx");
904
 
                o->_Jzz = (double*) MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_Jzz");
905
 
                o->_Jxz = (double*) MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_Jxz");
906
 
 
907
 
                o->_Jxx_plan = fftw_plan_dft_c2r_2d(o->_M,o->_N, o->_fft_in_jxx, o->_Jxx, FFTW_ESTIMATE);
908
 
                o->_Jzz_plan = fftw_plan_dft_c2r_2d(o->_M,o->_N, o->_fft_in_jzz, o->_Jzz, FFTW_ESTIMATE);
909
 
                o->_Jxz_plan = fftw_plan_dft_c2r_2d(o->_M,o->_N, o->_fft_in_jxz, o->_Jxz, FFTW_ESTIMATE);
 
873
                o->_fft_in_jxx = (fftw_complex *)MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex),
 
874
                                                             "ocean_fft_in_jxx");
 
875
                o->_fft_in_jzz = (fftw_complex *)MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex),
 
876
                                                             "ocean_fft_in_jzz");
 
877
                o->_fft_in_jxz = (fftw_complex *)MEM_mallocN(o->_M * (1 + o->_N / 2) * sizeof(fftw_complex),
 
878
                                                             "ocean_fft_in_jxz");
 
879
 
 
880
                o->_Jxx = (double *)MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_Jxx");
 
881
                o->_Jzz = (double *)MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_Jzz");
 
882
                o->_Jxz = (double *)MEM_mallocN(o->_M * o->_N * sizeof(double), "ocean_Jxz");
 
883
 
 
884
                o->_Jxx_plan = fftw_plan_dft_c2r_2d(o->_M, o->_N, o->_fft_in_jxx, o->_Jxx, FFTW_ESTIMATE);
 
885
                o->_Jzz_plan = fftw_plan_dft_c2r_2d(o->_M, o->_N, o->_fft_in_jzz, o->_Jzz, FFTW_ESTIMATE);
 
886
                o->_Jxz_plan = fftw_plan_dft_c2r_2d(o->_M, o->_N, o->_fft_in_jxz, o->_Jxz, FFTW_ESTIMATE);
910
887
        }
911
888
 
912
889
        BLI_rw_mutex_unlock(&o->oceanmutex);
921
898
 
922
899
        BLI_rw_mutex_lock(&oc->oceanmutex, THREAD_LOCK_WRITE);
923
900
 
924
 
        if (oc->_do_disp_y)
925
 
        {
 
901
        if (oc->_do_disp_y) {
926
902
                fftw_destroy_plan(oc->_disp_y_plan);
927
903
                MEM_freeN(oc->_disp_y);
928
904
        }
929
905
 
930
 
        if (oc->_do_normals)
931
 
        {
 
906
        if (oc->_do_normals) {
932
907
                MEM_freeN(oc->_fft_in_nx);
933
908
                MEM_freeN(oc->_fft_in_nz);
934
909
                fftw_destroy_plan(oc->_N_x_plan);
938
913
                MEM_freeN(oc->_N_z);
939
914
        }
940
915
 
941
 
        if (oc->_do_chop)
942
 
        {
 
916
        if (oc->_do_chop) {
943
917
                MEM_freeN(oc->_fft_in_x);
944
918
                MEM_freeN(oc->_fft_in_z);
945
919
                fftw_destroy_plan(oc->_disp_x_plan);
948
922
                MEM_freeN(oc->_disp_z);
949
923
        }
950
924
 
951
 
        if (oc->_do_jacobian)
952
 
        {
 
925
        if (oc->_do_jacobian) {
953
926
                MEM_freeN(oc->_fft_in_jxx);
954
927
                MEM_freeN(oc->_fft_in_jzz);
955
928
                MEM_freeN(oc->_fft_in_jxz);
993
966
/* ********* Baking/Caching ********* */
994
967
 
995
968
 
996
 
#define CACHE_TYPE_DISPLACE     1
997
 
#define CACHE_TYPE_FOAM         2
998
 
#define CACHE_TYPE_NORMAL       3
 
969
#define CACHE_TYPE_DISPLACE 1
 
970
#define CACHE_TYPE_FOAM     2
 
971
#define CACHE_TYPE_NORMAL   3
999
972
 
1000
973
static void cache_filename(char *string, const char *path, const char *relbase, int frame, int type)
1001
974
{
1002
975
        char cachepath[FILE_MAX];
1003
976
        const char *fname;
1004
977
 
1005
 
        switch(type) {
1006
 
        case CACHE_TYPE_FOAM:
1007
 
                fname= "foam_";
1008
 
                break;
1009
 
        case CACHE_TYPE_NORMAL:
1010
 
                fname= "normal_";
1011
 
                break;
1012
 
        case CACHE_TYPE_DISPLACE:
1013
 
        default:
1014
 
                fname= "disp_";
1015
 
                break;
 
978
        switch (type) {
 
979
                case CACHE_TYPE_FOAM:
 
980
                        fname = "foam_";
 
981
                        break;
 
982
                case CACHE_TYPE_NORMAL:
 
983
                        fname = "normal_";
 
984
                        break;
 
985
                case CACHE_TYPE_DISPLACE:
 
986
                default:
 
987
                        fname = "disp_";
 
988
                        break;
1016
989
        }
1017
990
 
1018
991
        BLI_join_dirfile(cachepath, sizeof(cachepath), path, fname);
1019
992
 
1020
 
        BKE_makepicstring(string, cachepath, relbase, frame, R_IMF_IMTYPE_OPENEXR, 1, TRUE);
 
993
        BKE_makepicstring_from_type(string, cachepath, relbase, frame, R_IMF_IMTYPE_OPENEXR, 1, TRUE);
1021
994
}
1022
995
 
1023
996
/* silly functions but useful to inline when the args do a lot of indirections */
1024
997
MINLINE void rgb_to_rgba_unit_alpha(float r_rgba[4], const float rgb[3])
1025
998
{
1026
 
        r_rgba[0]= rgb[0];
1027
 
        r_rgba[1]= rgb[1];
1028
 
        r_rgba[2]= rgb[2];
1029
 
        r_rgba[3]= 1.0f;
 
999
        r_rgba[0] = rgb[0];
 
1000
        r_rgba[1] = rgb[1];
 
1001
        r_rgba[2] = rgb[2];
 
1002
        r_rgba[3] = 1.0f;
1030
1003
}
1031
1004
MINLINE void value_to_rgba_unit_alpha(float r_rgba[4], const float value)
1032
1005
{
1033
 
        r_rgba[0]= value;
1034
 
        r_rgba[1]= value;
1035
 
        r_rgba[2]= value;
1036
 
        r_rgba[3]= 1.0f;
 
1006
        r_rgba[0] = value;
 
1007
        r_rgba[1] = value;
 
1008
        r_rgba[2] = value;
 
1009
        r_rgba[3] = 1.0f;
1037
1010
}
1038
1011
 
1039
1012
void BKE_free_ocean_cache(struct OceanCache *och)
1040
1013
{
1041
 
        int i, f=0;
 
1014
        int i, f = 0;
1042
1015
 
1043
1016
        if (!och) return;
1044
1017
 
1045
1018
        if (och->ibufs_disp) {
1046
 
                for (i=och->start, f=0; i<=och->end; i++, f++)
1047
 
                {
 
1019
                for (i = och->start, f = 0; i <= och->end; i++, f++) {
1048
1020
                        if (och->ibufs_disp[f]) {
1049
1021
                                IMB_freeImBuf(och->ibufs_disp[f]);
1050
1022
                        }
1053
1025
        }
1054
1026
 
1055
1027
        if (och->ibufs_foam) {
1056
 
                for (i=och->start, f=0; i<=och->end; i++, f++)
1057
 
                {
 
1028
                for (i = och->start, f = 0; i <= och->end; i++, f++) {
1058
1029
                        if (och->ibufs_foam[f]) {
1059
1030
                                IMB_freeImBuf(och->ibufs_foam[f]);
1060
1031
                        }
1063
1034
        }
1064
1035
 
1065
1036
        if (och->ibufs_norm) {
1066
 
                for (i=och->start, f=0; i<=och->end; i++, f++)
1067
 
                {
 
1037
                for (i = och->start, f = 0; i <= och->end; i++, f++) {
1068
1038
                        if (och->ibufs_norm[f]) {
1069
1039
                                IMB_freeImBuf(och->ibufs_norm[f]);
1070
1040
                        }
1090
1060
        if (v < 0) v += 1.0f;
1091
1061
 
1092
1062
        if (och->ibufs_disp[f]) {
1093
 
                ibuf_sample(och->ibufs_disp[f], u, v, (1.0f/(float)res_x), (1.0f/(float)res_y), result);
 
1063
                ibuf_sample(och->ibufs_disp[f], u, v, (1.0f / (float)res_x), (1.0f / (float)res_y), result);
1094
1064
                copy_v3_v3(ocr->disp, result);
1095
1065
        }
1096
1066
 
1097
1067
        if (och->ibufs_foam[f]) {
1098
 
                ibuf_sample(och->ibufs_foam[f], u, v, (1.0f/(float)res_x), (1.0f/(float)res_y), result);
 
1068
                ibuf_sample(och->ibufs_foam[f], u, v, (1.0f / (float)res_x), (1.0f / (float)res_y), result);
1099
1069
                ocr->foam = result[0];
1100
1070
        }
1101
1071
 
1102
1072
        if (och->ibufs_norm[f]) {
1103
 
                ibuf_sample(och->ibufs_norm[f], u, v, (1.0f/(float)res_x), (1.0f/(float)res_y), result);
 
1073
                ibuf_sample(och->ibufs_norm[f], u, v, (1.0f / (float)res_x), (1.0f / (float)res_y), result);
1104
1074
                copy_v3_v3(ocr->normal, result);
1105
1075
        }
1106
1076
}
1110
1080
        const int res_x = och->resolution_x;
1111
1081
        const int res_y = och->resolution_y;
1112
1082
 
1113
 
        if (i < 0) i= -i;
1114
 
        if (j < 0) j= -j;
 
1083
        if (i < 0) i = -i;
 
1084
        if (j < 0) j = -j;
1115
1085
 
1116
1086
        i = i % res_x;
1117
1087
        j = j % res_y;
1118
1088
 
1119
1089
        if (och->ibufs_disp[f]) {
1120
 
                copy_v3_v3(ocr->disp, &och->ibufs_disp[f]->rect_float[4*(res_x*j + i)]);
 
1090
                copy_v3_v3(ocr->disp, &och->ibufs_disp[f]->rect_float[4 * (res_x * j + i)]);
1121
1091
        }
1122
1092
 
1123
1093
        if (och->ibufs_foam[f]) {
1124
 
                ocr->foam = och->ibufs_foam[f]->rect_float[4*(res_x*j + i)];
 
1094
                ocr->foam = och->ibufs_foam[f]->rect_float[4 * (res_x * j + i)];
1125
1095
        }
1126
1096
 
1127
1097
        if (och->ibufs_norm[f]) {
1128
 
                copy_v3_v3(ocr->normal, &och->ibufs_norm[f]->rect_float[4*(res_x*j + i)]);
 
1098
                copy_v3_v3(ocr->normal, &och->ibufs_norm[f]->rect_float[4 * (res_x * j + i)]);
1129
1099
        }
1130
1100
}
1131
1101
 
1132
 
struct OceanCache *BKE_init_ocean_cache(const char *bakepath, const char *relbase,
1133
 
                                        int start, int end, float wave_scale,
 
1102
struct OceanCache *BKE_init_ocean_cache(const char *bakepath, const char *relbase, int start, int end, float wave_scale,
1134
1103
                                        float chop_amount, float foam_coverage, float foam_fade, int resolution)
1135
1104
{
1136
1105
        OceanCache *och = MEM_callocN(sizeof(OceanCache), "ocean cache data");
1145
1114
        och->chop_amount = chop_amount;
1146
1115
        och->foam_coverage = foam_coverage;
1147
1116
        och->foam_fade = foam_fade;
1148
 
        och->resolution_x = resolution*resolution;
1149
 
        och->resolution_y = resolution*resolution;
 
1117
        och->resolution_x = resolution * resolution;
 
1118
        och->resolution_y = resolution * resolution;
1150
1119
 
1151
 
        och->ibufs_disp = MEM_callocN(sizeof(ImBuf *)*och->duration, "displacement imbuf pointer array");
1152
 
        och->ibufs_foam = MEM_callocN(sizeof(ImBuf *)*och->duration, "foam imbuf pointer array");
1153
 
        och->ibufs_norm = MEM_callocN(sizeof(ImBuf *)*och->duration, "normal imbuf pointer array");
 
1120
        och->ibufs_disp = MEM_callocN(sizeof(ImBuf *) * och->duration, "displacement imbuf pointer array");
 
1121
        och->ibufs_foam = MEM_callocN(sizeof(ImBuf *) * och->duration, "foam imbuf pointer array");
 
1122
        och->ibufs_norm = MEM_callocN(sizeof(ImBuf *) * och->duration, "normal imbuf pointer array");
1154
1123
 
1155
1124
        och->time = NULL;
1156
1125
 
1165
1134
        /* ibufs array is zero based, but filenames are based on frame numbers */
1166
1135
        /* still need to clamp frame numbers to valid range of images on disk though */
1167
1136
        CLAMP(frame, och->start, och->end);
1168
 
        f = frame - och->start; // shift to 0 based
 
1137
        f = frame - och->start; /* shift to 0 based */
1169
1138
 
1170
1139
        /* if image is already loaded in mem, return */
1171
 
        if (och->ibufs_disp[f] != NULL ) return;
 
1140
        if (och->ibufs_disp[f] != NULL) return;
1172
1141
 
 
1142
        /* use default color spaces since we know for sure cache files were saved with default settings too */
1173
1143
 
1174
1144
        cache_filename(string, och->bakepath, och->relbase, frame, CACHE_TYPE_DISPLACE);
1175
 
        och->ibufs_disp[f] = IMB_loadiffname(string, 0);
1176
 
        //if (och->ibufs_disp[f] == NULL) printf("error loading %s\n", string);
1177
 
        //else printf("loaded cache %s\n", string);
 
1145
        och->ibufs_disp[f] = IMB_loadiffname(string, 0, NULL);
 
1146
#if 0
 
1147
        if (och->ibufs_disp[f] == NULL)
 
1148
                printf("error loading %s\n", string);
 
1149
        else
 
1150
                printf("loaded cache %s\n", string);
 
1151
#endif
1178
1152
 
1179
1153
        cache_filename(string, och->bakepath, och->relbase, frame, CACHE_TYPE_FOAM);
1180
 
        och->ibufs_foam[f] = IMB_loadiffname(string, 0);
1181
 
        //if (och->ibufs_foam[f] == NULL) printf("error loading %s\n", string);
1182
 
        //else printf("loaded cache %s\n", string);
 
1154
        och->ibufs_foam[f] = IMB_loadiffname(string, 0, NULL);
 
1155
#if 0
 
1156
        if (och->ibufs_foam[f] == NULL)
 
1157
                printf("error loading %s\n", string);
 
1158
        else
 
1159
                printf("loaded cache %s\n", string);
 
1160
#endif
1183
1161
 
1184
1162
        cache_filename(string, och->bakepath, och->relbase, frame, CACHE_TYPE_NORMAL);
1185
 
        och->ibufs_norm[f] = IMB_loadiffname(string, 0);
1186
 
        //if (och->ibufs_norm[f] == NULL) printf("error loading %s\n", string);
1187
 
        //else printf("loaded cache %s\n", string);
 
1163
        och->ibufs_norm[f] = IMB_loadiffname(string, 0, NULL);
 
1164
#if 0
 
1165
        if (och->ibufs_norm[f] == NULL)
 
1166
                printf("error loading %s\n", string);
 
1167
        else
 
1168
                printf("loaded cache %s\n", string);
 
1169
#endif
1188
1170
}
1189
1171
 
1190
1172
 
1191
 
void BKE_bake_ocean(struct Ocean *o, struct OceanCache *och, void (*update_cb)(void *, float progress, int *cancel), void *update_cb_data)
 
1173
void BKE_bake_ocean(struct Ocean *o, struct OceanCache *och, void (*update_cb)(void *, float progress, int *cancel),
 
1174
                    void *update_cb_data)
1192
1175
{
1193
1176
        /* note: some of these values remain uninitialized unless certain options
1194
1177
         * are enabled, take care that BKE_ocean_eval_ij() initializes a member
1195
1178
         * before use - campbell */
1196
1179
        OceanResult ocr;
1197
1180
 
1198
 
        ImageFormatData imf= {0};
 
1181
        ImageFormatData imf = {0};
1199
1182
 
1200
 
        int f, i=0, x, y, cancel=0;
 
1183
        int f, i = 0, x, y, cancel = 0;
1201
1184
        float progress;
1202
1185
 
1203
1186
        ImBuf *ibuf_foam, *ibuf_disp, *ibuf_normal;
1208
1191
 
1209
1192
        if (!o) return;
1210
1193
 
1211
 
        if (o->_do_jacobian) prev_foam = MEM_callocN(res_x*res_y*sizeof(float), "previous frame foam bake data");
1212
 
        else                 prev_foam = NULL;
 
1194
        if (o->_do_jacobian) prev_foam = MEM_callocN(res_x * res_y * sizeof(float), "previous frame foam bake data");
 
1195
        else prev_foam = NULL;
1213
1196
 
1214
1197
        BLI_srand(0);
1215
1198
 
1216
1199
        /* setup image format */
1217
 
        imf.imtype= R_IMF_IMTYPE_OPENEXR;
1218
 
        imf.depth=  R_IMF_CHAN_DEPTH_16;
1219
 
        imf.exr_codec= R_IMF_EXR_CODEC_ZIP;
 
1200
        imf.imtype = R_IMF_IMTYPE_OPENEXR;
 
1201
        imf.depth =  R_IMF_CHAN_DEPTH_16;
 
1202
        imf.exr_codec = R_IMF_EXR_CODEC_ZIP;
1220
1203
 
1221
 
        for (f=och->start, i=0; f<=och->end; f++, i++) {
 
1204
        for (f = och->start, i = 0; f <= och->end; f++, i++) {
1222
1205
 
1223
1206
                /* create a new imbuf to store image for this frame */
1224
1207
                ibuf_foam = IMB_allocImBuf(res_x, res_y, 32, IB_rectfloat);
1225
1208
                ibuf_disp = IMB_allocImBuf(res_x, res_y, 32, IB_rectfloat);
1226
1209
                ibuf_normal = IMB_allocImBuf(res_x, res_y, 32, IB_rectfloat);
1227
1210
 
1228
 
                ibuf_disp->profile = ibuf_foam->profile = ibuf_normal->profile = IB_PROFILE_LINEAR_RGB;
1229
 
 
1230
1211
                BKE_simulate_ocean(o, och->time[i], och->wave_scale, och->chop_amount);
1231
1212
 
1232
1213
                /* add new foam */
1233
 
                for (y=0; y < res_y; y++) {
1234
 
                        for (x=0; x < res_x; x++) {
 
1214
                for (y = 0; y < res_y; y++) {
 
1215
                        for (x = 0; x < res_x; x++) {
1235
1216
 
1236
1217
                                BKE_ocean_eval_ij(o, &ocr, x, y);
1237
1218
 
1238
1219
                                /* add to the image */
1239
 
                                rgb_to_rgba_unit_alpha(&ibuf_disp->rect_float[4*(res_x*y + x)], ocr.disp);
 
1220
                                rgb_to_rgba_unit_alpha(&ibuf_disp->rect_float[4 * (res_x * y + x)], ocr.disp);
1240
1221
 
1241
1222
                                if (o->_do_jacobian) {
1242
1223
                                        /* TODO, cleanup unused code - campbell */
1243
1224
 
1244
 
                                        float /*r,*/ /* UNUSED */ pr=0.0f, foam_result;
 
1225
                                        float /*r, */ /* UNUSED */ pr = 0.0f, foam_result;
1245
1226
                                        float neg_disp, neg_eplus;
1246
1227
 
1247
1228
                                        ocr.foam = BKE_ocean_jminus_to_foam(ocr.Jminus, och->foam_coverage);
1248
1229
 
1249
1230
                                        /* accumulate previous value for this cell */
1250
1231
                                        if (i > 0) {
1251
 
                                                pr = prev_foam[res_x*y + x];
 
1232
                                                pr = prev_foam[res_x * y + x];
1252
1233
                                        }
1253
1234
 
1254
 
                                        /* r = BLI_frand(); */ /* UNUSED */ // randomly reduce foam
1255
 
 
1256
 
                                        //pr = pr * och->foam_fade;             // overall fade
1257
 
 
1258
 
                                        // remember ocean coord sys is Y up!
1259
 
                                        // break up the foam where height (Y) is low (wave valley),
1260
 
                                        // and X and Z displacement is greatest
 
1235
                                        /* r = BLI_frand(); */ /* UNUSED */ /* randomly reduce foam */
 
1236
 
 
1237
                                        /* pr = pr * och->foam_fade; */         /* overall fade */
 
1238
 
 
1239
                                        /* remember ocean coord sys is Y up!
 
1240
                                         * break up the foam where height (Y) is low (wave valley), and X and Z displacement is greatest
 
1241
                                         */
1261
1242
 
1262
1243
#if 0
1263
1244
                                        vec[0] = ocr.disp[0];
1266
1247
                                        CLAMP(hor_stretch, 0.0, 1.0);
1267
1248
#endif
1268
1249
 
1269
 
                                        neg_disp = ocr.disp[1] < 0.0f ? 1.0f+ocr.disp[1] : 1.0f;
 
1250
                                        neg_disp = ocr.disp[1] < 0.0f ? 1.0f + ocr.disp[1] : 1.0f;
1270
1251
                                        neg_disp = neg_disp < 0.0f ? 0.0f : neg_disp;
1271
1252
 
1272
1253
                                        /* foam, 'ocr.Eplus' only initialized with do_jacobian */
1273
 
                                        neg_eplus = ocr.Eplus[2] < 0.0f ? 1.0f + ocr.Eplus[2]:1.0f;
1274
 
                                        neg_eplus = neg_eplus<0.0f ? 0.0f : neg_eplus;
1275
 
 
1276
 
                                        //if (ocr.disp[1] < 0.0 || r > och->foam_fade)
1277
 
                                        //      pr *= och->foam_fade;
1278
 
 
1279
 
 
1280
 
                                        //pr = pr * (1.0 - hor_stretch) * ocr.disp[1];
1281
 
                                        //pr = pr * neg_disp * neg_eplus;
1282
 
 
1283
 
                                        if (pr < 1.0f) pr *=pr;
 
1254
                                        neg_eplus = ocr.Eplus[2] < 0.0f ? 1.0f + ocr.Eplus[2] : 1.0f;
 
1255
                                        neg_eplus = neg_eplus < 0.0f ? 0.0f : neg_eplus;
 
1256
 
 
1257
#if 0
 
1258
                                        if (ocr.disp[1] < 0.0 || r > och->foam_fade)
 
1259
                                                pr *= och->foam_fade;
 
1260
 
 
1261
 
 
1262
                                        pr = pr * (1.0 - hor_stretch) * ocr.disp[1];
 
1263
                                        pr = pr * neg_disp * neg_eplus;
 
1264
#endif
 
1265
 
 
1266
                                        if (pr < 1.0f)
 
1267
                                                pr *= pr;
1284
1268
 
1285
1269
                                        pr *= och->foam_fade * (0.75f + neg_eplus * 0.25f);
1286
1270
 
1287
 
 
1288
 
                                        foam_result = pr + ocr.foam;
1289
 
 
1290
 
                                        prev_foam[res_x*y + x] = foam_result;
1291
 
 
1292
 
                                        value_to_rgba_unit_alpha(&ibuf_foam->rect_float[4*(res_x*y + x)], foam_result);
 
1271
                                        /* A full clamping should not be needed! */
 
1272
                                        foam_result = min_ff(pr + ocr.foam, 1.0f);
 
1273
 
 
1274
                                        prev_foam[res_x * y + x] = foam_result;
 
1275
 
 
1276
                                        /*foam_result = min_ff(foam_result, 1.0f); */
 
1277
 
 
1278
                                        value_to_rgba_unit_alpha(&ibuf_foam->rect_float[4 * (res_x * y + x)], foam_result);
1293
1279
                                }
1294
1280
 
1295
1281
                                if (o->_do_normals) {
1296
 
                                        rgb_to_rgba_unit_alpha(&ibuf_normal->rect_float[4*(res_x*y + x)], ocr.normal);
 
1282
                                        rgb_to_rgba_unit_alpha(&ibuf_normal->rect_float[4 * (res_x * y + x)], ocr.normal);
1297
1283
                                }
1298
1284
                        }
1299
1285
                }
1300
1286
 
1301
1287
                /* write the images */
1302
1288
                cache_filename(string, och->bakepath, och->relbase, f, CACHE_TYPE_DISPLACE);
1303
 
                if (0 == BKE_write_ibuf(ibuf_disp, string, &imf))
 
1289
                if (0 == BKE_imbuf_write(ibuf_disp, string, &imf))
1304
1290
                        printf("Cannot save Displacement File Output to %s\n", string);
1305
1291
 
1306
1292
                if (o->_do_jacobian) {
1307
 
                        cache_filename(string, och->bakepath, och->relbase,  f, CACHE_TYPE_FOAM);
1308
 
                        if (0 == BKE_write_ibuf(ibuf_foam, string, &imf))
 
1293
                        cache_filename(string, och->bakepath, och->relbase, f, CACHE_TYPE_FOAM);
 
1294
                        if (0 == BKE_imbuf_write(ibuf_foam, string, &imf))
1309
1295
                                printf("Cannot save Foam File Output to %s\n", string);
1310
1296
                }
1311
1297
 
1312
1298
                if (o->_do_normals) {
1313
 
                        cache_filename(string, och->bakepath,  och->relbase, f, CACHE_TYPE_NORMAL);
1314
 
                        if (0 == BKE_write_ibuf(ibuf_normal, string, &imf))
 
1299
                        cache_filename(string, och->bakepath, och->relbase, f, CACHE_TYPE_NORMAL);
 
1300
                        if (0 == BKE_imbuf_write(ibuf_normal, string, &imf))
1315
1301
                                printf("Cannot save Normal File Output to %s\n", string);
1316
1302
                }
1317
1303
 
1333
1319
        och->baked = 1;
1334
1320
}
1335
1321
 
1336
 
#else // WITH_OCEANSIM
 
1322
#else /* WITH_OCEANSIM */
1337
1323
 
1338
1324
/* stub */
1339
1325
typedef struct Ocean {
1347
1333
        return 0.0f;
1348
1334
}
1349
1335
 
1350
 
void BKE_ocean_eval_uv(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), float UNUSED(u),float UNUSED(v))
1351
 
{
1352
 
}
1353
 
 
1354
 
// use catmullrom interpolation rather than linear
1355
 
void BKE_ocean_eval_uv_catrom(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), float UNUSED(u),float UNUSED(v))
1356
 
{
1357
 
}
1358
 
 
1359
 
void BKE_ocean_eval_xz(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), float UNUSED(x),float UNUSED(z))
1360
 
{
1361
 
}
1362
 
 
1363
 
void BKE_ocean_eval_xz_catrom(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), float UNUSED(x),float UNUSED(z))
1364
 
{
1365
 
}
1366
 
 
1367
 
void BKE_ocean_eval_ij(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), int UNUSED(i),int UNUSED(j))
 
1336
void BKE_ocean_eval_uv(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), float UNUSED(u), float UNUSED(v))
 
1337
{
 
1338
}
 
1339
 
 
1340
/* use catmullrom interpolation rather than linear */
 
1341
void BKE_ocean_eval_uv_catrom(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), float UNUSED(u),
 
1342
                              float UNUSED(v))
 
1343
{
 
1344
}
 
1345
 
 
1346
void BKE_ocean_eval_xz(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), float UNUSED(x), float UNUSED(z))
 
1347
{
 
1348
}
 
1349
 
 
1350
void BKE_ocean_eval_xz_catrom(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), float UNUSED(x),
 
1351
                              float UNUSED(z))
 
1352
{
 
1353
}
 
1354
 
 
1355
void BKE_ocean_eval_ij(struct Ocean *UNUSED(oc), struct OceanResult *UNUSED(ocr), int UNUSED(i), int UNUSED(j))
1368
1356
{
1369
1357
}
1370
1358
 
1379
1367
        return oc;
1380
1368
}
1381
1369
 
1382
 
void BKE_init_ocean(struct Ocean* UNUSED(o), int UNUSED(M),int UNUSED(N), float UNUSED(Lx), float UNUSED(Lz), float UNUSED(V), float UNUSED(l), float UNUSED(A), float UNUSED(w), float UNUSED(damp),
1383
 
                                           float UNUSED(alignment), float UNUSED(depth), float UNUSED(time), short UNUSED(do_height_field), short UNUSED(do_chop), short UNUSED(do_normals), short UNUSED(do_jacobian), int UNUSED(seed))
 
1370
void BKE_init_ocean(struct Ocean *UNUSED(o), int UNUSED(M), int UNUSED(N), float UNUSED(Lx), float UNUSED(Lz),
 
1371
                    float UNUSED(V), float UNUSED(l), float UNUSED(A), float UNUSED(w), float UNUSED(damp),
 
1372
                    float UNUSED(alignment), float UNUSED(depth), float UNUSED(time), short UNUSED(do_height_field),
 
1373
                    short UNUSED(do_chop), short UNUSED(do_normals), short UNUSED(do_jacobian), int UNUSED(seed))
1384
1374
{
1385
1375
}
1386
1376
 
1405
1395
        MEM_freeN(och);
1406
1396
}
1407
1397
 
1408
 
void BKE_ocean_cache_eval_uv(struct OceanCache *UNUSED(och), struct OceanResult *UNUSED(ocr), int UNUSED(f), float UNUSED(u), float UNUSED(v))
1409
 
{
1410
 
}
1411
 
 
1412
 
void BKE_ocean_cache_eval_ij(struct OceanCache *UNUSED(och), struct OceanResult *UNUSED(ocr), int UNUSED(f), int UNUSED(i), int UNUSED(j))
1413
 
{
1414
 
}
1415
 
 
1416
 
struct OceanCache *BKE_init_ocean_cache(const char *UNUSED(bakepath), const char *UNUSED(relbase),
1417
 
                                        int UNUSED(start), int UNUSED(end), float UNUSED(wave_scale),
1418
 
                                        float UNUSED(chop_amount), float UNUSED(foam_coverage), float UNUSED(foam_fade), int UNUSED(resolution))
 
1398
void BKE_ocean_cache_eval_uv(struct OceanCache *UNUSED(och), struct OceanResult *UNUSED(ocr), int UNUSED(f),
 
1399
                             float UNUSED(u), float UNUSED(v))
 
1400
{
 
1401
}
 
1402
 
 
1403
void BKE_ocean_cache_eval_ij(struct OceanCache *UNUSED(och), struct OceanResult *UNUSED(ocr), int UNUSED(f),
 
1404
                             int UNUSED(i), int UNUSED(j))
 
1405
{
 
1406
}
 
1407
 
 
1408
OceanCache *BKE_init_ocean_cache(const char *UNUSED(bakepath), const char *UNUSED(relbase), int UNUSED(start),
 
1409
                                 int UNUSED(end), float UNUSED(wave_scale), float UNUSED(chop_amount),
 
1410
                                 float UNUSED(foam_coverage), float UNUSED(foam_fade), int UNUSED(resolution))
1419
1411
{
1420
1412
        OceanCache *och = MEM_callocN(sizeof(OceanCache), "ocean cache data");
1421
1413
 
1426
1418
{
1427
1419
}
1428
1420
 
1429
 
void BKE_bake_ocean(struct Ocean *UNUSED(o), struct OceanCache *UNUSED(och), void (*update_cb)(void *, float progress, int *cancel), void *UNUSED(update_cb_data))
 
1421
void BKE_bake_ocean(struct Ocean *UNUSED(o), struct OceanCache *UNUSED(och),
 
1422
                    void (*update_cb)(void *, float progress, int *cancel), void *UNUSED(update_cb_data))
1430
1423
{
1431
1424
        /* unused */
1432
1425
        (void)update_cb;
1433
1426
}
1434
 
#endif // WITH_OCEANSIM
 
1427
#endif /* WITH_OCEANSIM */