| 1 | #define _FILE_OFFSET_BITS 64 | 
| 2 |  | 
| 3 | #include <stdio.h> | 
| 4 | #include <stdlib.h> | 
| 5 | #include <sys/types.h> | 
| 6 | #include <sys/stat.h> | 
| 7 | #include <string.h> | 
| 8 |  | 
| 9 | #include "params.h" | 
| 10 | #include "tcProps.h" | 
| 11 | #include "readWrite.h" | 
| 12 |  | 
| 13 |  | 
| 14 | #define BUFFER_SIZE 2000 | 
| 15 | int isScanned; | 
| 16 | int fileOpen; | 
| 17 | int nFrames; | 
| 18 | double *frameTimes; | 
| 19 | FILE* inFile; | 
| 20 | char* inName; | 
| 21 |  | 
| 22 | struct linkedPos{ | 
| 23 |  | 
| 24 | fpos_t *myPos; | 
| 25 | double timeStamp; | 
| 26 | double Hmat[3][3]; | 
| 27 | struct linkedPos* next; | 
| 28 | }; | 
| 29 |  | 
| 30 | struct staticPos{ | 
| 31 |  | 
| 32 | fpos_t *myPos; | 
| 33 | double Hmat[3][3]; | 
| 34 | }; | 
| 35 | struct staticPos* posArray; | 
| 36 |  | 
| 37 | void parseDumpLine( char readBuffer[BUFFER_SIZE], int index, | 
| 38 | struct atomCoord* atoms ); | 
| 39 | void setU( double q[4], double u[3] ); | 
| 40 | void normalizeQ( double q[4] ); | 
| 41 |  | 
| 42 | void openFile( void ){ | 
| 43 |  | 
| 44 | inFile = fopen(inName, "r"); | 
| 45 | if(inFile ==NULL){ | 
| 46 | fprintf(stderr, | 
| 47 | "Error opening file \"%s\"\n", | 
| 48 | inName); | 
| 49 | exit(0); | 
| 50 | } | 
| 51 |  | 
| 52 | fileOpen = 1; | 
| 53 | } | 
| 54 |  | 
| 55 | void closeFile( void ){ | 
| 56 |  | 
| 57 | fclose( inFile ); | 
| 58 |  | 
| 59 | fileOpen = 0; | 
| 60 | } | 
| 61 |  | 
| 62 | void setFrames( void ){ | 
| 63 |  | 
| 64 | int i,j,k; | 
| 65 | struct linkedPos* headPos; | 
| 66 | struct linkedPos* currPos; | 
| 67 | fpos_t *currPT; | 
| 68 | char readBuffer[BUFFER_SIZE]; | 
| 69 | char trash[BUFFER_SIZE]; | 
| 70 | char* foo; | 
| 71 | int lineNum = 0; | 
| 72 |  | 
| 73 |  | 
| 74 | if( !fileOpen ){ | 
| 75 | fprintf(stderr, | 
| 76 | "Attempt to scan the file without first opening the file.\n" ); | 
| 77 | exit(0); | 
| 78 | } | 
| 79 |  | 
| 80 |  | 
| 81 | nFrames = 0; | 
| 82 | headPos = (struct linkedPos*)malloc(sizeof(struct linkedPos)); | 
| 83 | currPos = headPos; | 
| 84 |  | 
| 85 | currPT = (fpos_t *)malloc(sizeof(fpos_t)); | 
| 86 | fgetpos(inFile, currPT); | 
| 87 |  | 
| 88 | fgets( readBuffer, sizeof( readBuffer ), inFile ); | 
| 89 | lineNum++; | 
| 90 | if( feof( inFile ) ){ | 
| 91 | fprintf( stderr, | 
| 92 | "File \"%s\" ended unexpectedly at line %d\n", | 
| 93 | inName, | 
| 94 | lineNum ); | 
| 95 | exit(0); | 
| 96 | } | 
| 97 |  | 
| 98 | while( !feof( inFile ) ){ | 
| 99 |  | 
| 100 | currPos->next = (struct linkedPos*)malloc(sizeof(struct linkedPos)); | 
| 101 | currPos = currPos->next; | 
| 102 | currPos->myPos = currPT; | 
| 103 | nFrames++; | 
| 104 |  | 
| 105 | i = atoi(readBuffer); | 
| 106 |  | 
| 107 | fgets( readBuffer, sizeof( readBuffer ), inFile ); | 
| 108 | lineNum++; | 
| 109 | if( feof( inFile ) ){ | 
| 110 | fprintf( stderr, | 
| 111 | "File \"%s\" ended unexpectedly at line %d\n", | 
| 112 | inName, | 
| 113 | lineNum ); | 
| 114 | exit(0); | 
| 115 | } | 
| 116 |  | 
| 117 |  | 
| 118 | // parse the comment line | 
| 119 |  | 
| 120 | foo = strtok(readBuffer, " ,;\t"); | 
| 121 |  | 
| 122 | if(foo == NULL){ | 
| 123 | fprintf(stderr, | 
| 124 | "error in reading time from %s at line %d\n", | 
| 125 | inName, lineNum ); | 
| 126 | exit(0); | 
| 127 | } | 
| 128 | currPos->timeStamp = atof( foo ); | 
| 129 |  | 
| 130 | // get the Hx vector | 
| 131 |  | 
| 132 | foo = strtok(NULL, " ,;\t"); | 
| 133 | if(foo == NULL){ | 
| 134 | fprintf(stderr, | 
| 135 | "error in reading Hx[0] from %s at line %d\n", | 
| 136 | inName, lineNum ); | 
| 137 | exit(0); | 
| 138 | } | 
| 139 | currPos->Hmat[0][0] = atof( foo ); | 
| 140 |  | 
| 141 | foo = strtok(NULL, " ,;\t"); | 
| 142 | if(foo == NULL){ | 
| 143 | fprintf(stderr, | 
| 144 | "error in reading Hx[1] from %s at line %d\n", | 
| 145 | inName, lineNum ); | 
| 146 | exit(0); | 
| 147 | } | 
| 148 | currPos->Hmat[1][0] = atof( foo ); | 
| 149 |  | 
| 150 | foo = strtok(NULL, " ,;\t"); | 
| 151 | if(foo == NULL){ | 
| 152 | fprintf(stderr, | 
| 153 | "error in reading Hx[2] from %s at line %d\n", | 
| 154 | inName, lineNum ); | 
| 155 | exit(0); | 
| 156 | } | 
| 157 | currPos->Hmat[2][0] = atof( foo ); | 
| 158 |  | 
| 159 | // get the Hy vector | 
| 160 |  | 
| 161 | foo = strtok(NULL, " ,;\t"); | 
| 162 | if(foo == NULL){ | 
| 163 | fprintf(stderr, | 
| 164 | "error in reading Hy[0] from %s at line %d\n", | 
| 165 | inName, lineNum ); | 
| 166 | exit(0); | 
| 167 | } | 
| 168 | currPos->Hmat[0][1] = atof( foo ); | 
| 169 |  | 
| 170 | foo = strtok(NULL, " ,;\t"); | 
| 171 | if(foo == NULL){ | 
| 172 | fprintf(stderr, | 
| 173 | "error in reading Hy[1] from %s at line %d\n", | 
| 174 | inName, lineNum ); | 
| 175 | exit(0); | 
| 176 | } | 
| 177 | currPos->Hmat[1][1] = atof( foo ); | 
| 178 |  | 
| 179 | foo = strtok(NULL, " ,;\t"); | 
| 180 | if(foo == NULL){ | 
| 181 | fprintf(stderr, | 
| 182 | "error in reading Hy[2] from %s at line %d\n", | 
| 183 | inName, lineNum ); | 
| 184 | exit(0); | 
| 185 | } | 
| 186 | currPos->Hmat[2][1] = atof( foo ); | 
| 187 |  | 
| 188 | // get the Hz vector | 
| 189 |  | 
| 190 | foo = strtok(NULL, " ,;\t"); | 
| 191 | if(foo == NULL){ | 
| 192 | fprintf(stderr, | 
| 193 | "error in reading Hz[0] from %s at line %d\n", | 
| 194 | inName, lineNum ); | 
| 195 | exit(0); | 
| 196 | } | 
| 197 | currPos->Hmat[0][2] = atof( foo ); | 
| 198 |  | 
| 199 | foo = strtok(NULL, " ,;\t"); | 
| 200 | if(foo == NULL){ | 
| 201 | fprintf(stderr, | 
| 202 | "error in reading Hz[1] from %s at line %d\n", | 
| 203 | inName, lineNum ); | 
| 204 | exit(0); | 
| 205 | } | 
| 206 | currPos->Hmat[1][2] = atof( foo ); | 
| 207 |  | 
| 208 | foo = strtok(NULL, " ,;\t"); | 
| 209 | if(foo == NULL){ | 
| 210 | fprintf(stderr, | 
| 211 | "error in reading Hz[2] from %s at line %d\n", | 
| 212 | inName, lineNum ); | 
| 213 | exit(0); | 
| 214 | } | 
| 215 | currPos->Hmat[2][2] = atof( foo ); | 
| 216 |  | 
| 217 | // skip the atoms | 
| 218 |  | 
| 219 | for(j=0; j<i; j++){ | 
| 220 |  | 
| 221 | fgets( readBuffer, sizeof( readBuffer ), inFile ); | 
| 222 | lineNum++; | 
| 223 | if( feof( inFile ) ){ | 
| 224 | fprintf( stderr, | 
| 225 | "File \"%s\" ended unexpectedly at line %d," | 
| 226 | " with atom %d\n", | 
| 227 | inName, | 
| 228 | lineNum, | 
| 229 | j ); | 
| 230 | exit(0); | 
| 231 | } | 
| 232 | } | 
| 233 |  | 
| 234 | currPT = (fpos_t *)malloc(sizeof(fpos_t)); | 
| 235 | fgetpos(inFile, currPT); | 
| 236 |  | 
| 237 | fgets( readBuffer, sizeof( readBuffer ), inFile ); | 
| 238 | lineNum++; | 
| 239 | } | 
| 240 |  | 
| 241 | // allocate the static position array | 
| 242 |  | 
| 243 | posArray = (struct staticPos*)calloc(nFrames, sizeof(struct staticPos)); | 
| 244 | frameTimes = (double *)calloc(nFrames, sizeof(double) ); | 
| 245 |  | 
| 246 | i=0; | 
| 247 | currPos = headPos->next; | 
| 248 | while( currPos != NULL ){ | 
| 249 |  | 
| 250 | if(i>=nFrames){ | 
| 251 |  | 
| 252 | fprintf( stderr, | 
| 253 | "The number of frame pointers exceeded the number of frames.\n" | 
| 254 | "  OOPS!\n" ); | 
| 255 | exit(0); | 
| 256 | } | 
| 257 |  | 
| 258 | frameTimes[i] = currPos->timeStamp; | 
| 259 | posArray[i].myPos = currPos->myPos; | 
| 260 | for(j=0;j<3;j++){ | 
| 261 | for(k=0;k<3;k++){ | 
| 262 | posArray[i].Hmat[j][k] = currPos->Hmat[j][k]; | 
| 263 | } | 
| 264 | } | 
| 265 |  | 
| 266 | i++; | 
| 267 | currPos = currPos->next; | 
| 268 | } | 
| 269 |  | 
| 270 | // clean up the linked list | 
| 271 |  | 
| 272 | currPos = headPos; | 
| 273 | while( currPos != NULL ){ | 
| 274 | headPos = currPos; | 
| 275 | currPos = headPos->next; | 
| 276 | free(headPos); | 
| 277 | } | 
| 278 |  | 
| 279 |  | 
| 280 | isScanned = 1; | 
| 281 | } | 
| 282 |  | 
| 283 |  | 
| 284 | void readFrame(int theFrame, struct atomCoord* atoms, double Hmat[3][3] ){ | 
| 285 |  | 
| 286 | // list of 'a priori' constants | 
| 287 |  | 
| 288 | const int nLipAtoms = NL_ATOMS; | 
| 289 | const int nBonds    = NBONDS; | 
| 290 | const int nLipids   = NLIPIDS; | 
| 291 | const int nSSD      = NSSD; | 
| 292 | const int nAtoms    = nLipAtoms * nLipids + nSSD; | 
| 293 |  | 
| 294 | fpos_t *framePos; | 
| 295 | char readBuffer[BUFFER_SIZE]; | 
| 296 |  | 
| 297 | int i, j, k; | 
| 298 |  | 
| 299 |  | 
| 300 | // quick checks | 
| 301 |  | 
| 302 | if (theFrame >= nFrames){ | 
| 303 |  | 
| 304 | fprintf(stderr, | 
| 305 | "frame %d, is out of range\n", | 
| 306 | theFrame ); | 
| 307 | exit(0); | 
| 308 | } | 
| 309 |  | 
| 310 | if( !fileOpen ){ | 
| 311 |  | 
| 312 | fprintf( stderr, | 
| 313 | "File is closed, cannot read frame %d.\n", | 
| 314 | theFrame ); | 
| 315 | exit(0); | 
| 316 | } | 
| 317 |  | 
| 318 | // access the frame | 
| 319 |  | 
| 320 | framePos = posArray[theFrame].myPos; | 
| 321 | fsetpos( inFile, framePos ); | 
| 322 |  | 
| 323 | for(j=0;j<3;j++){ | 
| 324 | for(k=0;k<3;k++){ | 
| 325 | Hmat[j][k] = posArray[theFrame].Hmat[j][k]; | 
| 326 | } | 
| 327 | } | 
| 328 |  | 
| 329 | fgets( readBuffer, sizeof( readBuffer ), inFile ); | 
| 330 | if( feof( inFile ) ){ | 
| 331 | fprintf( stderr, | 
| 332 | "File \"%s\" ended unexpectedly in Frame %d\n", | 
| 333 | inName, | 
| 334 | theFrame ); | 
| 335 | exit(0); | 
| 336 | } | 
| 337 |  | 
| 338 | i = atoi( readBuffer ); | 
| 339 | if( i != nAtoms ){ | 
| 340 | fprintf( stderr, | 
| 341 | "Error in frame %d: frame atoms, %d, does not params atoms, %d\n", | 
| 342 | theFrame, | 
| 343 | i, | 
| 344 | nAtoms ); | 
| 345 | exit(0); | 
| 346 | } | 
| 347 |  | 
| 348 | // read and toss the comment line | 
| 349 |  | 
| 350 | fgets( readBuffer, sizeof( readBuffer ), inFile ); | 
| 351 | if( feof( inFile ) ){ | 
| 352 | fprintf( stderr, | 
| 353 | "File \"%s\" ended unexpectedly in Frame %d\n", | 
| 354 | inName, | 
| 355 | theFrame ); | 
| 356 | exit(0); | 
| 357 | } | 
| 358 |  | 
| 359 | // read the atoms | 
| 360 |  | 
| 361 | for( i=0; i < nAtoms; i++){ | 
| 362 |  | 
| 363 | fgets( readBuffer, sizeof( readBuffer ), inFile ); | 
| 364 | if( feof( inFile ) ){ | 
| 365 | fprintf( stderr, | 
| 366 | "File \"%s\" ended unexpectedly in Frame %d at atom %d\n", | 
| 367 | inName, | 
| 368 | theFrame, | 
| 369 | i ); | 
| 370 | exit(0); | 
| 371 | } | 
| 372 |  | 
| 373 |  | 
| 374 | parseDumpLine( readBuffer, i, atoms ); | 
| 375 | } | 
| 376 | } | 
| 377 |  | 
| 378 | void parseDumpLine( char readLine[BUFFER_SIZE], int i, | 
| 379 | struct atomCoord* atoms ){ | 
| 380 |  | 
| 381 | const int nLipAtoms = NL_ATOMS; | 
| 382 | const int nBonds    = NBONDS; | 
| 383 | const int nLipids   = NLIPIDS; | 
| 384 | const int nSSD      = NSSD; | 
| 385 | const int nAtoms    = nLipAtoms * nLipids + nSSD; | 
| 386 |  | 
| 387 | char* foo; | 
| 388 | double q[4]; | 
| 389 |  | 
| 390 | // set the string tokenizer | 
| 391 |  | 
| 392 | foo = strtok(readLine, " ,;\t"); | 
| 393 |  | 
| 394 | // check the atom ID | 
| 395 |  | 
| 396 | switch( atoms[i].type ){ | 
| 397 |  | 
| 398 | case HEAD: | 
| 399 | if( strcmp(foo, "HEAD") ){ | 
| 400 | fprintf( stderr, | 
| 401 | "Atom %s does not match master array type \"HEAD\".\n", | 
| 402 | foo ); | 
| 403 | exit(0); | 
| 404 | } | 
| 405 | break; | 
| 406 |  | 
| 407 | case CH: | 
| 408 | if( strcmp(foo, "CH") ){ | 
| 409 | fprintf( stderr, | 
| 410 | "Atom %s does not match master array type \"CH\".\n", | 
| 411 | foo ); | 
| 412 | exit(0); | 
| 413 | } | 
| 414 | break; | 
| 415 |  | 
| 416 | case CH2: | 
| 417 | if( strcmp(foo, "CH2") ){ | 
| 418 | fprintf( stderr, | 
| 419 | "Atom %s does not match master array type \"CH2\".\n", | 
| 420 | foo ); | 
| 421 | exit(0); | 
| 422 | } | 
| 423 | break; | 
| 424 |  | 
| 425 | case CH3: | 
| 426 | if( strcmp(foo, "CH3") ){ | 
| 427 | fprintf( stderr, | 
| 428 | "Atom %s does not match master array type \"CH3\".\n", | 
| 429 | foo ); | 
| 430 | exit(0); | 
| 431 | } | 
| 432 | break; | 
| 433 |  | 
| 434 | case SSD: | 
| 435 | if( strcmp(foo, "SSD") ){ | 
| 436 | fprintf( stderr, | 
| 437 | "Atom %s does not match master array type \"SSD\".\n", | 
| 438 | foo ); | 
| 439 | exit(0); | 
| 440 | } | 
| 441 | break; | 
| 442 |  | 
| 443 | default: | 
| 444 | fprintf(stderr, | 
| 445 | "Atom %d is an unrecognized enumeration\n", | 
| 446 | i ); | 
| 447 | exit(0); | 
| 448 | } | 
| 449 |  | 
| 450 | // get the positions | 
| 451 |  | 
| 452 | foo = strtok(NULL, " ,;\t"); | 
| 453 | if(foo == NULL){ | 
| 454 | fprintf( stderr, | 
| 455 | "error in reading postition x from %s\n" | 
| 456 | "natoms  = %d, i = %d\n", | 
| 457 | inName, nAtoms, i ); | 
| 458 | exit(0); | 
| 459 | } | 
| 460 | atoms[i].pos[0] = atof( foo ); | 
| 461 |  | 
| 462 | foo = strtok(NULL, " ,;\t"); | 
| 463 | if(foo == NULL){ | 
| 464 | fprintf( stderr, | 
| 465 | "error in reading postition y from %s\n" | 
| 466 | "natoms  = %d, i = %d\n", | 
| 467 | inName, nAtoms, i ); | 
| 468 | exit(0); | 
| 469 | } | 
| 470 | atoms[i].pos[1] = atof( foo ); | 
| 471 |  | 
| 472 | foo = strtok(NULL, " ,;\t"); | 
| 473 | if(foo == NULL){ | 
| 474 | fprintf( stderr, | 
| 475 | "error in reading postition z from %s\n" | 
| 476 | "natoms  = %d, i = %d\n", | 
| 477 | inName, nAtoms, i ); | 
| 478 | exit(0); | 
| 479 | } | 
| 480 | atoms[i].pos[2] = atof( foo ); | 
| 481 |  | 
| 482 |  | 
| 483 | // get the velocities | 
| 484 |  | 
| 485 | foo = strtok(NULL, " ,;\t"); | 
| 486 | if(foo == NULL){ | 
| 487 | fprintf( stderr, | 
| 488 | "error in reading velocity x from %s\n" | 
| 489 | "natoms  = %d, i = %d\n", | 
| 490 | inName, nAtoms, i ); | 
| 491 | exit(0); | 
| 492 | } | 
| 493 | atoms[i].vel[0] = atof( foo ); | 
| 494 |  | 
| 495 | foo = strtok(NULL, " ,;\t"); | 
| 496 | if(foo == NULL){ | 
| 497 | fprintf( stderr, | 
| 498 | "error in reading velocity y from %s\n" | 
| 499 | "natoms  = %d, i = %d\n", | 
| 500 | inName, nAtoms, i ); | 
| 501 | exit(0); | 
| 502 | } | 
| 503 | atoms[i].vel[1] = atof( foo ); | 
| 504 |  | 
| 505 | foo = strtok(NULL, " ,;\t"); | 
| 506 | if(foo == NULL){ | 
| 507 | fprintf( stderr, | 
| 508 | "error in reading velocity z from %s\n" | 
| 509 | "natoms  = %d, i = %d\n", | 
| 510 | inName, nAtoms, i ); | 
| 511 | exit(0); | 
| 512 | } | 
| 513 | atoms[i].vel[2] = atof( foo ); | 
| 514 |  | 
| 515 |  | 
| 516 | // get the quaternions | 
| 517 |  | 
| 518 | if( (atoms[i].type == HEAD) || (atoms[i].type == SSD) ){ | 
| 519 |  | 
| 520 | foo = strtok(NULL, " ,;\t"); | 
| 521 | if(foo == NULL){ | 
| 522 | fprintf( stderr, | 
| 523 | "error in reading quaternion 0 from %s\n" | 
| 524 | "natoms  = %d, i = %d\n", | 
| 525 | inName, nAtoms, i ); | 
| 526 | exit(0); | 
| 527 | } | 
| 528 | q[0] = atof( foo ); | 
| 529 |  | 
| 530 | foo = strtok(NULL, " ,;\t"); | 
| 531 | if(foo == NULL){ | 
| 532 | fprintf( stderr, | 
| 533 | "error in reading quaternion 1 from %s\n" | 
| 534 | "natoms  = %d, i = %d\n", | 
| 535 | inName, nAtoms, i ); | 
| 536 | exit(0); | 
| 537 | } | 
| 538 | q[1] = atof( foo ); | 
| 539 |  | 
| 540 | foo = strtok(NULL, " ,;\t"); | 
| 541 | if(foo == NULL){ | 
| 542 | fprintf( stderr, | 
| 543 | "error in reading quaternion 2 from %s\n" | 
| 544 | "natoms  = %d, i = %d\n", | 
| 545 | inName, nAtoms, i ); | 
| 546 | exit(0); | 
| 547 | } | 
| 548 | q[2] = atof( foo ); | 
| 549 |  | 
| 550 | foo = strtok(NULL, " ,;\t"); | 
| 551 | if(foo == NULL){ | 
| 552 | fprintf( stderr, | 
| 553 | "error in reading quaternion 3 from %s\n" | 
| 554 | "natoms  = %d, i = %d\n", | 
| 555 | inName, nAtoms, i ); | 
| 556 | exit(0); | 
| 557 | } | 
| 558 | q[3] = atof( foo ); | 
| 559 |  | 
| 560 | normalizeQ( q ); | 
| 561 |  | 
| 562 | setU( q, atoms[i].u ); | 
| 563 | } | 
| 564 | } | 
| 565 |  | 
| 566 |  | 
| 567 | void setU( double the_q[4], double u[3] ){ | 
| 568 |  | 
| 569 | double q0Sqr, q1Sqr, q2Sqr, q3Sqr; | 
| 570 | double A[3][3]; | 
| 571 | double rb[3]; | 
| 572 | int i,j,k; | 
| 573 |  | 
| 574 | // initialize the axis | 
| 575 |  | 
| 576 | rb[0] = 0.0; | 
| 577 | rb[1] = 0.0; | 
| 578 | rb[2] = 1.0; | 
| 579 |  | 
| 580 | // set the rotation matrix | 
| 581 |  | 
| 582 | q0Sqr = the_q[0] * the_q[0]; | 
| 583 | q1Sqr = the_q[1] * the_q[1]; | 
| 584 | q2Sqr = the_q[2] * the_q[2]; | 
| 585 | q3Sqr = the_q[3] * the_q[3]; | 
| 586 |  | 
| 587 | A[0][0] = q0Sqr + q1Sqr - q2Sqr - q3Sqr; | 
| 588 | A[0][1] = 2.0 * ( the_q[1] * the_q[2] + the_q[0] * the_q[3] ); | 
| 589 | A[0][2] = 2.0 * ( the_q[1] * the_q[3] - the_q[0] * the_q[2] ); | 
| 590 |  | 
| 591 | A[1][0]  = 2.0 * ( the_q[1] * the_q[2] - the_q[0] * the_q[3] ); | 
| 592 | A[1][1] = q0Sqr - q1Sqr + q2Sqr - q3Sqr; | 
| 593 | A[1][2] = 2.0 * ( the_q[2] * the_q[3] + the_q[0] * the_q[1] ); | 
| 594 |  | 
| 595 | A[2][0] = 2.0 * ( the_q[1] * the_q[3] + the_q[0] * the_q[2] ); | 
| 596 | A[2][1] = 2.0 * ( the_q[2] * the_q[3] - the_q[0] * the_q[1] ); | 
| 597 | A[2][2] = q0Sqr - q1Sqr -q2Sqr +q3Sqr; | 
| 598 |  | 
| 599 | // perform the rotation | 
| 600 |  | 
| 601 | for( i=0; i<3; i++ ){ | 
| 602 | u[i] = 0.0; | 
| 603 |  | 
| 604 | for( j=0; j<3; j++ ){ | 
| 605 | u[i] += A[j][i] * rb[j]; | 
| 606 | } | 
| 607 | } | 
| 608 | } | 
| 609 |  | 
| 610 | void normalizeQ( double q[4] ){ | 
| 611 |  | 
| 612 | double qSqr, qL; | 
| 613 | int i; | 
| 614 |  | 
| 615 | qSqr = 0.0; | 
| 616 | for(i=0;i<4;i++) | 
| 617 | qSqr += q[i] * q[i]; | 
| 618 |  | 
| 619 | for(i=0;i<4;i++) | 
| 620 | q[i] /= qSqr; | 
| 621 | } | 
| 622 |  | 
| 623 | double scanSmallestZ(int startFrame){ | 
| 624 | double smallest; | 
| 625 | int i; | 
| 626 |  | 
| 627 | smallest = posArray[startFrame].Hmat[2][2]; | 
| 628 | for(i=startFrame;i<nFrames;i++) | 
| 629 | smallest = | 
| 630 | (smallest<posArray[i].Hmat[2][2])? smallest : posArray[i].Hmat[2][2]; | 
| 631 |  | 
| 632 | return smallest; | 
| 633 | } | 
| 634 |  | 
| 635 | double calcAverageXY( double startTime ){ | 
| 636 |  | 
| 637 | int i, count, startFrame, startFound; | 
| 638 | double avg; | 
| 639 |  | 
| 640 | startFound = 0; | 
| 641 | startFrame = -1; | 
| 642 | while( !startFound ){ | 
| 643 |  | 
| 644 | startFrame++; | 
| 645 |  | 
| 646 | if(startFrame >= nFrames){ | 
| 647 |  | 
| 648 | fprintf( stderr, | 
| 649 | "Start Time, %G, was not found in the dump file.\n", | 
| 650 | startTime ); | 
| 651 | exit(0); | 
| 652 | } | 
| 653 |  | 
| 654 | if(startTime <= frameTimes[startFrame]) | 
| 655 | startFound = 1; | 
| 656 |  | 
| 657 |  | 
| 658 | } | 
| 659 |  | 
| 660 | avg = 0.0; | 
| 661 | count = 0; | 
| 662 | for(i=startFrame;i<nFrames;i++){ | 
| 663 | avg += posArray[i].Hmat[0][0] * posArray[i].Hmat[1][1]; | 
| 664 | count++; | 
| 665 | } | 
| 666 | avg /= (double)count; | 
| 667 |  | 
| 668 | return avg; | 
| 669 | } |