1 /* 2 Teem: Tools to process and visualize scientific data and images . 3 Copyright (C) 2012, 2011, 2010, 2009 University of Chicago 4 Copyright (C) 2008, 2007, 2006, 2005 Gordon Kindlmann 5 Copyright (C) 2004, 2003, 2002, 2001, 2000, 1999, 1998 University of Utah 6 7 This library is free software; you can redistribute it and/or 8 modify it under the terms of the GNU Lesser General Public License 9 (LGPL) as published by the Free Software Foundation; either 10 version 2.1 of the License, or (at your option) any later version. 11 The terms of redistributing and/or modifying this software also 12 include exceptions to the LGPL that facilitate static linking. 13 14 This library is distributed in the hope that it will be useful, 15 but WITHOUT ANY WARRANTY; without even the implied warranty of 16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 17 Lesser General Public License for more details. 18 19 You should have received a copy of the GNU Lesser General Public License 20 along with this library; if not, write to Free Software Foundation, Inc., 21 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA 22 */ 23 24 #ifndef PULL_HAS_BEEN_INCLUDED 25 #define PULL_HAS_BEEN_INCLUDED 26 27 #include <teem/air.h> 28 #include <teem/hest.h> 29 #include <teem/biff.h> 30 #include <teem/nrrd.h> 31 #include <teem/ell.h> 32 #include <teem/gage.h> 33 #include <teem/limn.h> 34 #include <teem/ten.h> 35 36 /* 37 ** This library was created to implement the research and methods of: 38 ** 39 ** Gordon L. Kindlmann, Ra{\'u}l San Jos{\'e} Est{\'e}par, 40 ** Stephen M. Smith, Carl-Fredrik Westin. 41 ** Sampling and Visualizing Creases with Scale-Space Particles. 42 ** IEEE Trans. on Visualization and Computer Graphics, 43 ** 15(6):1415-1424 (2009) 44 ** 45 ** Further information and usage examples: 46 ** http://people.cs.uchicago.edu/~glk/ssp/ 47 ** 48 ** The library is still being actively developed to support research on 49 ** particle systems for general feature sampling. At some point, for 50 ** example, it should subsume the "push" library for glyph packing in 51 ** tensor fields. 52 */ 53 54 #if defined(_WIN32) && !defined(__CYGWIN__) && !defined(TEEM_STATIC) 55 # if defined(TEEM_BUILD) || defined(pull_EXPORTS) || defined(teem_EXPORTS) 56 # define PULL_EXPORT extern __declspec(dllexport) 57 # else 58 # define PULL_EXPORT extern __declspec(dllimport) 59 # endif 60 #else /* TEEM_STATIC || UNIX */ 61 # define PULL_EXPORT extern 62 #endif 63 64 #ifdef __cplusplus 65 extern "C" { 66 #endif 67 68 #define PULL pullBiffKey 69 #define PULL_THREAD_MAXNUM 512 70 #define PULL_VOLUME_MAXNUM 4 71 #define PULL_POINT_NEIGH_INCR 16 72 #define PULL_BIN_MAXNUM 40000000 /* sanity check on max number bins */ 73 #define PULL_PHIST 0 74 #define PULL_HINTER 0 75 #define PULL_TANCOVAR 1 76 77 /* 78 ******** pullInfo enum 79 ** 80 ** all the things that might be *learned* about the local neighborhood 81 ** that are used as ingredients in the computation of particle motion. 82 ** This info was originally learned only from gage, but now (according 83 ** to value of pullSource) can come from other kinds of information. 84 ** 85 ** There are multiple scalars (and associated) derivatives that can 86 ** be used for dynamics: 87 ** - Inside: just for nudging things to stay inside a mask 88 ** - Height: value for computer-vision-y features of ridges, valleys, 89 ** and edges. Setting pullInfoHeight as a constraint does valley 90 ** sampling (flip the sign to get ridges), based on the various 91 ** "tangents" Setting pullInfoHeightLaplacian as a constraint 92 * does zero-crossing edge detection. 93 ** - Isovalue: just for implicit surfaces f=0 94 ** - Strength: some measure of feature strength, with the assumption 95 ** that it can't be analytically differentiated in space or scale. 96 */ 97 enum { 98 pullInfoUnknown, /* 0 */ 99 pullInfoTensor, /* 1: [7] tensor here */ 100 pullInfoTensorInverse, /* 2: [7] inverse of tensor here */ 101 pullInfoHessian, /* 3: [9] hessian used for force distortion */ 102 pullInfoInside, /* 4: [1] containment scalar */ 103 pullInfoInsideGradient, /* 5: [3] containment vector */ 104 pullInfoHeight, /* 6: [1] for gravity, and edge and crease 105 feature detection */ 106 pullInfoHeightGradient, /* 7: [3] */ 107 pullInfoHeightHessian, /* 8: [9] */ 108 pullInfoHeightLaplacian, /* 9: [1] for zero-crossing edge detection */ 109 pullInfoSeedPreThresh, /* 10: [1] scalar for pre-thresholding seeding, 110 so that points can be quickly eliminated 111 (e.g. prior to constraint satisfaction) */ 112 pullInfoSeedThresh, /* 11: [1] scalar for thresholding seeding */ 113 pullInfoLiveThresh, /* 12: [1] scalar for thresholding extant 114 particles, AND for future additions from 115 population control */ 116 pullInfoLiveThresh2, /* 13: [1] another pullInfoLiveThresh */ 117 pullInfoLiveThresh3, /* 14: [1] yet another pullInfoLiveThresh */ 118 pullInfoTangent1, /* 15: [3] first tangent to motion allowed 119 for constraint satisfaction */ 120 pullInfoTangent2, /* 16: [3] second tangent to motion allowed 121 for constraint satisfaction */ 122 pullInfoNegativeTangent1, /* 17: [3] like the tangents, but with a negated 123 objective function */ 124 pullInfoNegativeTangent2, /* 18: [3] second negative tangent */ 125 pullInfoIsovalue, /* 19: [1] for isosurface extraction */ 126 pullInfoIsovalueGradient, /* 20: [3] */ 127 pullInfoIsovalueHessian, /* 21: [9] */ 128 pullInfoStrength, /* 22: [1] */ 129 pullInfoQuality, /* 23: [1] */ 130 pullInfoLast 131 }; 132 #define PULL_INFO_MAX 23 133 134 /* 135 ********* pullProp* enum: the various properties of particles in the system 136 ** 137 ** These are things that are not learned from the image data, but are 138 ** descriptions of a particle, its neighborhood of particles in the 139 ** system, and the current state of system computation. As of revision 140 ** 5080, these may refer to information *computed* from the per-particle 141 ** info, but which may not itself be saved as a field in pullPoint. 142 ** 143 ** consider adding: dot between normalized directions of force and movmt 144 */ 145 enum { 146 pullPropUnknown, /* 0: nobody knows */ 147 pullPropIdtag, /* 1: [1] idtag (unsigned int) */ 148 pullPropIdCC, /* 2: [1] idcc (unsigned int) */ 149 pullPropEnergy, /* 3: [1] energy from last iteration */ 150 pullPropStepEnergy, /* 4: [1] step size for minimizing energy */ 151 pullPropStepConstr, /* 5: [1] step size for constraint satis. */ 152 pullPropStuck, /* 6: [1] how many iters its been stuck */ 153 pullPropPosition, /* 7: [4] position */ 154 pullPropForce, /* 8: [4] force accumulation */ 155 pullPropNeighDistMean, /* 9: [1] "mean distance" to neighbors */ 156 pullPropScale, /* 10: [1] scale position */ 157 pullPropNeighCovar, /* 11: [10] unique coeffs of covariance matrix 158 of offsets to all interacting neighbors, 159 in rs-normalized space, updated only during 160 pullProcessModeNeighLearn. The layout is: 161 0:xx 1:xy 2:xz 3:xs 162 . 4:yy 5:yz 6:ys 163 . 7:zz 8:zs 164 . 9:ss */ 165 pullPropNeighCovar7Ten, /* 12: [7] spatial 3x3 submatrix of covariance, 166 formatted as ten-compatible 7-tensor */ 167 pullPropNeighTanCovar, /* 13: [6] covariance of "tangents" of neighbors, 168 (e.g. pullInfoTangent1 for crease surfaces) 169 including point itself */ 170 pullPropNeighInterNum, /* 14: [1] # neighbors last iter */ 171 pullPropNeighCovarTrace, /* 15: [1] trace of NeighCovar */ 172 pullPropNeighCovarDet, /* 16: [1] det of NeighCovar */ 173 pullPropStability, /* 17: [1] some measure of NeighCovar */ 174 pullPropLast 175 }; 176 #define PULL_PROP_MAX 17 177 178 /* 179 ** the components of a point's status that are set as a bitflag 180 ** in point->status 181 */ 182 enum { 183 pullStatusUnknown, /* 0: nobody knows */ 184 pullStatusStuck, /* 1: couldn't move to decrease energy */ 185 #define PULL_STATUS_STUCK_BIT (1<< 1) 186 pullStatusNewbie, /* 2: not binned, testing if the system 187 would be better with me in it */ 188 #define PULL_STATUS_NEWBIE_BIT (1<< 2) 189 pullStatusNixMe, /* 3: nix me at the *end* of this iter, 190 and don't look at me for energy 191 during this iteration */ 192 #define PULL_STATUS_NIXME_BIT (1<< 3) 193 pullStatusEdge, /* 4: at the spatial edge of one of the 194 volumes: gage had to invent values for 195 some samples in the kernel support */ 196 #define PULL_STATUS_EDGE_BIT (1<< 4) 197 pullStatusLast 198 }; 199 200 /* 201 ******** pullInterType* enum 202 ** 203 ** the different types of scale-space interaction that can happen 204 ** in scale-space. The descriptions here overlook the normalization 205 ** by radiusScale and radiusSpace 206 */ 207 enum { 208 pullInterTypeUnknown, /* 0 */ 209 pullInterTypeJustR, /* 1: phi(r,s) = phi_r(r) */ 210 pullInterTypeUnivariate, /* 2: phi(r,s) = phi_r(u); u = sqrt(r*r+s*s) */ 211 pullInterTypeSeparable, /* 3: phi(r,s) = phi_r(r)*phi_s(s) */ 212 pullInterTypeAdditive, /* 4: phi(r,s) = beta*phi_r(r)*win(s) 213 + (1-beta)*win(r)*phi_s(s) */ 214 pullInterTypeLast 215 }; 216 #define PULL_INTER_TYPE_MAX 4 217 218 /* 219 ******** pullEnergyType* enum 220 ** 221 ** the different shapes of potential energy profiles that can be used 222 */ 223 enum { 224 pullEnergyTypeUnknown, /* 0 */ 225 pullEnergyTypeSpring, /* 1 */ 226 pullEnergyTypeGauss, /* 2 */ 227 pullEnergyTypeBspln, /* 3 */ 228 pullEnergyTypeButterworth, /* 4 */ 229 pullEnergyTypeCotan, /* 5 */ 230 pullEnergyTypeCubic, /* 6 */ 231 pullEnergyTypeQuartic, /* 7 */ 232 pullEnergyTypeCubicWell, /* 8 */ 233 pullEnergyTypeBetterCubicWell, /* 9 */ 234 pullEnergyTypeQuarticWell, /* 10 */ 235 pullEnergyTypeHepticWell, /* 11 */ 236 pullEnergyTypeZero, /* 12 */ 237 pullEnergyTypeButterworthParabola, /* 13 */ 238 pullEnergyTypeLast 239 }; 240 #define PULL_ENERGY_TYPE_MAX 13 241 #define PULL_ENERGY_PARM_NUM 3 242 243 enum { 244 pullProcessModeUnknown, /* 0 */ 245 pullProcessModeDescent, /* 1 */ 246 pullProcessModeNeighLearn, /* 2 */ 247 pullProcessModeAdding, /* 3 */ 248 pullProcessModeNixing, /* 4 */ 249 pullProcessModeLast 250 }; 251 #define PULL_PROCESS_MODE_MAX 4 252 253 /* 254 ** the conditions under which a point may find itself at some position 255 */ 256 enum { 257 pullCondUnknown, /* 0 */ 258 pullCondOld, /* 1 */ 259 pullCondConstraintSatA, /* 2 */ 260 pullCondConstraintSatB, /* 3 */ 261 pullCondEnergyTry, /* 4 */ 262 pullCondConstraintFail, /* 5 */ 263 pullCondEnergyBad, /* 6 */ 264 pullCondNew, /* 7 */ 265 pullCondLast 266 }; 267 268 /* 269 ** the places that "info" can be learned from. Originally this was 270 ** strictly from gage and no where else; it can be handy to allow it 271 ** to originate from different places as well, but still flow through 272 ** the channels now organized around pullInfo 273 */ 274 enum { 275 pullSourceUnknown, /* 0 */ 276 pullSourceGage, /* 1: measured from gage */ 277 pullSourceProp, /* 2: copied from a pullProp */ 278 pullSourceLast 279 }; 280 #define PULL_SOURCE_MAX 2 281 282 /* 283 ** the different kinds of computations and entities that one can 284 ** count, for book-keeping and meta-optimization purposes 285 */ 286 enum { 287 pullCountUnknown, /* 0 */ 288 pullCountDescent, /* 1 */ 289 pullCountTestStep, /* 2 */ 290 pullCountEnergyFromImage, /* 3 */ 291 pullCountForceFromImage, /* 4 */ 292 pullCountEnergyFromPoints, /* 5 */ 293 pullCountForceFromPoints, /* 6 */ 294 pullCountProbe, /* 7 */ 295 pullCountConstraintSatisfy, /* 8 */ 296 pullCountAdding, /* 9 */ 297 pullCountNixing, /* 10 */ 298 pullCountPointsStuck, /* 11 */ 299 pullCountPoints, /* 12 */ 300 pullCountCC, /* 13 */ 301 pullCountIteration, /* 14 */ 302 pullCountLast 303 }; 304 #define PULL_COUNT_MAX 14 305 306 /* 307 ** reasons for pullTraceSet to stop (or go nowhere) 308 */ 309 enum { 310 pullTraceStopUnknown, /* 0 */ 311 pullTraceStopSpeeding, /* 1 */ 312 pullTraceStopConstrFail, /* 2 */ 313 pullTraceStopBounds, /* 3 */ 314 pullTraceStopLength, /* 4 */ 315 pullTraceStopStub, /* 5 */ 316 pullTraceStopLast 317 }; 318 #define PULL_TRACE_STOP_MAX 5 319 320 /* 321 ** Defines how par-particle information can be learned. This is 322 ** typically via measurements in the image by gage, but other sources 323 ** are possible (as indicated by the source field). 324 ** 325 ** There is a basic design decision in whether to make pullInfos the 326 ** most general representation of information (including referring to 327 ** pullProps), or, whether pullProps should be the general thing, 328 ** which can then refer to some pullInfo. In the long run, these 329 ** should probably be merged. In the short term, there is better 330 ** infrastructure for representing and parsing pullInfos, so that has 331 ** been leveraged to include other things as well. 332 */ 333 typedef struct pullInfoSpec_t { 334 /* ------ INPUT ------ */ 335 int info, /* from the pullInfo* enum: what is the 336 "info" that this infospec defines */ 337 source; /* from the pullSource* enum: where does 338 this information come from */ 339 char *volName; /* volume name */ 340 int item, /* which gage item (for pullSourceGage) */ 341 prop; /* which pull property (for pullSourceProp) */ 342 double scale, /* scaling factor (including sign) */ 343 zero; /* for height and inside: where is zero, 344 for seedThresh, threshold value */ 345 int constraint; /* (for scalar items) minimizing this 346 is a constraint to enforce per-point 347 per-iteration, not merely a contribution 348 to the point's energy */ 349 /* ------ INTERNAL ------ */ 350 unsigned int volIdx; /* which volume (for pullSourceGage) */ 351 } pullInfoSpec; 352 353 /* 354 ******** pullPoint 355 ** 356 */ 357 typedef struct pullPoint_t { 358 unsigned int idtag, /* unique point ID */ 359 idCC; /* id for connected component analysis */ 360 struct pullPoint_t **neighPoint; /* list of neighboring points */ 361 unsigned int neighPointNum; 362 airArray *neighPointArr; /* airArray around neighPoint and neighNum 363 (no callbacks used here) */ 364 double neighDistMean; /* average of distance to neighboring 365 points with whom this point interacted, 366 in rs-normalized space */ 367 float neighCovar[10], /* unique coeffs in 4x4 covariance matrix of 368 neighbors with whom this point interacted */ 369 #if PULL_TANCOVAR 370 neighTanCovar[6], /* covariance of "tangent" info of neighbors */ 371 #endif 372 stability; /* the scalar stability measure */ 373 unsigned int neighInterNum, /* number of particles with which I had some 374 non-zero interaction on last iteration */ 375 stuckIterNum; /* how many iterations I've been stuck */ 376 #if PULL_PHIST 377 double *phist; /* history of positions tried in the last iter, 378 in sets of 5 doubles: (x,y,z,t,info) */ 379 unsigned int phistNum; /* number of positions stored */ 380 airArray *phistArr; /* airArray around phist */ 381 #endif 382 int status; /* bit-flag of status info */ 383 double pos[4], /* position in space and scale */ 384 energy, /* energy accumulator for this iteration */ 385 force[4], /* force accumulator for this iteration */ 386 stepEnergy, /* step size for energy minimization */ 387 stepConstr, /* step size for constraint satisfaction; 388 HEY: this doens't seem to be really used? */ 389 info[1]; /* all information learned from gage that matters 390 for particle dynamics. This is sneakily 391 allocated for *more*, depending on needs, 392 so this MUST be last field */ 393 } pullPoint; 394 395 /* 396 ******** pullBin 397 ** 398 ** the data structure for doing spatial binning. 399 */ 400 typedef struct pullBin_t { 401 pullPoint **point; /* dyn. alloc. array of point pointers */ 402 unsigned int pointNum; /* # of points in this bin */ 403 airArray *pointArr; /* airArray around point and pointNum 404 (no callbacks used here) */ 405 struct pullBin_t **neighBin; /* NULL-terminated list of all 406 neighboring bins, including myself */ 407 } pullBin; 408 409 /* 410 ******** pullEnergy 411 ** 412 ** the functions which determine inter-point forces 413 ** 414 ** NOTE: the eval() function probably does NOT check to see it was passed 415 ** a non-NULL pointer into which to store the derivative of energy ("denr") 416 ** 417 ** Thu Apr 10 12:40:08 EDT 2008: nixed the "support" function, since it 418 ** was annoying to deal with variable support potentials. Now everything 419 ** cuts off at dist=1. You can still use the parm vector to change the 420 ** shape inside the support. 421 */ 422 typedef struct { 423 char name[AIR_STRLEN_SMALL]; 424 unsigned int parmNum; 425 double (*well)(double *wx, const double parm[PULL_ENERGY_PARM_NUM]); 426 double (*eval)(double *denr, double dist, 427 const double parm[PULL_ENERGY_PARM_NUM]); 428 } pullEnergy; 429 430 typedef struct { 431 const pullEnergy *energy; 432 double parm[PULL_ENERGY_PARM_NUM]; 433 } pullEnergySpec; 434 435 /* 436 ** In the interests of simplicity (and with the cost of some redundancy), 437 ** this is going to copied per-task, which is why it contains the gageContext 438 ** The idea is that the first of these is somehow set up by the user 439 ** or something, and the rest of them are created within pull per-task. 440 */ 441 typedef struct { 442 int verbose; /* blah blah blah */ 443 char *name; /* how the volume will be identified 444 (like its a variable name) */ 445 const gageKind *kind; 446 const Nrrd *ninSingle; /* don't own */ 447 const Nrrd *const *ninScale; /* don't own; 448 NOTE: only one of ninSingle and ninScale 449 can be non-NULL */ 450 unsigned int scaleNum; /* number of scale-space samples (volumes) */ 451 double *scalePos; /* location of all samples in scale */ 452 int scaleDerivNorm; /* normalize derivatives based on scale */ 453 double scaleDerivNormBias; /* bias on derivative normalization by scale */ 454 NrrdKernelSpec *ksp00, /* for sampling tensor field */ 455 *ksp11, /* for gradient of mask, other 1st derivs */ 456 *ksp22, /* for 2nd derivatives */ 457 *kspSS; /* for reconstructing from scale-space 458 samples */ 459 gageQuery pullValQuery; /* if this is a pullValGageKind volume, 460 then we don't have a real gageContext, 461 and we have to manage our own query */ 462 gageContext *gctx; /* do own, and set based on info here */ 463 gagePerVolume *gpvl, /* stupid gage API . . . */ 464 **gpvlSS; /* stupid gage API . . . */ 465 int seedOnly, /* volume only required for seeding, for 466 either pullInfoSeedThresh or 467 pullInfoSeedPreThresh */ 468 forSeedPreThresh; /* we learn pullInfoSeedPreThresh from this */ 469 } pullVolume; 470 471 /* 472 ******** pullTask 473 ** 474 ** The information specific for a thread. 475 */ 476 typedef struct pullTask_t { 477 struct pullContext_t 478 *pctx; /* parent's context; not const because the 479 tasks assign themselves bins to do work */ 480 pullVolume 481 *vol[PULL_VOLUME_MAXNUM]; /* volumes copied from parent */ 482 const double 483 *ans[PULL_INFO_MAX+1]; /* answer *pointers* for all possible infos, 484 pointing into per-task per-volume gctxs 485 (or into above per-task pullValAnswer), 486 OR: NULL if that info is not being used */ 487 int processMode, /* what kind of point processing is being 488 done by this task right now */ 489 probeSeedPreThreshOnly; /* hack-ish flag to communicate to pullProbe 490 that we only care about SeedPreThresh */ 491 airThread *thread; /* my thread */ 492 unsigned int threadIdx; /* which thread am I */ 493 airRandMTState *rng; /* state for my RNG */ 494 pullPoint *pointBuffer, /* place for copying point into during 495 strength ascent computation; can't be 496 statically allocated because pullPoint 497 size is known only at run-time */ 498 **neighPoint; /* array of point pointers, either all 499 possible points from neighbor bins, or 500 last learned interacting neighbors */ 501 pullPoint **addPoint; /* points to add before next iter */ 502 unsigned int addPointNum; /* # of points to add */ 503 airArray *addPointArr; /* airArray around addPoint, addPointNum */ 504 pullPoint **nixPoint; /* points to nix before next iter */ 505 unsigned int nixPointNum; /* # of points to nix */ 506 airArray *nixPointArr; /* airArray around nixPoint, nixPointNum */ 507 void *returnPtr; /* for airThreadJoin */ 508 unsigned int stuckNum; /* # stuck particles seen by this task */ 509 } pullTask; 510 511 /* 512 ******** pullInitMethod* enum 513 ** 514 ** the different ways pull can be initialized 515 */ 516 enum { 517 pullInitMethodUnknown, /* 0 */ 518 pullInitMethodRandom, /* 1 */ 519 pullInitMethodHalton, /* 2 */ 520 pullInitMethodPointPerVoxel, /* 3 */ 521 pullInitMethodGivenPos, /* 4 */ 522 pullInitMethodLast 523 }; 524 525 /* 526 ******** pullInitParm 527 ** 528 ** none of this is directly user-set; set with pullInit*Set function 529 ** (note that there is no pullInit* enum; these values are too diverse) 530 */ 531 typedef struct { 532 int method; /* from pullInitMethod* enum */ 533 int liveThreshUse, /* use the liveThresh info as a criterion 534 for seeding */ 535 unequalShapesAllow; /* when taking in volumes, allow their shapes 536 to be unuequal. Some confusing usage problems 537 are due to using differently shaped volumes */ 538 double jitter; /* w/ PointPerVoxel, 539 how much to jitter index space positions */ 540 unsigned int numInitial, /* w/ Random OR Halton, # points to start with */ 541 haltonStartIndex, /* w/ Halton, first index to use */ 542 samplesAlongScaleNum, /* w/ PointPerVoxel, 543 # of samples along scale (distributed 544 uniformly in scale's *index* space*/ 545 ppvZRange[2]; /* (hack to permit seeding only in part of 546 volume, when initialization is painfully 547 the main bottleneck) w/ PointPerVoxel, 548 range of indices along Z to do seeding 549 by pointPerVoxel, or, {1,0} to do the 550 whole volume as normal */ 551 int pointPerVoxel; /* w/ PointPerVoxel, 552 if this (ppv) is > 0, then use 553 ppv points per voxel. If ppv < 0, then jitter 554 point point every -ppv'th voxel 555 (so ppv=-1 is same as ppv=1) */ 556 const Nrrd *npos; /* positions (4xN array) to start with */ 557 } pullInitParm; 558 559 /* 560 ******** pullIterParm* enum 561 ** 562 ** parameters related to iterations and their periods 563 */ 564 enum { 565 pullIterParmUnknown, 566 567 /* if non-zero, minimum number of iterations for whole system. */ 568 pullIterParmMin, 569 570 /* if non-zero, max number of iterations for whole system. 571 if zero: no explicit limit on the number of iterations */ 572 pullIterParmMax, 573 574 /* if non-zero, max number of iterations we allow something to be 575 continuously stuck before nixing it */ 576 pullIterParmStuckMax, 577 578 /* if non-zero, max number of iterations for enforcing each constraint */ 579 pullIterParmConstraintMax, 580 581 /* how many intervals to wait between attemps at population control, 582 or, 0 to say: "no pop cntl" */ 583 pullIterParmPopCntlPeriod, 584 585 /* how many iterations for which to run descent on tentative new 586 points during population control, so that they end up at a good 587 location, at which we can meaningfully test whether adding the 588 point will reducing system energy */ 589 pullIterParmAddDescent, 590 591 /* periodicity with which to call the pullContext->iter_cb */ 592 pullIterParmCallback, 593 594 /* if non-zero, interval between iters at which output snapshots are saved */ 595 pullIterParmSnap, 596 597 /* the half-life of energyIncreasePermit, in terms of iterations, or 598 0 if no such decay is wanted */ 599 pullIterParmEnergyIncreasePermitHalfLife, 600 601 pullIterParmLast 602 }; 603 604 typedef struct { 605 unsigned int min, max, 606 popCntlPeriod, 607 addDescent, 608 constraintMax, 609 stuckMax, 610 callback, 611 snap, 612 energyIncreasePermitHalfLife; 613 } pullIterParm; 614 615 /* 616 ******** pullSysParm* enum 617 ** 618 ** various continuous parameters. It should be understood that the number 619 ** of these is a reflection of the experimentation and exploration that 620 ** went into the creation of these particles systems, rather than a need 621 ** for onerous parameter tweaking in order to get something useful done. 622 ** Except for alpha, beta, and gamma, these have reasonable defaults. 623 */ 624 enum { 625 pullSysParmUnknown, 626 627 /* balance between particle-image energy E_i and inter-particle energy E_ij, 628 the most important parameter in the governing equation of the system. 629 alpha = 0: only particle-image; alpha = 1: only inter-particle */ 630 pullSysParmAlpha, 631 632 /* when using pullInterAdditive ("Phi_2") inter-particle energy, beta 633 sets the balance between spatial repulsion and scale attraction. 634 (if not using this energy, this is moot) 635 beta = 0: only spatial repulsion; beta = 1: only scale attraction */ 636 pullSysParmBeta, 637 638 /* when energyFromStrength is non-zero: scaling factor on energy 639 from strength */ 640 pullSysParmGamma, 641 642 /* when learning a suitable gamma from data (via pullGammaLearn) 643 in the context of separable energy functions, how much to 644 scale the learned gamma. Making this greater than 1.0 helps 645 let the energy from strength overpower the slight repulsion 646 along scale that might exist by the separable phi construction. */ 647 pullSysParmSeparableGammaLearnRescale, 648 649 /* to be more selective for pullInfoSeedThresh and 650 pullInfoLiveThresh at higher scales, set theta > 0, and the 651 effective threshold will be base threshold + theta*scale. 652 HOWEVER, the way this is implemented is a hack: 653 the probed strength value is decremented by theta*scale */ 654 pullSysParmTheta, 655 656 /* initial (time) step for dynamics */ 657 pullSysParmStepInitial, 658 659 /* radius/scaling of inter-particle interactions in the spatial domain */ 660 pullSysParmRadiusSpace, 661 662 /* radius/scaling of inter-particle interactions in the scale domain */ 663 pullSysParmRadiusScale, 664 665 /* spatial width of bin, as multiple of pullSysParmRadiusSpace (the 666 width along scale is set to pullSysParmRadiusScale). Can't be 667 lower than 1, but may be usefully set greater that 1 to reduce 668 the total number of bins, especially if caching neighbor lists */ 669 pullSysParmBinWidthSpace, 670 671 /* probability that we find the true neighbors of the particle, as 672 opposed to using a cached list */ 673 pullSysParmNeighborTrueProb, 674 675 /* probability that we do image probing to find out what's really going on */ 676 pullSysParmProbeProb, 677 678 /* (>= 1.0) how much to opportunistically scale up step size (for 679 energy minimization) with every iteration */ 680 pullSysParmOpporStepScale, 681 682 /* (< 1.0) when energy goes up instead of down, or when constraint 683 satisfaction seems to be going the wrong way, how to scale (down) 684 step size */ 685 pullSysParmBackStepScale, 686 687 /* pseudo-convergence threshold that controls when population control is 688 activated (has to be higher than (less strict) energyDecreaseMin */ 689 pullSysParmEnergyDecreasePopCntlMin, 690 691 /* epsilon amount by which its okay for particle energy to increase, 692 in the context of gradient descent */ 693 pullSysParmEnergyIncreasePermit, 694 695 /* convergence threshold: stop when fractional improvement 696 (decrease) in total system energy dips below this */ 697 pullSysParmEnergyDecreaseMin, 698 699 /* convergence threshold for constraint satisfaction: finished if 700 stepsize goes below this times pctx->voxelSizeSpace */ 701 pullSysParmConstraintStepMin, 702 703 /* spring constant on bbox wall */ 704 pullSysParmWall, 705 706 /* when doing population control nixing, don't nix a particle if this 707 fraction of its neighbors have already been nixed (Section 3.5 of 708 paper implies that this value should be 0.5; lower values also work) */ 709 pullSysParmFracNeighNixedMax, 710 711 pullSysParmLast 712 }; 713 714 typedef struct { 715 double alpha, beta, gamma, separableGammaLearnRescale, theta, wall, 716 radiusSpace, radiusScale, binWidthSpace, 717 neighborTrueProb, probeProb, 718 stepInitial, opporStepScale, backStepScale, constraintStepMin, 719 energyDecreaseMin, 720 energyDecreasePopCntlMin, 721 energyIncreasePermit, 722 fracNeighNixedMax; 723 } pullSysParm; 724 725 /* 726 ******** pullConstraintFail* enum 727 ** 728 ** the various ways constriant satisfaction can fail 729 */ 730 enum { 731 pullConstraintFailUnknown, /* 0 */ 732 pullConstraintFailProjGradZeroA, /* 1 */ 733 pullConstraintFailProjGradZeroB, /* 2 */ 734 pullConstraintFailIterMaxed, /* 3 */ 735 pullConstraintFailTravel, /* 4 */ 736 pullConstraintFailLast 737 }; 738 #define PULL_CONSTRAINT_FAIL_MAX 4 739 740 /* 741 ******** pullFlag* enum 742 ** 743 ** the various booleans 744 */ 745 enum { 746 pullFlagUnknown, 747 748 /* permute points during rebinning between iters, to randomize their 749 ordering within their bin (bins are still processed in scanline 750 order) */ 751 pullFlagPermuteOnRebin, 752 753 /* when alpha is exactly zero, probably with the purpose of migrating 754 particles along scale towards the scale of maximal strength while 755 avoiding all inter-particle interaction, don't try any population 756 control */ 757 pullFlagNoPopCntlWithZeroAlpha, 758 759 /* when learning gamma via Eq. 26 of the paper, also use beta to 760 more accurately reflect the 2nd derivative of energy wrt scale */ 761 pullFlagUseBetaForGammaLearn, 762 763 /* whether or not to deny adding points to bins where there are 764 close points already */ 765 pullFlagRestrictiveAddToBins, 766 767 /* if non-zero, strength is a particle-image energy term that is 768 minimized by motion along scale, which in turn requires extra 769 probing to determine the strength gradient along scale. */ 770 pullFlagEnergyFromStrength, 771 772 /* if non-zero, nix points that got near enough to the volume edge 773 that gage had to invent values for the kernel support */ 774 pullFlagNixAtVolumeEdgeSpace, 775 776 /* if non-zero, during initialization, try constraint satisfaction 777 (if there is a constraint) before testing whether the seedThresh 778 is met. Doing the constraint will take longer, but a point is 779 more likely to meet a threshold based on feature strength */ 780 pullFlagConstraintBeforeSeedThresh, 781 782 /* do no adding during population control */ 783 pullFlagNoAdd, 784 785 /* use the targetDim-based "enough" heuristic to guess whether 786 adding a point could usefully reduce system energy. On by 787 default; turn this off when using a large-support energy that 788 involves more neighbors than the single-neighbor-deep energies 789 that are normally used (like cwell:0.66,x or qwell:0.64) */ 790 pullFlagPopCntlEnoughTest, 791 792 /* no binning: all particles can potentially interact (for debugging) */ 793 pullFlagBinSingle, 794 795 /* whether to allow codimension-3 (point) constraint manifolds. 796 Typical uses of constraints are for extracting lines and surfaces */ 797 pullFlagAllowCodimension3Constraints, 798 799 /* what we call scale is not sigma but rather the tau of gage */ 800 pullFlagScaleIsTau, 801 802 /* pullStart should skip initializing the points */ 803 pullFlagStartSkipsPoints, 804 805 pullFlagLast 806 }; 807 808 typedef struct { 809 int permuteOnRebin, 810 noPopCntlWithZeroAlpha, 811 useBetaForGammaLearn, 812 restrictiveAddToBins, 813 energyFromStrength, 814 nixAtVolumeEdgeSpace, 815 constraintBeforeSeedThresh, 816 popCntlEnoughTest, 817 noAdd, 818 binSingle, 819 allowCodimension3Constraints, 820 scaleIsTau, 821 startSkipsPoints; 822 } pullFlag; 823 824 /* 825 ******** pullContext 826 ** 827 ** everything for doing one computation 828 ** 829 */ 830 typedef struct pullContext_t { 831 /* INPUT ----------------------------- */ 832 833 pullInitParm initParm; /* parms for initialization, set with 834 the pullInit*Set() functions */ 835 pullIterParm iterParm; /* parms about iterations and periods, set 836 with pullIterParmSet() */ 837 pullSysParm sysParm; /* continuous parameters for system, 838 set with pullSysParmSet() */ 839 pullFlag flag; /* all flags, set with pullFlagSet() */ 840 int verbose; /* verbosity level, set with 841 pullVerboseSet() */ 842 unsigned int threadNum, /* number of threads to use, set with 843 pullThreadNumSet() */ 844 rngSeed, /* seed value for random number generator, 845 NOT directly related to seed point 846 placement (as name might suggest), 847 set with pullRngSeedSet() */ 848 progressBinMod; /* if non-zero, progress indication by 849 printing "." is given when the bin index 850 is a multiple of this; higher numbers 851 give less feedback, set with 852 pullProgressBinModSet() */ 853 void (*iter_cb)(void *data_cb); /* callback to call from pullRun() every 854 iterParm.callback iterations. This and 855 data_cb are set w/ pullCallbackSet() */ 856 void *data_cb; /* data to pass to callback */ 857 pullVolume 858 *vol[PULL_VOLUME_MAXNUM]; /* the volumes we analyze (we DO OWN), 859 set by either pullVolumeSingleAdd() 860 or pullVolumeStackAdd() */ 861 unsigned int volNum; /* actual length of vol[] used */ 862 pullInfoSpec 863 *ispec[PULL_INFO_MAX+1]; /* info ii is in effect if ispec[ii] is 864 non-NULL (and we DO OWN ispec[ii]), 865 set by pullInfoSpecAdd() */ 866 int interType; /* from the pullInterType* enum. This 867 and the energy specs below are set by 868 pullInterEnergySet() */ 869 pullEnergySpec *energySpecR, /* starting point for radial potential 870 energy function, phi_r */ 871 *energySpecS, /* second energy potential function, for 872 scale-space behavior, phi_s */ 873 *energySpecWin; /* function used to window phi_r along s, 874 and phi_s along r, for use with 875 pullInterTypeAdditive */ 876 877 /* INTERNAL -------------------------- */ 878 879 unsigned int haltonOffset; /* with pullInitMethodHalton, add this to 880 the index to sequence generation, to 881 account for the points previously 882 generated (which did not meet the 883 constraint satisfaction) */ 884 double bboxMin[4], bboxMax[4]; /* scale-space bounding box of all volumes: 885 region over which binning is defined. 886 In 3-D space, the bbox is axis aligned, 887 even when the volume is not so aligned, 888 which means that some bins might be 889 under- or un- utilized, oh well. 890 bboxMin[3] and bboxMax[3] are the 891 bounds of the volume in *scale* (sigma), 892 not t, or tau */ 893 unsigned int infoTotalLen, /* total length of the info buffers needed, 894 which determines size of allocated 895 binPoint */ 896 infoIdx[PULL_INFO_MAX+1], /* index of answer within pullPoint->info */ 897 idtagNext; /* next per-point igtag value */ 898 int haveScale, /* non-zero iff one of the volumes is in 899 scale-space */ 900 constraint, /* if non-zero, we have a constraint to 901 satisfy, and this is its info number */ 902 constraintDim, /* dimension of *spatial* constraint 903 manifold we're working on; or 904 -1 if unknown/unset */ 905 targetDim, /* dimension of total constraint manifold 906 which can be different than constraintDim 907 because of scale-space, and either 908 repulsive (+1) or attractive (+0) 909 behavior along scale; or 910 -1 if unknown/unset */ 911 finished; /* used to signal all threads to return */ 912 double maxDistSpace, /* max dist of point-point interaction in 913 the spatial axes.*/ 914 maxDistScale, /* max dist of point-point interaction 915 along scale */ 916 voxelSizeSpace, /* mean spatial voxel edge length, for 917 limiting travel distance for descent 918 and constraint satisfaction */ 919 voxelSizeScale, /* mean voxel edge length in space, for 920 limiting travel (along scale) distance 921 during descent */ 922 eipScale; /* how to scale energyIncreasePermit 923 at each iteration, in accordance with 924 energyIncreasePermitHalfLife */ 925 pullBin *bin; /* volume of bins (see binsEdge, binNum) */ 926 unsigned int binsEdge[4], /* # bins along each volume edge, 927 determined by maxEval and scale */ 928 binNum, /* total # bins in grid */ 929 binNextIdx; /* next bin of points to be processed, 930 we're done when binNextIdx == binNum */ 931 unsigned int *tmpPointPerm; /* storing points during rebinning */ 932 pullPoint **tmpPointPtr; 933 unsigned int tmpPointNum; 934 935 airThreadMutex *binMutex; /* mutex around bin, needed because bins 936 are the unit of work for the tasks */ 937 pullTask **task; /* dynamically allocated array of tasks */ 938 airThreadBarrier *iterBarrierA; /* barriers between iterations */ 939 airThreadBarrier *iterBarrierB; /* barriers between iterations */ 940 #if PULL_HINTER 941 Nrrd *nhinter; /* 2-D histogram of (r,s)-space relative 942 locations of interacting particles 943 (NOT thread safe) */ 944 #endif 945 FILE *logAdd; /* text-file record of all the particles 946 that have been added 947 (NOT thread-safe) */ 948 949 /* OUTPUT ---------------------------- */ 950 951 double timeIteration, /* time needed for last (single) iter */ 952 timeRun, /* total time spent in pullRun() */ 953 energy; /* final energy of system */ 954 unsigned int addNum, /* # prtls added by PopCntl in last iter */ 955 nixNum, /* # prtls nixed by PopCntl in last iter */ 956 stuckNum, /* # stuck particles in last iter */ 957 pointNum, /* total # particles */ 958 CCNum, /* # connected components */ 959 iter, /* how many iterations were needed */ 960 /* HEY: this should really be per-task, to be thread-safe!! */ 961 count[PULL_COUNT_MAX+1]; /* all possible kinds of counts */ 962 } pullContext; 963 964 /* 965 ******** pullTrace 966 */ 967 typedef struct { 968 double seedPos[4]; /* where was the seed point */ 969 /* ------- output ------- */ 970 Nrrd *nvert, /* locations of tract vertices */ 971 *nstrn, /* if non-NULL, 1-D array of strengths */ 972 *nvelo; /* 1-D list of velocities */ 973 unsigned int seedIdx; /* which index in nvert is for seedpoint */ 974 int whyStop[2], /* why backward/forward (0/1) tracing stopped 975 (from pullTraceStop* enum) */ 976 whyNowhere; /* why trace never started (from pullTraceStop*) */ 977 } pullTrace; 978 979 /* 980 ******** pullTraceMulti 981 */ 982 typedef struct { 983 pullTrace **trace; 984 unsigned int traceNum; 985 airArray *traceArr; 986 } pullTraceMulti; 987 988 /* 989 ******** pullPtrPtrUnion 990 ** 991 ** deal with "dereferencing type-punned pointer will 992 ** break strict-aliasing rules" 993 */ 994 typedef union { 995 pullPoint ***points; /* address of array of point pointers */ 996 void **v; 997 } pullPtrPtrUnion; 998 999 /* defaultsPull.c */ 1000 PULL_EXPORT const int pullPresent; 1001 PULL_EXPORT const int pullPhistEnabled; 1002 PULL_EXPORT const char *pullBiffKey; 1003 1004 /* initPull.c */ 1005 PULL_EXPORT int pullInitRandomSet(pullContext *pctx, unsigned int numInitial); 1006 PULL_EXPORT int pullInitHaltonSet(pullContext *pctx, unsigned int numInitial, 1007 unsigned int start); 1008 PULL_EXPORT int pullInitPointPerVoxelSet(pullContext *pctx, int pointPerVoxel, 1009 unsigned int zSlcMin, 1010 unsigned int zSlcMax, 1011 unsigned int alongScaleNum, 1012 double jitter); 1013 PULL_EXPORT int pullInitGivenPosSet(pullContext *pctx, const Nrrd *npos); 1014 PULL_EXPORT int pullInitLiveThreshUseSet(pullContext *pctx, int liveThreshUse); 1015 PULL_EXPORT int pullInitUnequalShapesAllowSet(pullContext *pctx, int allow); 1016 1017 /* parmPull.c */ 1018 PULL_EXPORT int pullIterParmSet(pullContext *pctx, int which, 1019 unsigned int pval); 1020 PULL_EXPORT int pullSysParmSet(pullContext *pctx, int which, 1021 double pval); 1022 PULL_EXPORT int pullFlagSet(pullContext *pctx, int which, int flag); 1023 PULL_EXPORT int pullVerboseSet(pullContext *pctx, int verbose); 1024 PULL_EXPORT int pullThreadNumSet(pullContext *pctx, unsigned int threadNum); 1025 PULL_EXPORT int pullRngSeedSet(pullContext *pctx, unsigned int rngSeed); 1026 PULL_EXPORT int pullProgressBinModSet(pullContext *pctx, unsigned int bmod); 1027 PULL_EXPORT int pullCallbackSet(pullContext *pctx, 1028 void (*iter_cb)(void *data_cb), 1029 void *data_cb); 1030 PULL_EXPORT int pullInterEnergySet(pullContext *pctx, int interType, 1031 const pullEnergySpec *enspR, 1032 const pullEnergySpec *enspS, 1033 const pullEnergySpec *enspWin); 1034 PULL_EXPORT int pullLogAddSet(pullContext *pctx, FILE *log); 1035 1036 /* energy.c */ 1037 PULL_EXPORT const airEnum *const pullInterType; 1038 PULL_EXPORT const airEnum *const pullEnergyType; 1039 PULL_EXPORT const pullEnergy *const pullEnergyUnknown; 1040 PULL_EXPORT const pullEnergy *const pullEnergySpring; 1041 PULL_EXPORT const pullEnergy *const pullEnergyGauss; 1042 PULL_EXPORT const pullEnergy *const pullEnergyBspln; 1043 PULL_EXPORT const pullEnergy *const pullEnergyButterworth; 1044 PULL_EXPORT const pullEnergy *const pullEnergyCotan; 1045 PULL_EXPORT const pullEnergy *const pullEnergyCubic; 1046 PULL_EXPORT const pullEnergy *const pullEnergyQuartic; 1047 PULL_EXPORT const pullEnergy *const pullEnergyCubicWell; 1048 PULL_EXPORT const pullEnergy *const pullEnergyBetterCubicWell; 1049 PULL_EXPORT const pullEnergy *const pullEnergyQuarticWell; 1050 PULL_EXPORT const pullEnergy *const pullEnergyHepticWell; 1051 PULL_EXPORT const pullEnergy *const pullEnergyZero; 1052 PULL_EXPORT const pullEnergy *const pullEnergyButterworthParabola; 1053 PULL_EXPORT const pullEnergy *const pullEnergyAll[PULL_ENERGY_TYPE_MAX+1]; 1054 PULL_EXPORT pullEnergySpec *pullEnergySpecNew(void); 1055 PULL_EXPORT void pullEnergySpecSet(pullEnergySpec *ensp, 1056 const pullEnergy *energy, 1057 const double parm[PULL_ENERGY_PARM_NUM]); 1058 PULL_EXPORT void pullEnergySpecCopy(pullEnergySpec *esDst, 1059 const pullEnergySpec *esSrc); 1060 PULL_EXPORT pullEnergySpec *pullEnergySpecNix(pullEnergySpec *ensp); 1061 PULL_EXPORT int pullEnergySpecParse(pullEnergySpec *ensp, const char *str); 1062 PULL_EXPORT hestCB *pullHestEnergySpec; 1063 1064 /* volumePull.c */ 1065 PULL_EXPORT pullVolume *pullVolumeNew(void); 1066 PULL_EXPORT pullVolume *pullVolumeNix(pullVolume *vol); 1067 PULL_EXPORT int pullVolumeSingleAdd(pullContext *pctx, 1068 const gageKind *kind, 1069 char *name, const Nrrd *nin, 1070 const NrrdKernelSpec *ksp00, 1071 const NrrdKernelSpec *ksp11, 1072 const NrrdKernelSpec *ksp22); 1073 PULL_EXPORT int pullVolumeStackAdd(pullContext *pctx, 1074 const gageKind *kind, 1075 char *name, 1076 const Nrrd *nin, 1077 const Nrrd *const *ninSS, 1078 double *scalePos, 1079 unsigned int ninNum, 1080 int scaleDerivNorm, 1081 double scaleDerivNormBias, 1082 const NrrdKernelSpec *ksp00, 1083 const NrrdKernelSpec *ksp11, 1084 const NrrdKernelSpec *ksp22, 1085 const NrrdKernelSpec *kspSS); 1086 PULL_EXPORT const pullVolume *pullVolumeLookup(const pullContext *pctx, 1087 const char *volName); 1088 PULL_EXPORT int pullConstraintScaleRange(pullContext *pctx, 1089 double ssrange[2]); 1090 1091 /* enumsPull.c */ 1092 PULL_EXPORT const airEnum *const pullInfo; 1093 PULL_EXPORT const airEnum *const pullSource; 1094 PULL_EXPORT const airEnum *const pullProp; 1095 PULL_EXPORT const airEnum *const pullProcessMode; 1096 PULL_EXPORT const airEnum *const pullTraceStop; 1097 PULL_EXPORT const airEnum *const pullCount; 1098 PULL_EXPORT const airEnum *const pullConstraintFail; 1099 1100 /* infoPull.c */ 1101 PULL_EXPORT unsigned int pullPropLen(int prop); 1102 PULL_EXPORT unsigned int pullInfoLen(int info); 1103 PULL_EXPORT pullInfoSpec *pullInfoSpecNew(void); 1104 PULL_EXPORT pullInfoSpec *pullInfoSpecNix(pullInfoSpec *ispec); 1105 PULL_EXPORT int pullInfoSpecAdd(pullContext *pctx, pullInfoSpec *ispec); 1106 PULL_EXPORT int pullInfoGet(Nrrd *ninfo, int info, pullContext *pctx); 1107 PULL_EXPORT int pullInfoSpecSprint(char str[AIR_STRLEN_LARGE], 1108 const pullContext *pctx, 1109 const pullInfoSpec *ispec); 1110 1111 /* contextPull.c */ 1112 PULL_EXPORT pullContext *pullContextNew(void); 1113 PULL_EXPORT pullContext *pullContextNix(pullContext *pctx); 1114 PULL_EXPORT int pullOutputGet(Nrrd *nPosOut, Nrrd *nTenOut, 1115 Nrrd *nStrengthOut, 1116 const double scaleVec[3], double scaleRad, 1117 pullContext *pctx); 1118 PULL_EXPORT int pullOutputGetFilter(Nrrd *nPosOut, Nrrd *nTenOut, 1119 Nrrd *nStrengthOut, 1120 const double scaleVec[3], double scaleRad, 1121 pullContext *pctx, 1122 unsigned int idtagMin, 1123 unsigned int idtagMax); 1124 PULL_EXPORT int pullPositionHistoryGet(limnPolyData *pld, pullContext *pctx); 1125 PULL_EXPORT int pullPropGet(Nrrd *nprop, int prop, pullContext *pctx); 1126 1127 /* pointPull.c */ 1128 PULL_EXPORT int pullPointInitializePerVoxel(const pullContext *pctx, 1129 const unsigned int pointIdx, 1130 pullPoint *point, 1131 pullVolume *scaleVol, 1132 int *createFailP); 1133 PULL_EXPORT int pullPointInitializeRandomOrHalton(pullContext *pctx, 1134 const unsigned int pointIdx, 1135 pullPoint *point, 1136 pullVolume *scaleVol); 1137 PULL_EXPORT int pullPointInitializeGivenPos(pullContext *pctx, 1138 const double *posData, 1139 const unsigned int pointIdx, 1140 pullPoint *point, 1141 int *createFailP); 1142 PULL_EXPORT double pullPointScalar(const pullContext *pctx, 1143 const pullPoint *point, int sclInfo, 1144 double grad[4], double hess[9]); 1145 PULL_EXPORT unsigned int pullPointNumber(const pullContext *pctx); 1146 PULL_EXPORT unsigned int pullPointNumberFilter(const pullContext *pctx, 1147 unsigned int idtagMin, 1148 unsigned int idtagMax); 1149 PULL_EXPORT pullPoint *pullPointNew(pullContext *pctx); 1150 PULL_EXPORT pullPoint *pullPointNix(pullPoint *pnt); 1151 PULL_EXPORT int pullProbe(pullTask *task, pullPoint *point); 1152 1153 /* binningPull.c */ 1154 PULL_EXPORT int pullBinsPointAdd(pullContext *pctx, pullPoint *point, 1155 /* output */ 1156 pullBin **binUsed); 1157 PULL_EXPORT int pullBinsPointMaybeAdd(pullContext *pctx, pullPoint *point, 1158 /* output */ 1159 pullBin **binUsed, int *added); 1160 1161 /* trace.c */ 1162 PULL_EXPORT pullTrace *pullTraceNew(void); 1163 PULL_EXPORT pullTrace *pullTraceNix(pullTrace *pts); 1164 PULL_EXPORT size_t pullTraceMultiSizeof(const pullTraceMulti *mtrc); 1165 PULL_EXPORT int pullTraceSet(pullContext *pctx, pullTrace *trc, 1166 int recordStrength, 1167 double scaleDelta, double halfScaleWin, 1168 double velocityMax, unsigned int arrIncr, 1169 const double seedPos[4]); 1170 PULL_EXPORT pullTraceMulti *pullTraceMultiNew(void); 1171 PULL_EXPORT pullTraceMulti *pullTraceMultiNix(pullTraceMulti *mtrc); 1172 PULL_EXPORT int pullTraceMultiAdd(pullTraceMulti *mtrc, pullTrace *trc, 1173 int *addedP); 1174 PULL_EXPORT int pullTraceMultiFilterConcaveDown(Nrrd *nfilt, 1175 const pullTraceMulti *mtrc, 1176 double winLenFrac); 1177 PULL_EXPORT int pullTraceMultiPlotAdd(Nrrd *nplot, 1178 const pullTraceMulti *mtrc, 1179 const Nrrd *nfilt, 1180 unsigned int trcIdxMin, 1181 unsigned int trcNum); 1182 PULL_EXPORT int pullTraceMultiWrite(FILE *file, const pullTraceMulti *mtrc); 1183 PULL_EXPORT int pullTraceMultiRead(pullTraceMulti *mtrc, FILE *file); 1184 1185 /* actionPull.c */ 1186 PULL_EXPORT int pullEnergyPlot(pullContext *pctx, Nrrd *nplot, 1187 double xx, double yy, double zz, 1188 unsigned int res); 1189 PULL_EXPORT int pullBinProcess(pullTask *task, unsigned int myBinIdx); 1190 PULL_EXPORT int pullGammaLearn(pullContext *pctx); 1191 1192 /* corePull.c */ 1193 PULL_EXPORT int pullStart(pullContext *pctx); 1194 PULL_EXPORT int pullRun(pullContext *pctx); 1195 PULL_EXPORT int pullFinish(pullContext *pctx); 1196 1197 /* ccPull.c */ 1198 PULL_EXPORT int pullCCFind(pullContext *pctx); 1199 PULL_EXPORT int pullCCMeasure(pullContext *pctx, Nrrd *nmeas, 1200 int measrInfo, double rho); 1201 PULL_EXPORT int pullCCSort(pullContext *pctx, int measrInfo, double rho); 1202 1203 #ifdef __cplusplus 1204 } 1205 #endif 1206 1207 #endif /* PULL_HAS_BEEN_INCLUDED */ 1208