12 |
|
#include "ReadWrite.hpp" |
13 |
|
|
14 |
|
|
15 |
< |
void map( double *x, double *y, double *z, double centerX, double centerY, |
16 |
< |
double centerZ, double boxX, double boxY, double boxZ ); |
15 |
> |
void map( double &x, double &y, double &z, |
16 |
> |
double boxX, double boxY, double boxZ ); |
17 |
|
|
18 |
+ |
void rotate( double &x, double &y, double &z, |
19 |
+ |
double theta, double phi, double psi ); |
20 |
+ |
|
21 |
|
char* program_name; |
22 |
|
using namespace std; |
23 |
|
|
48 |
|
|
49 |
|
double water_shell = 10.0; |
50 |
|
double water_padding = 2.5; |
51 |
< |
double lipid_spaceing = 4.0; |
51 |
> |
double lipid_spaceing = 2.5; |
52 |
|
|
53 |
|
srand48( 1337 ); // initialize the random number generator. |
54 |
|
|
155 |
|
// create an fcc lattice in the water box. |
156 |
|
|
157 |
|
|
158 |
< |
index = 0; |
158 |
> |
int ndx = 0; |
159 |
|
for( i=0; i < n_cellX; i++ ){ |
160 |
|
for( j=0; j < n_cellY; j++ ){ |
161 |
|
for( k=0; k < n_cellZ; k++ ){ |
162 |
|
|
163 |
< |
waterX[index] = i * water_cell + x0; |
164 |
< |
waterY[index] = j * water_cell + y0; |
165 |
< |
waterZ[index] = k * water_cell + z0; |
166 |
< |
index++; |
163 |
> |
waterX[ndx] = i * water_cell + x0; |
164 |
> |
waterY[ndx] = j * water_cell + y0; |
165 |
> |
waterZ[ndx] = k * water_cell + z0; |
166 |
> |
ndx++; |
167 |
|
|
168 |
< |
waterX[index] = i * water_cell + 0.5 * water_cell + x0; |
169 |
< |
waterY[index] = j * water_cell + 0.5 * water_cell + y0; |
170 |
< |
waterZ[index] = k * water_cell + z0; |
171 |
< |
index++; |
168 |
> |
waterX[ndx] = i * water_cell + 0.5 * water_cell + x0; |
169 |
> |
waterY[ndx] = j * water_cell + 0.5 * water_cell + y0; |
170 |
> |
waterZ[ndx] = k * water_cell + z0; |
171 |
> |
ndx++; |
172 |
|
|
173 |
< |
waterX[index] = i * water_cell + x0; |
174 |
< |
waterY[index] = j * water_cell + 0.5 * water_cell + y0; |
175 |
< |
waterZ[index] = k * water_cell + 0.5 * water_cell + z0; |
176 |
< |
index++; |
173 |
> |
waterX[ndx] = i * water_cell + x0; |
174 |
> |
waterY[ndx] = j * water_cell + 0.5 * water_cell + y0; |
175 |
> |
waterZ[ndx] = k * water_cell + 0.5 * water_cell + z0; |
176 |
> |
ndx++; |
177 |
|
|
178 |
< |
waterX[index] = i * water_cell + 0.5 * water_cell + x0; |
179 |
< |
waterY[index] = j * water_cell + y0; |
180 |
< |
waterZ[index] = k * water_cell + 0.5 * water_cell + z0; |
181 |
< |
index++; |
178 |
> |
waterX[ndx] = i * water_cell + 0.5 * water_cell + x0; |
179 |
> |
waterY[ndx] = j * water_cell + y0; |
180 |
> |
waterZ[ndx] = k * water_cell + 0.5 * water_cell + z0; |
181 |
> |
ndx++; |
182 |
|
} |
183 |
|
} |
184 |
|
} |
201 |
|
dy = waterY[i] - lipidAtoms[j]->getY(); |
202 |
|
dz = waterZ[i] - lipidAtoms[j]->getZ(); |
203 |
|
|
204 |
< |
map( &dx, &dy, &dz, cx, cy, cz, box_x, box_y, box_z ); |
204 |
> |
map( dx, dy, dz, box_x, box_y, box_z ); |
205 |
|
|
206 |
|
dx2 = dx * dx; |
207 |
|
dy2 = dy * dy; |
214 |
|
} |
215 |
|
} |
216 |
|
} |
217 |
< |
|
217 |
> |
|
218 |
> |
std::cerr << "nTarget before: " << n_h2o_target; |
219 |
> |
|
220 |
|
n_h2o_target += n_deleted * n_lipids; |
221 |
+ |
|
222 |
+ |
std::cerr << ", after: " << n_h2o_target << ", n_deleted: " << n_deleted |
223 |
+ |
<< "\n"; |
224 |
|
|
225 |
|
// find a box size that best suits the number of waters we need. |
226 |
|
|
227 |
|
int done = 0; |
228 |
|
|
229 |
< |
if( n_waters < n_h2o_target ){ |
229 |
> |
if( n_water < n_h2o_target ){ |
230 |
|
|
231 |
|
int n_generated = n_cellX; |
232 |
|
int n_test, nx, ny, nz; |
234 |
|
ny = n_cellY; |
235 |
|
nz = n_cellZ; |
236 |
|
|
237 |
+ |
n_test = 4 * nx * ny * nz; |
238 |
+ |
|
239 |
|
while( n_test < n_h2o_target ){ |
240 |
|
|
241 |
|
nz++; |
244 |
|
|
245 |
|
int n_diff, goodX, goodY, goodZ; |
246 |
|
|
247 |
< |
n_diff = ntest - n_h2o_target; |
247 |
> |
n_diff = n_test - n_h2o_target; |
248 |
|
goodX = nx; |
249 |
|
goodY = ny; |
250 |
|
goodZ = nz; |
251 |
|
|
252 |
|
int test_diff; |
253 |
< |
int n_limit = n_z; |
254 |
< |
n_z = n_cellZ; |
253 |
> |
int n_limit = nz; |
254 |
> |
nz = n_cellZ; |
255 |
|
|
256 |
|
for( i=n_generated; i<=n_limit; i++ ){ |
257 |
|
for( j=i; j<=n_limit; j++ ){ |
280 |
|
n_cellZ = goodZ; |
281 |
|
} |
282 |
|
|
283 |
< |
// we now have the best box size for the simulation. |
284 |
< |
|
285 |
< |
|
283 |
> |
// we now have the best box size for the simulation. Next we |
284 |
> |
// recreate the water box to the new specifications. |
285 |
> |
|
286 |
> |
n_water = n_cellX * n_cellY * n_cellZ * 4; |
287 |
> |
|
288 |
> |
std::cerr << "new waters = " << n_water << "\n"; |
289 |
> |
|
290 |
> |
delete[] waterX; |
291 |
> |
delete[] waterY; |
292 |
> |
delete[] waterZ; |
293 |
> |
|
294 |
> |
waterX = new double[n_water]; |
295 |
> |
waterY = new double[n_water]; |
296 |
> |
waterZ = new double[n_water]; |
297 |
> |
|
298 |
> |
box_x = water_cell * n_cellX; |
299 |
> |
box_y = water_cell * n_cellY; |
300 |
> |
box_z = water_cell * n_cellZ; |
301 |
> |
|
302 |
> |
x0 = 0.0; |
303 |
> |
y0 = 0.0; |
304 |
> |
z0 = 0.0; |
305 |
> |
|
306 |
> |
cx = ( box_x * 0.5 ); |
307 |
> |
cy = ( box_y * 0.5 ); |
308 |
> |
cz = ( box_z * 0.5 ); |
309 |
> |
|
310 |
> |
// create an fcc lattice in the water box. |
311 |
> |
|
312 |
> |
ndx = 0; |
313 |
> |
for( i=0; i < n_cellX; i++ ){ |
314 |
> |
for( j=0; j < n_cellY; j++ ){ |
315 |
> |
for( k=0; k < n_cellZ; k++ ){ |
316 |
> |
|
317 |
> |
waterX[ndx] = i * water_cell + x0; |
318 |
> |
waterY[ndx] = j * water_cell + y0; |
319 |
> |
waterZ[ndx] = k * water_cell + z0; |
320 |
> |
ndx++; |
321 |
> |
|
322 |
> |
waterX[ndx] = i * water_cell + 0.5 * water_cell + x0; |
323 |
> |
waterY[ndx] = j * water_cell + 0.5 * water_cell + y0; |
324 |
> |
waterZ[ndx] = k * water_cell + z0; |
325 |
> |
ndx++; |
326 |
> |
|
327 |
> |
waterX[ndx] = i * water_cell + x0; |
328 |
> |
waterY[ndx] = j * water_cell + 0.5 * water_cell + y0; |
329 |
> |
waterZ[ndx] = k * water_cell + 0.5 * water_cell + z0; |
330 |
> |
ndx++; |
331 |
> |
|
332 |
> |
waterX[ndx] = i * water_cell + 0.5 * water_cell + x0; |
333 |
> |
waterY[ndx] = j * water_cell + y0; |
334 |
> |
waterZ[ndx] = k * water_cell + 0.5 * water_cell + z0; |
335 |
> |
ndx++; |
336 |
> |
} |
337 |
> |
} |
338 |
> |
} |
339 |
> |
|
340 |
> |
// ************************************************************** |
341 |
> |
|
342 |
> |
|
343 |
> |
|
344 |
> |
// start a 3D RSA for the for the lipid placements |
345 |
> |
|
346 |
> |
srand48( 1337 ); |
347 |
> |
|
348 |
> |
int rsaNAtoms = n_lipids * lipidNAtoms; |
349 |
> |
Atom** rsaAtoms = new Atom*[rsaNAtoms]; |
350 |
> |
|
351 |
> |
DirectionalAtom* dAtom; |
352 |
> |
DirectionalAtom* dAtomNew; |
353 |
> |
|
354 |
> |
double rotMat[3][3]; |
355 |
> |
double unitRotMat[3][3]; |
356 |
> |
|
357 |
> |
unitRotMat[0][0] = 1.0; |
358 |
> |
unitRotMat[0][1] = 0.0; |
359 |
> |
unitRotMat[0][2] = 0.0; |
360 |
> |
|
361 |
> |
unitRotMat[1][0] = 0.0; |
362 |
> |
unitRotMat[1][1] = 1.0; |
363 |
> |
unitRotMat[1][2] = 0.0; |
364 |
> |
|
365 |
> |
unitRotMat[2][0] = 0.0; |
366 |
> |
unitRotMat[2][1] = 0.0; |
367 |
> |
unitRotMat[2][2] = 1.0; |
368 |
> |
|
369 |
> |
ndx = 0; |
370 |
> |
for(i=0; i<n_lipids; i++ ){ |
371 |
> |
for(j=0; j<lipidNAtoms; j++){ |
372 |
> |
|
373 |
> |
if( lipidAtoms[j]->isDirectional() ){ |
374 |
> |
dAtom = (DirectionalAtom *)lipidAtoms[j]; |
375 |
> |
|
376 |
> |
dAtomNew = new DirectionalAtom(); |
377 |
> |
dAtomNew->setSUx( dAtom->getSUx() ); |
378 |
> |
dAtomNew->setSUx( dAtom->getSUx() ); |
379 |
> |
dAtomNew->setSUx( dAtom->getSUx() ); |
380 |
> |
|
381 |
> |
dAtom->getA( rotMat ); |
382 |
> |
dAtomNew->setA( rotMat ); |
383 |
> |
|
384 |
> |
rsaAtoms[ndx] = dAtomNew; |
385 |
> |
} |
386 |
> |
else{ |
387 |
> |
|
388 |
> |
rsaAtoms[ndx] = new GeneralAtom(); |
389 |
> |
} |
390 |
> |
|
391 |
> |
rsaAtoms[ndx]->setType( lipidAtoms[j]->getType() ); |
392 |
> |
|
393 |
> |
ndx++; |
394 |
> |
} |
395 |
> |
} |
396 |
> |
|
397 |
> |
double testX, testY, testZ; |
398 |
> |
double theta, phi, psi; |
399 |
> |
double tempX, tempY, tempZ; |
400 |
> |
int reject; |
401 |
> |
int testDX, acceptedDX; |
402 |
> |
|
403 |
> |
rCutSqr = lipid_spaceing * lipid_spaceing; |
404 |
> |
|
405 |
> |
for(i=0; i<n_lipids; i++ ){ |
406 |
> |
done = 0; |
407 |
> |
while( !done ){ |
408 |
> |
|
409 |
> |
testX = drand48() * box_x; |
410 |
> |
testY = drand48() * box_y; |
411 |
> |
testZ = drand48() * box_z; |
412 |
> |
|
413 |
> |
theta = drand48() * 2.0 * M_PI; |
414 |
> |
phi = drand48() * 2.0 * M_PI; |
415 |
> |
psi = drand48() * 2.0 * M_PI; |
416 |
> |
|
417 |
> |
ndx = i * lipidNAtoms; |
418 |
> |
for(j=0; j<lipidNAtoms; j++){ |
419 |
> |
|
420 |
> |
tempX = lipidAtoms[j]->getX(); |
421 |
> |
tempY = lipidAtoms[j]->getY(); |
422 |
> |
tempZ = lipidAtoms[j]->getZ(); |
423 |
> |
|
424 |
> |
rotate( tempX, tempY, tempZ, theta, phi, psi ); |
425 |
> |
|
426 |
> |
rsaAtoms[ndx + j]->setX( tempX + testX ); |
427 |
> |
rsaAtoms[ndx + j]->setY( tempY + testY ); |
428 |
> |
rsaAtoms[ndx + j]->setZ( tempZ + testZ ); |
429 |
> |
} |
430 |
> |
|
431 |
> |
reject = 0; |
432 |
> |
for( j=0; !reject && j<i; j++){ |
433 |
> |
for(k=0; !reject && k<lipidNAtoms; k++){ |
434 |
|
|
435 |
+ |
acceptedDX = j*lipidNAtoms + k; |
436 |
+ |
for(l=0; !reject && l<lipidNAtoms; l++){ |
437 |
+ |
|
438 |
+ |
testDX = ndx + l; |
439 |
|
|
440 |
+ |
dx = rsaAtoms[testDX]->getX() - rsaAtoms[acceptedDX]->getX(); |
441 |
+ |
dy = rsaAtoms[testDX]->getY() - rsaAtoms[acceptedDX]->getY(); |
442 |
+ |
dz = rsaAtoms[testDX]->getZ() - rsaAtoms[acceptedDX]->getZ(); |
443 |
+ |
|
444 |
+ |
map( dx, dy, dz, box_x, box_y, box_z ); |
445 |
+ |
|
446 |
+ |
dx2 = dx * dx; |
447 |
+ |
dy2 = dy * dy; |
448 |
+ |
dz2 = dz * dz; |
449 |
+ |
|
450 |
+ |
dSqr = dx2 + dy2 + dz2; |
451 |
+ |
if( dSqr < rCutSqr ) reject = 1; |
452 |
+ |
} |
453 |
+ |
} |
454 |
+ |
} |
455 |
|
|
456 |
+ |
if( !reject ){ |
457 |
+ |
done = 1; |
458 |
+ |
std::cerr << i << " has been accepted\n"; |
459 |
+ |
} |
460 |
+ |
} |
461 |
+ |
} |
462 |
+ |
|
463 |
+ |
// cut out the waters that overlap with the lipids. |
464 |
+ |
|
465 |
+ |
delete[] isActive; |
466 |
+ |
isActive = new int[n_water]; |
467 |
+ |
for(i=0; i<n_water; i++) isActive[i] = 1; |
468 |
+ |
int n_active = n_water; |
469 |
+ |
rCutSqr = water_padding * water_padding; |
470 |
+ |
|
471 |
+ |
for(i=0; ( (i<n_water) && isActive[i] ); i++){ |
472 |
+ |
for(j=0; ( (j<rsaNAtoms) && isActive[i] ); j++){ |
473 |
|
|
474 |
+ |
dx = waterX[i] - rsaAtoms[j]->getX(); |
475 |
+ |
dy = waterY[i] - rsaAtoms[j]->getY(); |
476 |
+ |
dz = waterZ[i] - rsaAtoms[j]->getZ(); |
477 |
|
|
478 |
+ |
map( dx, dy, dz, box_x, box_y, box_z ); |
479 |
|
|
480 |
+ |
dx2 = dx * dx; |
481 |
+ |
dy2 = dy * dy; |
482 |
+ |
dz2 = dz * dz; |
483 |
+ |
|
484 |
+ |
dSqr = dx2 + dy2 + dz2; |
485 |
+ |
if( dSqr < rCutSqr ){ |
486 |
+ |
isActive[i] = 0; |
487 |
+ |
n_active--; |
488 |
+ |
} |
489 |
+ |
} |
490 |
+ |
} |
491 |
|
|
492 |
+ |
std::cerr << "final n_waters = " << n_active << "\n"; |
493 |
|
|
494 |
< |
int new_nAtoms = group_nAtoms + n_active; |
494 |
> |
// place all of the waters and lipids into one new array |
495 |
> |
|
496 |
> |
int new_nAtoms = rsaNAtoms + n_active; |
497 |
|
Atom** new_atoms = new Atom*[new_nAtoms]; |
498 |
|
|
499 |
< |
index = 0; |
500 |
< |
for(i=0; i<group_nAtoms; i++ ){ |
499 |
> |
ndx = 0; |
500 |
> |
for(i=0; i<rsaNAtoms; i++ ){ |
501 |
|
|
502 |
< |
if( group_atoms[i]->isDirectional() ){ |
503 |
< |
dAtom = (DirectionalAtom *)group_atoms[i]; |
502 |
> |
if( rsaAtoms[i]->isDirectional() ){ |
503 |
> |
dAtom = (DirectionalAtom *)rsaAtoms[i]; |
504 |
|
|
505 |
|
dAtomNew = new DirectionalAtom(); |
506 |
|
dAtomNew->setSUx( dAtom->getSUx() ); |
507 |
|
dAtomNew->setSUx( dAtom->getSUx() ); |
508 |
|
dAtomNew->setSUx( dAtom->getSUx() ); |
509 |
|
|
510 |
+ |
dAtom->getA( rotMat ); |
511 |
|
dAtomNew->setA( rotMat ); |
512 |
|
|
513 |
< |
new_atoms[index] = dAtomNew; |
513 |
> |
new_atoms[ndx] = dAtomNew; |
514 |
|
} |
515 |
|
else{ |
516 |
|
|
517 |
< |
new_atoms[index] = new GeneralAtom(); |
517 |
> |
new_atoms[ndx] = new GeneralAtom(); |
518 |
|
} |
519 |
|
|
520 |
< |
new_atoms[index]->setType( group_atoms[i]->getType() ); |
520 |
> |
new_atoms[ndx]->setType( rsaAtoms[i]->getType() ); |
521 |
|
|
522 |
< |
new_atoms[index]->setX( group_atoms[i]->getX() ); |
523 |
< |
new_atoms[index]->setY( group_atoms[i]->getY() ); |
524 |
< |
new_atoms[index]->setZ( group_atoms[i]->getZ() ); |
522 |
> |
new_atoms[ndx]->setX( rsaAtoms[i]->getX() ); |
523 |
> |
new_atoms[ndx]->setY( rsaAtoms[i]->getY() ); |
524 |
> |
new_atoms[ndx]->setZ( rsaAtoms[i]->getZ() ); |
525 |
|
|
526 |
< |
new_atoms[index]->set_vx( 0.0 ); |
527 |
< |
new_atoms[index]->set_vy( 0.0 ); |
528 |
< |
new_atoms[index]->set_vz( 0.0 ); |
526 |
> |
new_atoms[ndx]->set_vx( 0.0 ); |
527 |
> |
new_atoms[ndx]->set_vy( 0.0 ); |
528 |
> |
new_atoms[ndx]->set_vz( 0.0 ); |
529 |
|
|
530 |
< |
index++; |
530 |
> |
ndx++; |
531 |
|
} |
532 |
|
|
320 |
– |
|
321 |
– |
|
322 |
– |
|
533 |
|
for(i=0; i<n_water; i++){ |
534 |
|
if(isActive[i]){ |
535 |
|
|
536 |
< |
new_atoms[index] = new DirectionalAtom(); |
537 |
< |
new_atoms[index]->setType( "SSD" ); |
536 |
> |
new_atoms[ndx] = new DirectionalAtom(); |
537 |
> |
new_atoms[ndx]->setType( "SSD" ); |
538 |
|
|
539 |
< |
new_atoms[index]->setX( waterX[i] ); |
540 |
< |
new_atoms[index]->setY( waterY[i] ); |
541 |
< |
new_atoms[index]->setZ( waterZ[i] ); |
539 |
> |
new_atoms[ndx]->setX( waterX[i] ); |
540 |
> |
new_atoms[ndx]->setY( waterY[i] ); |
541 |
> |
new_atoms[ndx]->setZ( waterZ[i] ); |
542 |
|
|
543 |
< |
new_atoms[index]->set_vx( 0.0 ); |
544 |
< |
new_atoms[index]->set_vy( 0.0 ); |
545 |
< |
new_atoms[index]->set_vz( 0.0 ); |
543 |
> |
new_atoms[ndx]->set_vx( 0.0 ); |
544 |
> |
new_atoms[ndx]->set_vy( 0.0 ); |
545 |
> |
new_atoms[ndx]->set_vz( 0.0 ); |
546 |
|
|
547 |
< |
dAtom = (DirectionalAtom *) new_atoms[index]; |
547 |
> |
dAtom = (DirectionalAtom *) new_atoms[ndx]; |
548 |
|
|
549 |
|
dAtom->setSUx( 0.0 ); |
550 |
|
dAtom->setSUy( 0.0 ); |
551 |
|
dAtom->setSUz( 1.0 ); |
552 |
|
|
553 |
< |
dAtom->setA( rotMat ); |
553 |
> |
dAtom->setA( unitRotMat ); |
554 |
|
|
555 |
< |
index++; |
555 |
> |
ndx++; |
556 |
|
} |
557 |
|
} |
558 |
|
|
600 |
|
} |
601 |
|
|
602 |
|
|
603 |
< |
void map( x, y, z, centerX, centerY, centerZ, boxX, boxY, boxZ ) |
604 |
< |
double *x, *y, *z; |
605 |
< |
double centerX, centerY, centerZ; |
606 |
< |
double boxX, boxY, boxZ; |
607 |
< |
{ |
603 |
> |
void map( double &x, double &y, double &z, |
604 |
> |
double boxX, double boxY, double boxZ ){ |
605 |
> |
|
606 |
> |
if(x < 0) x -= boxX * (double)( (int)( (x / boxX) - 0.5 ) ); |
607 |
> |
else x -= boxX * (double)( (int)( (x / boxX ) + 0.5)); |
608 |
|
|
609 |
< |
*x -= centerX; |
610 |
< |
*y -= centerY; |
611 |
< |
*z -= centerZ; |
609 |
> |
if(y < 0) y -= boxY * (double)( (int)( (y / boxY) - 0.5 ) ); |
610 |
> |
else y -= boxY * (double)( (int)( (y / boxY ) + 0.5)); |
611 |
> |
|
612 |
> |
if(z < 0) z -= boxZ * (double)( (int)( (z / boxZ) - 0.5 ) ); |
613 |
> |
else z -= boxZ * (double)( (int)( (z / boxZ ) + 0.5)); |
614 |
> |
} |
615 |
|
|
403 |
– |
if(*x < 0) *x -= boxX * (double)( (int)( (*x / boxX) - 0.5 ) ); |
404 |
– |
else *x -= boxX * (double)( (int)( (*x / boxX ) + 0.5)); |
616 |
|
|
617 |
< |
if(*y < 0) *y -= boxY * (double)( (int)( (*y / boxY) - 0.5 ) ); |
618 |
< |
else *y -= boxY * (double)( (int)( (*y / boxY ) + 0.5)); |
617 |
> |
void rotate( double &x, double &y, double &z, |
618 |
> |
double theta, double phi, double psi ){ |
619 |
> |
|
620 |
> |
double newX, newY, newZ; |
621 |
> |
|
622 |
> |
double A[3][3]; |
623 |
> |
|
624 |
> |
A[0][0] = (cos(phi) * cos(psi)) - (sin(phi) * cos(theta) * sin(psi)); |
625 |
> |
A[0][1] = (sin(phi) * cos(psi)) + (cos(phi) * cos(theta) * sin(psi)); |
626 |
> |
A[0][2] = sin(theta) * sin(psi); |
627 |
|
|
628 |
< |
if(*z < 0) *z -= boxZ * (double)( (int)( (*z / boxZ) - 0.5 ) ); |
629 |
< |
else *z -= boxZ * (double)( (int)( (*z / boxZ ) + 0.5)); |
628 |
> |
A[1][0] = -(cos(phi) * sin(psi)) - (sin(phi) * cos(theta) * cos(psi)); |
629 |
> |
A[1][1] = -(sin(phi) * sin(psi)) + (cos(phi) * cos(theta) * cos(psi)); |
630 |
> |
A[1][2] = sin(theta) * cos(psi); |
631 |
> |
|
632 |
> |
A[2][0] = sin(phi) * sin(theta); |
633 |
> |
A[2][1] = -cos(phi) * sin(theta); |
634 |
> |
A[2][2] = cos(theta); |
635 |
> |
|
636 |
> |
newX = (x * A[0][0]) + (y * A[0][1]) + (z * A[0][2]); |
637 |
> |
newY = (x * A[1][0]) + (y * A[1][1]) + (z * A[1][2]); |
638 |
> |
newZ = (x * A[2][0]) + (y * A[2][1]) + (z * A[2][2]); |
639 |
> |
|
640 |
> |
x = newX; |
641 |
> |
y = newY; |
642 |
> |
z = newZ; |
643 |
|
} |