00186 {
00187     static v2s_param_t   P;
00188     v2s_opts_t         * sopt, sopt_def;
00189     v2s_results        * res;
00190 
00191 ENTRY("afni_vol2surf");
00192 
00193     if ( !gpar )                 RETURN(NULL);
00194 
00195     if (       check_SUMA_surface(sA) ) RETURN(NULL);
00196     if ( sB && check_SUMA_surface(sB) ) RETURN(NULL);
00197 
00198     if ( use_defaults )
00199     {
00200         sopt = &sopt_def;
00201         fill_sopt_default(sopt, sB ? 2 : 1);  /* 1 or 2 surfaces */
00202 
00203         /* but apply any debug options */
00204         sopt->debug = gv2s_plug_opts.sopt.debug;
00205         sopt->dnode = gv2s_plug_opts.sopt.dnode;
00206     }
00207     else 
00208         sopt = &gv2s_plug_opts.sopt;
00209 
00210     sopt->gp_index = gp_index;
00211 
00212     /* now fill the param struct based on the inputs */
00213     memset(&P, 0, sizeof(P));
00214     P.gpar         = gpar;
00215     P.cmask        = mask;
00216     P.nvox         = DSET_NVOX(gpar);
00217     P.over_steps   = v2s_vals_over_steps(sopt->map);
00218     P.nsurf        = sB ? 2 : 1;
00219     P.surf[0]      = *sA;
00220 
00221     /* verify steps, in case the user has not selected 2 surfaces */
00222     if ( P.nsurf == 1 && ! sopt->use_norms )
00223         sopt->f_steps = 1;
00224 
00225     if ( sB ) P.surf[1] = *sB;
00226 
00227     if ( gv2s_plug_opts.sopt.debug > 1 )
00228         disp_v2s_opts_t("   surf options: ", sopt);
00229 
00230     /* fire it up */
00231 
00232     res = vol2surf(sopt, &P);
00233 
00234     /* if the user wants output files, here they are (don't error check) */
00235     if (res && sopt->outfile_1D )
00236     {
00237         if ( THD_is_file(sopt->outfile_1D) )
00238             fprintf(stderr,"** over-writing 1D output file '%s'\n",
00239                     sopt->outfile_1D);
00240         v2s_write_outfile_1D(sopt, res, P.surf[0].label);
00241     }
00242 
00243     if (res && sopt->outfile_niml )
00244     {
00245         if ( THD_is_file(sopt->outfile_niml) )
00246             fprintf(stderr,"** over-writing niml output file '%s'\n",
00247                     sopt->outfile_niml);
00248         v2s_write_outfile_niml(sopt, res, 0);
00249     }
00250 
00251     RETURN(res);
00252 }

01493 {
01494 ENTRY("alloc_ints");
01495 
01496     *ptr = (int *)malloc(length * sizeof(int));
01497     if ( ! *ptr )
01498     {
01499         fprintf(stderr,"** ai: failed to alloc %d ints for '%s'\n",length,dstr);
01500         RETURN(1);
01501     }
01502     if ( debug > 1 )
01503         fprintf(stderr,"-d ai: alloc'd %d ints for '%s'\n", length, dstr);
01504 
01505     RETURN(0);
01506 }

01319 {
01320     v2s_results * sd;
01321     int           rv = 0, mem, nnodes;
01322 
01323 ENTRY("allocate_output_mem");
01324 
01325     /* if last <= 0, it will be set to nodes-1 */
01326     if ( sopt->first_node > sopt->last_node && sopt->last_node > 0 )
01327     {
01328         fprintf(stderr,"** error : first_node (%d) > last_node (%d)\n",
01329                 sopt->first_node, sopt->last_node);
01330         RETURN(NULL);
01331     }
01332 
01333     nnodes = p->surf[0].num_ixyz;               /* just for typing ease */
01334 
01335     sd = calloc(1, sizeof(v2s_results));
01336     if ( ! sd )
01337     {
01338         fprintf(stderr,"** aom: failed to allocate v2s_results struct\n");
01339         RETURN(NULL);
01340     }
01341 
01342     /* explicitly initialize all pointers with NULL */
01343     sd->nodes = sd->volind = sd->i = sd->j = sd->k = sd->nvals = NULL;
01344     sd->vals  = NULL;
01345  
01346     /* rcr - eventually, this may not apply */
01347     if ( sopt->first_node <  0      ) sopt->first_node = 0;
01348     if ( sopt->first_node >= nnodes ) sopt->first_node = nnodes-1;
01349     if ( sopt->last_node  <= 0      ) sopt->last_node  = nnodes-1;
01350     if ( sopt->last_node  >= nnodes ) sopt->last_node  = nnodes-1;
01351 
01352     sd->nalloc   = sopt->last_node - sopt->first_node + 1;
01353     sd->nused    = 0;
01354 
01355     /* decide the maximum number of entries per row */
01356     if ( p->over_steps )            sd->max_vals = sopt->f_steps;
01357     else if ( sopt->gp_index >= 0 ) sd->max_vals = 1;
01358     else                            sd->max_vals = DSET_NVALS(p->gpar);
01359 
01360     if ( sopt->skip_cols & V2S_SKIP_VALS ) sd->max_vals = 1;
01361 
01362     if ( p->over_steps && (sopt->skip_cols & V2S_SKIP_VALS) )
01363     {
01364         fprintf(stderr,"** if SKIP_VALS, cannot compute results over steps\n");
01365         free(sd);
01366         RETURN(NULL);
01367     }
01368 
01369     sd->memory   = 0;
01370 
01371     /* first, compute the memory needed for one row */
01372     mem = 0;
01373     if ( !(sopt->skip_cols & V2S_SKIP_NODES ) ) mem += sizeof(int);
01374     if ( !(sopt->skip_cols & V2S_SKIP_VOLIND) ) mem += sizeof(int);
01375     if ( !(sopt->skip_cols & V2S_SKIP_I     ) ) mem += sizeof(int);
01376     if ( !(sopt->skip_cols & V2S_SKIP_J     ) ) mem += sizeof(int);
01377     if ( !(sopt->skip_cols & V2S_SKIP_K     ) ) mem += sizeof(int);
01378     if ( !(sopt->skip_cols & V2S_SKIP_NVALS ) ) mem += sizeof(int);
01379 
01380     /* now note the actual output values */
01381     if ( !(sopt->skip_cols & V2S_SKIP_VALS  ) ) mem += sizeof(float);
01382     else                         mem += sd->max_vals * sizeof(float);
01383 
01384     mem *= sd->nalloc;  /* and multiply by the height of each column */
01385     sd->memory = mem;
01386 
01387     if ( sopt->debug > 1 )
01388         fprintf(stderr,"+d alloc result mem: %d bytes, max_vals = %d\n",
01389                 mem, sd->max_vals);
01390 
01391     /* okay, this time let's allocate something... */
01392 
01393     if ( ! (sopt->skip_cols & V2S_SKIP_NODES) )
01394         rv = alloc_ints(&sd->nodes, sd->nalloc, "nodes", sopt->debug);
01395 
01396     if ( ! rv && ! (sopt->skip_cols & V2S_SKIP_VOLIND) )
01397         rv = alloc_ints(&sd->volind, sd->nalloc, "volind", sopt->debug);
01398 
01399     if ( ! rv && ! (sopt->skip_cols & V2S_SKIP_I) )
01400         rv = alloc_ints(&sd->i, sd->nalloc, "i", sopt->debug);
01401 
01402     if ( ! rv && ! (sopt->skip_cols & V2S_SKIP_J) )
01403         rv = alloc_ints(&sd->j, sd->nalloc, "j", sopt->debug);
01404 
01405     if ( ! rv && ! (sopt->skip_cols & V2S_SKIP_K) )
01406         rv = alloc_ints(&sd->k, sd->nalloc, "k", sopt->debug);
01407 
01408     if ( ! rv && ! (sopt->skip_cols & V2S_SKIP_NVALS) )
01409         rv = alloc_ints(&sd->nvals, sd->nalloc, "nvals", sopt->debug);
01410 
01411     if ( ! rv && ! (sopt->skip_cols & V2S_SKIP_VALS) )
01412         rv = alloc_vals_list(&sd->vals, sd->nalloc, sd->max_vals, sopt->debug);
01413     else if ( ! rv )
01414         rv = alloc_vals_list(&sd->vals, sd->nalloc, 1, sopt->debug);
01415 
01416     if ( rv )
01417     {
01418         fprintf(stderr,"** failed v2s_results allocation\n");
01419         free_v2s_results(sd);
01420         sd = NULL;
01421     }
01422 
01423     RETURN(sd);
01424 }

01432 {
01433     int c;
01434 
01435 ENTRY("alloc_vals_list");
01436 
01437     *ptr = (float **)malloc(width * sizeof(float *));
01438     if ( !*ptr )
01439         fprintf(stderr,"** avl: failed to alloc %d floats pointers\n", width);
01440 
01441     for ( c = 0; c < width; c++ )
01442     {
01443         (*ptr)[c] = (float *)malloc(length * sizeof(float));
01444         if ( (*ptr)[c] == NULL )
01445             fprintf(stderr,"** avl: failed to alloc %d floats (# %d of %d)\n",
01446                     length, c, width);
01447     }
01448 
01449     if ( debug > 1 )
01450         fprintf(stderr,"-d alloc'd %d x %d floats for surf data\n",
01451                 width, length);
01452 
01453     RETURN(0);
01454 }

01960 {
01961     int rv = 0;
01962 ENTRY("is_valid_SUMA_surface");
01963 
01964     if ( !s ) { fprintf(stderr,"** ivSs: no surface\n");  RETURN(0); }
01965 
01966     if ( s->type != SUMA_SURFACE_TYPE )
01967     {
01968         fprintf(stderr,"** surf '%s': invalid type %d\n",
01969                 CHECK_EMPTY_STR(s->label), s->type);
01970         rv++;
01971     }
01972 
01973     if ( s->num_ixyz < 0 || s->nall_ixyz < 0 || s->num_ixyz < s->nall_ixyz )
01974     {
01975         fprintf(stderr,"** surf '%s': invalid num_ixyz = %d, nall_ixyz = %d\n",
01976                 CHECK_EMPTY_STR(s->label), s->num_ixyz, s->nall_ixyz);
01977         rv++;
01978     }
01979 
01980     if ( s->seq == 0 || s->sorted == 0 || s->seqbase != 0 )
01981     {
01982         fprintf(stderr,"** surf '%s': invalid seq, sort or base (%d, %d, %d)\n",
01983                 CHECK_EMPTY_STR(s->label), s->seq, s->sorted, s->seqbase);
01984         rv++;
01985     }
01986 
01987     if ( !s->ixyz )
01988     {
01989         fprintf(stderr,"** surf '%s': invalid, missing nodes\n",
01990                 CHECK_EMPTY_STR(s->label));
01991         rv++;
01992     }
01993 
01994     RETURN(rv);
01995 }

00309 {
00310     int rv = 0, mem = 0;
00311 ENTRY("compact_results");
00312 
00313     if ( sd->nused > sd->nalloc )  /* should not happen, of course */
00314     {
00315         fprintf(stderr,"** cr: nused (%d) > nalloc (%d) !!\n",
00316                 sd->nused, sd->nalloc);
00317         RETURN(-1);
00318     }
00319 
00320     if ( sd->nused == sd->nalloc ) RETURN(0);   /* we're good */
00321 
00322     /* otherwise, realloc everthing */
00323 
00324     sd->nalloc = sd->nused;
00325 
00326     if ( sd->nodes )
00327     {
00328         mem += sizeof(int);
00329         rv = realloc_ints(&sd->nodes, sd->nalloc, "nodes", debug);
00330     }
00331 
00332     if ( ! rv && sd->volind )
00333     {
00334         mem += sizeof(int);
00335         rv = realloc_ints(&sd->volind, sd->nalloc, "volind", debug);
00336     }
00337 
00338     if ( ! rv && sd->i )
00339     {
00340         mem += sizeof(int);
00341         rv = realloc_ints(&sd->i, sd->nalloc, "i", debug);
00342     }
00343 
00344     if ( ! rv && sd->j )
00345     {
00346         mem += sizeof(int);
00347         rv = realloc_ints(&sd->j, sd->nalloc, "j", debug);
00348     }
00349 
00350     if ( ! rv && sd->k )
00351     {
00352         mem += sizeof(int);
00353         rv = realloc_ints(&sd->k, sd->nalloc, "k", debug);
00354     }
00355 
00356     if ( ! rv && sd->nvals )
00357     {
00358         mem += sizeof(int);
00359         rv = realloc_ints(&sd->nvals, sd->nalloc, "nvals", debug);
00360     }
00361 
00362     if ( ! rv )
00363     {
00364         mem += (sizeof(float) * sd->max_vals);
00365         rv = realloc_vals_list(sd->vals, sd->nalloc, sd->max_vals, debug);
00366     }
00367 
00368     if ( rv ) RETURN(rv);       /* if there was a failure, just leave */
00369 
00370     mem *= sd->nalloc;          /* now multiply be the array length */
00371 
00372     if ( debug > 1 )
00373         fprintf(stderr,"+d compact results: reallocated %d bytes down to %d\n",
00374                 sd->memory, mem);
00375 
00376     sd->memory = mem;
00377 
00378     RETURN(rv);
00379 }

01195 {
01196     double mag, ratio;
01197     int    c;
01198 
01199 ENTRY("directed_dist");
01200 
01201     mag = magnitude_f(dir, 3);
01202 
01203     if ( mag < V2S_EPSILON )    /* can't be negative */
01204         ratio = 0.0;
01205     else
01206         ratio = dist/mag;
01207 
01208     for ( c = 0; c < 3; c++ )
01209         pnew[c] = pold[c] + dir[c]*ratio;
01210 
01211     RETURN(dist);
01212 }

01727 {
01728     float * fp;
01729     int     cr, cc;
01730 
01731 ENTRY("disp_mri_imarr");
01732 
01733     if ( info )
01734         fputs( info, stderr );
01735 
01736     if ( dp == NULL )
01737     {
01738         fprintf(stderr, "disp_mri_imarr: dp == NULL\n");
01739         RETURN(-1);
01740     }
01741 
01742     fprintf(stderr,
01743             "mri_imarr struct at %p :\n"
01744             "    num, nall = %d, %d\n",
01745             dp, dp->num, dp->nall );
01746 
01747     for ( cr = 0; cr < dp->num; cr++ )
01748     {
01749         fp = MRI_FLOAT_PTR(dp->imarr[cr]);
01750         fprintf(stderr, "    %3d: ", cr);
01751         for ( cc = 0; cc < dp->imarr[cr]->nx; cc++, fp++ )
01752             fprintf(stderr, "%f  ", *fp );
01753         fputc( '\n', stderr );
01754     }
01755 
01756     RETURN(0);
01757 }

01650 {
01651 ENTRY("disp_range_3dmm");
01652 
01653     if ( info )
01654         fputs( info, stderr );
01655 
01656     if ( dp == NULL )
01657     {
01658         fprintf(stderr, "disp_range_3dmm: dp == NULL\n");
01659         RETURN(-1);
01660     }
01661 
01662     fprintf(stderr,
01663             "range_3dmm struct at %p :\n"
01664             "    dset             = %p : %s\n"
01665             "    p1               = (%f, %f, %f)\n"
01666             "    pn               = (%f, %f, %f)\n"
01667             "    dset_min         = (%f, %f, %f)\n"
01668             "    dset_max         = (%f, %f, %f)\n"
01669             "    oob_check, debug = (%d, %d)\n",
01670             dp, dp->dset, ISVALID_DSET(dp->dset) ? "valid" : "invalid",
01671             dp->p1.xyz[0], dp->p1.xyz[1], dp->p1.xyz[2],
01672             dp->pn.xyz[0], dp->pn.xyz[1], dp->pn.xyz[2],
01673             dp->dset_min.xyz[0], dp->dset_min.xyz[1], dp->dset_min.xyz[2],
01674             dp->dset_max.xyz[0], dp->dset_max.xyz[1], dp->dset_max.xyz[2],
01675             dp->oob_check, dp->debug
01676         );
01677 
01678     RETURN(0);
01679 }

01687 {
01688 ENTRY("disp_range_3dmm_res");
01689 
01690     if ( info )
01691         fputs( info, stderr );
01692 
01693     if ( dp == NULL )
01694     {
01695         fprintf(stderr, "disp_range_3dmm: dp == NULL\n");
01696         RETURN(-1);
01697     }
01698 
01699     fprintf(stderr,
01700             "range_3dmm_res struct at %p :\n"
01701             "    ims.imarr         = %p\n"
01702             "    ims.num, ims.nall = %d, %d\n"
01703             "    repeats, masked   = %d, %d\n"
01704             "    oob, ifirst       = %d, %d\n"
01705             "    i3first[0,1,2]    = %d, %d, %d\n"
01706             "    i3arr             = %p\n"
01707             , dp,
01708             dp->ims.imarr, dp->ims.num, dp->ims.nall,
01709             dp->repeats, dp->masked, dp->oob, dp->ifirst,
01710             dp->i3first.ijk[0], dp->i3first.ijk[1], dp->i3first.ijk[2],
01711             dp->i3arr );
01712 
01713     if ( dp->i3arr )
01714         fprintf(stderr,
01715             "    i3arr[0].ijk       = %d, %d, %d\n",
01716             dp->i3arr[0].ijk[0], dp->i3arr[0].ijk[1], dp->i3arr[0].ijk[2] );
01717 
01718     RETURN(0);
01719 }

01765 {
01766     int index, c;
01767 
01768 ENTRY("disp_surf_vals");
01769 
01770     fprintf(stderr, "-------------------------------------------------\n");
01771     if ( mesg ) fputs( mesg, stderr );
01772     if ( sd->nused < 1 )
01773     {
01774         fprintf(stderr,"** no surf nodes defined\n");
01775         RETURN(-1);
01776     }
01777 
01778     index = (node >= 0) ? node : sd->nused - 1;
01779 
01780     fprintf(stderr, "v2s_results vals for sd_index %d, node %d :\n"
01781                     "    volind, (i, j, k) = %d, (%d, %d, %d)\n"
01782                     "    nvals: values...  = %d:  ", index,
01783                 sd->nodes  ? sd->nodes[index]  : 0,
01784                 sd->volind ? sd->volind[index] : 0,
01785                 sd->i      ? sd->i[index]      : 0,
01786                 sd->j      ? sd->j[index]      : 0,
01787                 sd->k      ? sd->k[index]      : 0,
01788                 sd->nvals  ? sd->nvals[index]  : 0 );
01789 
01790     for ( c = 0; c < sd->max_vals; c++ )
01791         fprintf(stderr,"%s  ", MV_format_fval(sd->vals[c][index]));
01792     fputc('\n', stderr);
01793 
01794     RETURN(0);
01795 }

01567 {
01568 ENTRY("disp_v2s_opts_t");
01569 
01570     if ( info )
01571         fputs( info, stderr );
01572 
01573     if ( sopt == NULL )
01574     {
01575         fprintf(stderr, "disp_v2s_opts_t: sopt == NULL\n");
01576         RETURN(-1);
01577     }
01578 
01579     fprintf(stderr,
01580             "v2s_opts_t struct at %p  :\n"
01581             "    map, gp_index         = %d, %d\n"
01582             "    debug, dnode          = %d, %d\n"
01583             "    no_head, skip_cols    = %d, %d\n"
01584             "    first_node, last_node = %d, %d\n"
01585             "    use_norms, norm_len   = %d, %f\n"
01586             "    norm_dir              = %d\n"
01587             "    f_index, f_steps      = %d, %d\n"
01588             "    f_p1_fr, f_pn_fr      = %f, %f\n"
01589             "    f_p1_mm, f_pn_mm      = %f, %f\n"
01590             "    outfile_1D            = %s\n"
01591             "    outfile_niml          = %s\n"
01592             , sopt,
01593             sopt->map, sopt->gp_index, sopt->debug, sopt->dnode,
01594             sopt->no_head, sopt->skip_cols,
01595             sopt->first_node, sopt->last_node,
01596             sopt->use_norms, sopt->norm_len, sopt->norm_dir,
01597             sopt->f_index, sopt->f_steps,
01598             sopt->f_p1_fr, sopt->f_pn_fr, sopt->f_p1_mm, sopt->f_pn_mm,
01599             CHECK_NULL_STR(sopt->outfile_1D),
01600             CHECK_NULL_STR(sopt->outfile_niml)
01601             );
01602 
01603     RETURN(0);
01604 }

01535 {
01536 ENTRY("disp_v2s_param_t");
01537 
01538     if ( info )
01539         fputs( info, stderr );
01540 
01541     if ( p == NULL )
01542     {
01543         fprintf(stderr, "disp_v2s_param_t: p == NULL\n" );
01544         RETURN(-1);
01545     }
01546 
01547     fprintf(stderr,
01548             "v2s_param_t struct at     %p :\n"
01549             "    gpar  : vcheck      = %p : %s\n"
01550             "    cmask               = %p\n"
01551             "    nvox, over_steps    = %d, %d\n"
01552             "    nsurf               = %d\n"
01553             , p,
01554             p->gpar, ISVALID_DSET(p->gpar) ? "valid" : "invalid",
01555             p->cmask, p->nvox, p->over_steps, p->nsurf
01556             );
01557 
01558     RETURN(0);
01559 }

01826 {
01827 ENTRY("disp_v2s_plugin_opts");
01828 
01829     if ( mesg ) fputs( mesg, stderr );
01830 
01831     fprintf(stderr, "v2s_plugin_opts @ %p\n"
01832                     "    ready      = %d\n"
01833                     "    use0, use1 = %d, %d\n"
01834                     "    s0A, s0B   = %d, %d\n"
01835                     "    s1A, s1B   = %d, %d\n",
01836                     d, d->ready, d->use0, d->use1,
01837                     d->s0A, d->s0B, d->s1A, d->s1B );
01838     RETURN(0);
01839 }

01803 {
01804 ENTRY("disp_v2s_results");
01805 
01806     if ( mesg ) fputs( mesg, stderr );
01807 
01808     fprintf(stderr, "v2s_results @ %p\n"
01809                     "    nalloc, nused    = %d, %d\n"
01810                     "    max_vals, memory = %d, %d\n"
01811                     "    nodes, volind    = %p, %p\n"
01812                     "    i, j, k          = %p, %p, %p\n"
01813                     "    nvals, vals      = %p, %p\n",
01814                     d, d->nalloc, d->nused, d->max_vals, d->memory,
01815                     d->nodes, d->volind, d->i, d->j, d->k, d->nvals, d->vals);
01816 
01817     RETURN(0);
01818 }

01628 {
01629     double d0, d1, d2;
01630 
01631     if ( p1 == NULL || p2 == NULL )
01632     {
01633         fprintf( stderr, "** dist_f3mm: invalid params (%p,%p)\n", p1, p2 );
01634         return(0.0);
01635     }
01636 
01637     d0 = p1->xyz[0] - p2->xyz[0];
01638     d1 = p1->xyz[1] - p2->xyz[1];
01639     d2 = p1->xyz[2] - p2->xyz[2];
01640 
01641     return(sqrt(d0*d0 + d1*d1 + d2*d2));
01642 }

00388 {
00389     range_3dmm_res r3mm_res;
00390     range_3dmm     r3mm;
00391     float          dist, min_dist, max_dist;
00392     int            sub, subs, nindex, findex = 0;
00393     int            oobc, oomc, max_index;
00394     int            oob1, oob2;
00395 
00396 ENTRY("dump_surf_3dt");
00397 
00398     if ( ! sopt || ! p || ! sd )
00399     {
00400         fprintf(stderr, "** ds3 : bad params (%p,%p,%p)\n", sopt, p, sd );
00401         RETURN(-1);
00402     }
00403 
00404     /* note the number of sub-bricks, unless the user has given just one */
00405     subs = sopt->gp_index >= 0 ? 1 : DSET_NVALS(p->gpar);
00406     init_range_structs( &r3mm, &r3mm_res );                 /* to empty */
00407     r3mm.dset = p->gpar;
00408     set_3dmm_bounds( p->gpar, &r3mm.dset_min, &r3mm.dset_max );
00409 
00410     if ( sopt->debug > 1 )
00411         fprintf(stderr, "-d dset bounding box: (%f, %f, %f)\n"
00412                         "                      (%f, %f, %f)\n",
00413                 r3mm.dset_min.xyz[0],r3mm.dset_min.xyz[1],r3mm.dset_min.xyz[2], 
00414                 r3mm.dset_max.xyz[0],r3mm.dset_max.xyz[1],r3mm.dset_max.xyz[2]);
00415 
00416     min_dist = 9999.9;                                          /* v2.3 */
00417     max_dist = -1.0;
00418     oobc     = 0;                         /* init out-of-bounds counter */
00419     oomc     = 0;                         /* init out-of-mask counter   */
00420 
00421     /* note, NodeList elements are in dicomm mm orientation */
00422 
00423     for ( nindex = sopt->first_node; nindex <= sopt->last_node; nindex++ )
00424     {
00425         /* init default max for oob and oom cases */
00426         max_index = p->over_steps ? sopt->f_steps : subs;
00427 
00428         init_seg_endpoints(sopt, p, &r3mm, nindex);    /* segment endpoints */
00429         v2s_adjust_endpts( sopt, &r3mm.p1, &r3mm.pn );
00430 
00431         if ( r3mm.debug )
00432             r3mm.debug = 0;
00433 
00434         if ( nindex == sopt->dnode )      /* if we have dnode, forget debug */
00435             r3mm.debug = sopt->debug > 0 ? sopt->debug : 1;
00436 
00437         /* if both points are outside our dataset, skip the pair   v2.3 */
00438         oob1 = f3mm_out_of_bounds( &r3mm.p1, &r3mm.dset_min, &r3mm.dset_max );
00439         oob2 = f3mm_out_of_bounds( &r3mm.pn, &r3mm.dset_min, &r3mm.dset_max );
00440         if ( oob1 && oob2 )
00441         {
00442             oobc++;
00443             if ( sopt->oob.show )
00444                 if ( set_all_surf_vals( sd, nindex, sopt->oob.index,
00445                                         sopt->oob.index, sopt->oob.index,
00446                                         sopt->oob.index, sopt->oob.value) )
00447                     RETURN(1);
00448             if ( nindex == sopt->dnode )
00449             {
00450                 disp_surf_vals("-d debug node, out-of-bounds : ", sd, -1);
00451                 fprintf(stderr,"-d dnode coords: (%f, %f, %f)\n",
00452                         r3mm.p1.xyz[0], r3mm.p1.xyz[1], r3mm.p1.xyz[2]);
00453             }
00454             continue;
00455         }
00456         else
00457             r3mm.oob_check = oob1 || oob2;
00458 
00459         dist = dist_f3mm( &r3mm.p1, &r3mm.pn );
00460         if ( dist < min_dist ) min_dist = dist;
00461         if ( dist > max_dist ) max_dist = dist;
00462 
00463         if ( segment_imarr( &r3mm_res, &r3mm, sopt, p->cmask ) != 0 )
00464             continue;
00465 
00466         if ( r3mm_res.ims.num == 0 )    /* any good voxels in the bunch? */
00467         {
00468             /* oob or oom? */
00469             if ( r3mm_res.oob == sopt->f_steps ) /* out of bounds */
00470             {
00471                 oobc++;
00472                 if ( sopt->oob.show )
00473                     if ( set_all_surf_vals( sd, nindex, sopt->oob.index,
00474                                             sopt->oob.index, sopt->oob.index,
00475                                             sopt->oob.index, sopt->oob.value) )
00476                         RETURN(1);
00477                 if ( nindex == sopt->dnode )
00478                     disp_surf_vals("-d debug node, out-of-bounds : ", sd, -1);
00479             }
00480             else   /* then we consider it out of mask */
00481             {
00482                 oomc++;
00483                 if ( sopt->oom.show )
00484                     if ( set_all_surf_vals( sd, nindex, r3mm_res.ifirst,
00485                             r3mm_res.i3first.ijk[0], r3mm_res.i3first.ijk[1],
00486                             r3mm_res.i3first.ijk[2], sopt->oom.value ) )
00487                         RETURN(1);
00488                 if ( nindex == sopt->dnode )
00489                     disp_surf_vals("-d debug node, out-of-mask : ", sd, -1);
00490             }
00491 
00492             if ( nindex == sopt->dnode )
00493                 fprintf(stderr,"-d dnode coords: (%f, %f, %f)\n",
00494                         r3mm.p1.xyz[0], r3mm.p1.xyz[1], r3mm.p1.xyz[2]);
00495 
00496             continue;   /* in any case */
00497         }
00498 
00499         /* get element 0, just for the findex */
00500         (void)v2s_apply_filter(&r3mm_res, sopt, 0, &findex);
00501 
00502         if ( set_surf_results(p, sopt, sd, &r3mm_res, nindex, findex) )
00503             RETURN(-1);
00504 
00505         /* clean up the MRI_IMARR struct, but don't free imarr */
00506         for ( sub = 0; sub < r3mm_res.ims.num; sub++ )
00507         {
00508             mri_free(r3mm_res.ims.imarr[sub]);
00509             r3mm_res.ims.imarr[sub] = NULL;
00510         }
00511         r3mm_res.ims.num = 0;
00512     }
00513 
00514     if ( sopt->debug > 0 )                                      /* v2.3 */
00515     {
00516         fprintf( stderr, "-d node pair dist (min,max) = (%f,%f)\n",
00517                  min_dist, max_dist );
00518         fprintf( stderr, "-d out-of-bounds, o-o-mask counts : %d, %d (of %d)\n",
00519                  oobc, oomc, sopt->last_node - sopt->first_node + 1);
00520     }
00521 
00522     /* now we can free the imarr and voxel lists */
00523     if ( r3mm_res.ims.nall > 0 )
00524     {
00525         free(r3mm_res.ims.imarr);
00526         free(r3mm_res.i3arr);
00527         r3mm_res.ims.nall = 0;
00528     }
00529 
00530     RETURN(0);
00531 }

00793 {
00794     int c;
00795 
00796     if ( !cp || !min || !max )
00797         return(-1);
00798 
00799     for ( c = 0; c < 3; c++ )
00800     {
00801         if ( ( cp->xyz[c] < min->xyz[c] ) ||
00802              ( cp->xyz[c] > max->xyz[c] ) )
00803             return(-1);
00804     }
00805 
00806     return(0);
00807 }

02198 {
02199 
02200 ENTRY("fill_sopt_default");
02201 
02202     if ( !sopt || nsurf < 1 || nsurf > 2 )
02203     {
02204         fprintf(stderr,"** fill_sopt_default: bad params (%p,%d)\n",sopt,nsurf);
02205         RETURN(1);
02206     }
02207 
02208     /* first set any zeros */
02209     memset(sopt, 0, sizeof(*sopt));
02210 
02211     if ( nsurf == 2 )
02212         sopt->map = E_SMAP_MIDPT;
02213     else
02214         sopt->map = E_SMAP_MASK;
02215 
02216     sopt->gp_index     = -1;
02217     sopt->no_head      = 1;
02218     sopt->skip_cols    = V2S_SKIP_ALL ^ V2S_SKIP_NODES;  /* nodes and 1 val */
02219     sopt->f_steps      = 1;
02220     sopt->outfile_1D   = NULL;
02221     sopt->outfile_niml = NULL;
02222 
02223     RETURN(0);
02224 }

01157 {
01158 ENTRY("float_list_alloc");
01159 
01160     if ( (f->nalloc < size) ||
01161          (trunc && (f->nalloc > size)) )
01162     {
01163         f->list = (float *)realloc(f->list, size * sizeof(float));
01164         if ( ! f->list )
01165         {
01166             fprintf(stderr,"** float_list_alloc: failed for %d floats\n", size);
01167             RETURN(-2);
01168         }
01169         f->nalloc = size;
01170 
01171         if ( ilist )     /* then allocate accompanying ilist */
01172         {
01173             *ilist = (int *)realloc(*ilist, size * sizeof(int));
01174             if ( ! *ilist )
01175             {
01176                 fprintf(stderr,"** float_list_alloc: failed for %d ints\n",
01177                         size);
01178                 RETURN(-2);
01179             }
01180         }
01181 
01182         if ( trunc && (f->nused > f->nalloc) )
01183             f->nused = f->nalloc;
01184     }
01185 
01186     RETURN(0);
01187 }

01056 {
01057     float fcur;
01058     int   ncur, c;
01059 
01060 ENTRY("float_list_comp_mode");
01061 
01062     /* init default results */
01063     *nvals = ncur = 1;
01064     *mode  = fcur = f->list[0];
01065     *index = 0;
01066 
01067     for ( c = 1; c < f->nused; c++ )
01068     {
01069         if ( f->list[c] == fcur )
01070             ncur ++;
01071         else                        /* found a new entry to count   */
01072         {
01073             if ( ncur > *nvals )     /* keep track of any new winner */
01074             {
01075                 *mode  = fcur;
01076                 *nvals = ncur;
01077                 *index = c;
01078             }
01079 
01080             fcur = f->list[c];
01081             ncur = 1;
01082         }
01083     }
01084 
01085     if ( ncur > *nvals )     /* keep track of any new winner */
01086     {
01087         *mode  = fcur;
01088         *nvals = ncur;
01089         *index = c;
01090     }
01091 
01092     RETURN(0);
01093 }

01106 {
01107     float * list, save;
01108     int     c0, c1, sindex;
01109 
01110 ENTRY("float_list_slow_sort");
01111 
01112     list = f->list;  /* for any little speed gain */
01113 
01114     for ( c0 = 0; c0 < f->nused-1; c0++ )
01115     {
01116         sindex = c0;
01117         save   = list[c0];
01118 
01119         /* find smallest remaining */
01120         for ( c1 = c0+1; c1 < f->nused; c1++ )
01121             if ( list[c1] < save )
01122             {
01123                 sindex = c1;
01124                 save   = list[sindex];
01125             }
01126 
01127         /* swap if smaller was found */
01128         if ( sindex > c0 )
01129         {
01130             list[sindex] = list[c0];
01131             list[c0]     = save;
01132 
01133             if ( ilist )        /* make same swap of indices */
01134             {
01135                 c1            = ilist[sindex];
01136                 ilist[sindex] = ilist[c0];
01137                 ilist[c0]     = c1;
01138             }
01139         }
01140     }
01141 
01142     RETURN(0);
01143 }

01847 {
01848     int c;
01849                                                                                 
01850 ENTRY("free_v2s_results");
01851 
01852     if( ! sd ) RETURN(0);
01853                                                                                 
01854     if (sd->nodes)  { free(sd->nodes);   sd->nodes  = NULL; }
01855     if (sd->volind) { free(sd->volind);  sd->volind = NULL; }
01856     if (sd->i)      { free(sd->i);       sd->i      = NULL; }
01857     if (sd->j)      { free(sd->j);       sd->j      = NULL; }
01858     if (sd->k)      { free(sd->k);       sd->k      = NULL; }
01859     if (sd->nvals)  { free(sd->nvals);   sd->nvals  = NULL; }
01860                                                                                 
01861     if (sd->vals)
01862     {
01863         for ( c = 0; c < sd->max_vals; c++ )
01864             if ( sd->vals[c] ) { free(sd->vals[c]);  sd->vals[c] = NULL; }
01865                                                                                 
01866         free(sd->vals);
01867         sd->vals = NULL;
01868     }
01869 
01870     free(sd);
01871                                                                                 
01872     RETURN(0);
01873 }

01257 {
01258 ENTRY("init_range_structs");
01259 
01260     r3->dset        = NULL;
01261     r3->debug       = 0;
01262 
01263     res3->ims.num   = 0;
01264     res3->ims.nall  = 0;
01265     res3->ims.imarr = NULL;
01266     res3->i3arr     = NULL;
01267 
01268     RETURN(0);
01269 }

01221 {
01222     SUMA_surface * sp;
01223 ENTRY("init_seg_endpoints");
01224 
01225     /* get node from first surface */
01226     sp = p->surf;
01227     R->p1.xyz[0] = sp->ixyz[node].x;
01228     R->p1.xyz[1] = sp->ixyz[node].y;
01229     R->p1.xyz[2] = sp->ixyz[node].z;
01230 
01231     /* set pn based on other parameters */
01232     if ( sopt->use_norms )
01233         directed_dist(R->pn.xyz, R->p1.xyz, sp->norm[node].xyz, sopt->norm_len);
01234     else if ( p->nsurf > 1 )
01235     {
01236         /* get node from second surface */
01237         sp = p->surf + 1;
01238         R->pn.xyz[0] = sp->ixyz[node].x;
01239         R->pn.xyz[1] = sp->ixyz[node].y;
01240         R->pn.xyz[2] = sp->ixyz[node].z;
01241     }
01242     else
01243         R->pn = R->p1;
01244 
01245     R->p1 = THD_dicomm_to_3dmm(R->dset, R->p1);
01246     R->pn = THD_dicomm_to_3dmm(R->dset, R->pn);
01247 
01248     RETURN(0);
01249 }

01612 {
01613     int     c;
01614     double  sums = 0.0;
01615 
01616     for ( c = 0; c < length; c++, p++ )
01617         sums += (*p) * (*p);
01618 
01619     return(sqrt(sums));
01620 }

02328 {
02329     int val;
02330                                                                                 
02331 ENTRY("print_header");
02332                                                                                 
02333     fprintf( outfp, "# --------------------------------------------------\n" );
02334     fprintf( outfp, "# surface '%s', '%s' :\n", surf, map );
02335     fprintf( outfp, "#\n" );
02336                                                                                 
02337     /* keep old style, but don't presume all columns get used (v 6.0) :
02338      *     fprintf( outfp, "#    node     1dindex    i    j    k     vals" );
02339      *     fprintf( outfp, "#   ------    -------   ---  ---  ---    ----" );
02340      */
02341                                                                                 
02342     /* output column headers */
02343     fputc( '#', outfp );        /* still comment line */
02344                                                                                 
02345     if ( sd->nodes  ) fprintf(outfp, "    node ");
02346     if ( sd->volind ) fprintf(outfp, "    1dindex ");
02347     if ( sd->i      ) fprintf(outfp, "   i ");
02348     if ( sd->j      ) fprintf(outfp, "   j ");
02349     if ( sd->k      ) fprintf(outfp, "   k ");
02350     if ( sd->nvals  ) fprintf(outfp, "    vals");
02351                                                                                 
02352     for ( val = 0; val < sd->max_vals; val++ )
02353         fprintf( outfp, "       v%-2d  ", val );
02354     fputc( '\n', outfp );
02355                                                                                 
02356     fputc( '#', outfp );
02357     /* underline the column headers */
02358     if ( sd->nodes  ) fprintf(outfp, "   ------");
02359     if ( sd->volind ) fprintf(outfp, "    ------- ");
02360     if ( sd->i      ) fprintf(outfp, "  ---");
02361     if ( sd->j      ) fprintf(outfp, "  ---");
02362     if ( sd->k      ) fprintf(outfp, "  ---");
02363     if ( sd->nvals  ) fprintf(outfp, "    ----");
02364                                                                                 
02365     fputs( "   ", outfp );
02366     for ( val = 0; val < sd->max_vals; val++ )
02367         fprintf( outfp, " --------   " );
02368     fputc( '\n', outfp );
02369                                                                                 
02370     RETURN(0);
02371 }

01514 {
01515 ENTRY("realloc_ints");
01516 
01517     *ptr = (int *)realloc(*ptr, length * sizeof(int));
01518     if ( ! *ptr )
01519     {
01520         fprintf(stderr,"** ri: failed to alloc %d ints for '%s'\n",length,dstr);
01521         RETURN(1);
01522     }
01523     if ( debug > 1 )
01524         fprintf(stderr,"-d ri: alloc'd %d ints for '%s'\n", length, dstr);
01525 
01526     RETURN(0);
01527 }

01462 {
01463     int c, count;
01464 
01465 ENTRY("realloc_vals_list");
01466 
01467     count = 0;
01468     for ( c = 0; c < width; c++ )
01469     {
01470         if ( ptr[c] )
01471         {
01472             ptr[c] = (float *)realloc(ptr[c], length * sizeof(float));
01473             if ( ptr[c] == NULL )
01474                 fprintf(stderr,"** rvl: fail to realloc %d floats (%d of %d)\n",
01475                     length, c, width);
01476             count++;
01477         }
01478     }
01479 
01480     if ( debug > 1 )
01481         fprintf(stderr,"-d realloc'd %d x %d floats for surf data\n",
01482                 count, length);
01483 
01484     RETURN(0);
01485 }

00656 {
00657     THD_fvec3 f3mm;
00658     THD_ivec3 i3ind;
00659     float     rat1, ratn;
00660     int       nx, ny;
00661     int       step, vindex, prev_ind;
00662 
00663 ENTRY("segment_imarr");
00664 
00665     /* check params for validity */
00666     if ( !R || !sopt || !res || !R->dset )
00667     {
00668         fprintf(stderr, "** seg_imarr: invalid params (%p,%p,%p)\n",R,sopt,res);
00669         if ( R ) disp_range_3dmm("segment_imarr: bad inputs:", R );
00670         RETURN(-1);
00671     }
00672 
00673     if ( R->debug > 1 )
00674         disp_range_3dmm("-d segment_imarr: ", R );
00675 
00676     /* handle this as an acceptable, trivial case */
00677     if ( sopt->f_steps < 1 )
00678     {
00679         res->ims.num = 0;
00680         RETURN(0);
00681     }
00682 
00683     nx = DSET_NX(R->dset);
00684     ny = DSET_NY(R->dset);
00685 
00686     /* if we don't have enough memory for results, (re)allocate some */
00687     if ( res->ims.nall < sopt->f_steps )
00688     {
00689         if ( R->debug > 1 )
00690             fprintf(stderr,"+d realloc of imarr (from %d to %d pointers)\n",
00691                     res->ims.nall,sopt->f_steps);
00692 
00693         res->ims.nall  = sopt->f_steps;
00694         res->ims.imarr = realloc(res->ims.imarr,
00695                                  sopt->f_steps*sizeof(MRI_IMAGE *));
00696         res->i3arr     = realloc(res->i3arr, sopt->f_steps*sizeof(THD_ivec3));
00697         if ( !res->ims.imarr || !res->i3arr )
00698         {
00699             fprintf(stderr,"** seg_imarr: no memory for %d MRI_IMAGE ptrs\n",
00700                     sopt->f_steps);
00701             res->ims.nall = res->ims.num = 0;
00702             /* one might be good */
00703             if ( res->ims.imarr ) free(res->ims.imarr);
00704             if ( res->i3arr )     free(res->i3arr);
00705             RETURN(-1);
00706         }
00707     }
00708 
00709     /* init return structure */
00710     res->ims.num = 0;
00711     res->repeats = 0;
00712     res->masked  = 0;
00713     res->oob     = 0;
00714     res->i3first = THD_3dmm_to_3dind_no_wod( R->dset, R->p1 );
00715     res->ifirst  = res->i3first.ijk[0] +
00716                    nx * (res->i3first.ijk[1] + ny * res->i3first.ijk[2] );
00717 
00718     prev_ind = -1;                      /* in case we want unique voxels */
00719 
00720     for ( step = 0; step < sopt->f_steps; step++ )
00721     {
00722         /* set our endpoint ratios */
00723         if ( sopt->f_steps < 2 )     /* if this is based on a single point */
00724             ratn = 0.0;
00725         else
00726             ratn = (float)step / (sopt->f_steps - 1);
00727         rat1 = 1.0 - ratn;
00728 
00729         f3mm.xyz[0] = rat1 * R->p1.xyz[0] + ratn * R->pn.xyz[0];
00730         f3mm.xyz[1] = rat1 * R->p1.xyz[1] + ratn * R->pn.xyz[1];
00731         f3mm.xyz[2] = rat1 * R->p1.xyz[2] + ratn * R->pn.xyz[2];
00732 
00733         /* accept part being oob                30 Sep 2004 [rickr] */
00734         if ( R->oob_check && 
00735              f3mm_out_of_bounds( &f3mm, &R->dset_min, &R->dset_max ) )
00736         {
00737             res->oob++;
00738             continue;
00739         }
00740 
00741         i3ind = THD_3dmm_to_3dind_no_wod( R->dset, f3mm );
00742         vindex = i3ind.ijk[0] + nx * (i3ind.ijk[1] + ny * i3ind.ijk[2] );
00743 
00744         /* is this voxel masked out? */
00745         if ( cmask && !cmask[vindex] )
00746         {
00747             res->masked++;
00748             continue;
00749         }
00750 
00751         /* is this voxel repeated, and if so, do we skip it? */
00752         if ( sopt->f_index == V2S_INDEX_VOXEL && vindex == prev_ind )
00753         {
00754             res->repeats++;
00755             continue;
00756         }
00757 
00758         /* Huston, we have a good voxel... */
00759 
00760         prev_ind = vindex;              /* note this for next time  */
00761 
00762         /* now for the big finish, get and insert the actual data */
00763 
00764         res->i3arr    [res->ims.num] = i3ind;   /* store the 3-D indices */
00765         res->ims.imarr[res->ims.num] = THD_extract_series( vindex, R->dset, 0 );
00766         res->ims.num++;
00767 
00768         if ( ! res->ims.imarr[res->ims.num-1] )
00769         {
00770             fprintf(stderr,"** failed THD_extract_series, vox %d\n", vindex);
00771             RETURN(-1);
00772         }
00773 
00774         if ( R->debug > 2 )
00775             fprintf(stderr, "-d seg step %2d, vindex %d, coords %f %f %f\n",
00776                     step,vindex,f3mm.xyz[0],f3mm.xyz[1],f3mm.xyz[2]);
00777     }
00778 
00779     if ( R->debug > 1 )
00780         disp_range_3dmm_res( "+d i3mm_seg_imarr results: ", res );
00781 
00782     RETURN(0);
00783 }

01279 {
01280     float tmp;
01281     int   c;
01282                                                                                 
01283 ENTRY("set_3dmm_bounds");
01284                                                                                 
01285     if ( !dset || !min || !max )
01286     {
01287         fprintf(stderr, "** invalid params to set_3dmm_bounds: (%p,%p,%p)\n",
01288                 dset, min, max );
01289         RETURN(-1);
01290     }
01291                                                                                 
01292     /* get undirected bounds */
01293     min->xyz[0] = DSET_XORG(dset) - 0.5 * DSET_DX(dset);
01294     max->xyz[0] = min->xyz[0] + DSET_NX(dset) * DSET_DX(dset);
01295                                                                                 
01296     min->xyz[1] = DSET_YORG(dset) - 0.5 * DSET_DY(dset);
01297     max->xyz[1] = min->xyz[1] + DSET_NY(dset) * DSET_DY(dset);
01298                                                                                 
01299     min->xyz[2] = DSET_ZORG(dset) - 0.5 * DSET_DZ(dset);
01300     max->xyz[2] = min->xyz[2] + DSET_NZ(dset) * DSET_DZ(dset);
01301                                                                                 
01302     for ( c = 0; c < 3; c++ )
01303         if ( min->xyz[c] > max->xyz[c] )
01304         {
01305             tmp = min->xyz[c];
01306             min->xyz[c] = max->xyz[c];
01307             max->xyz[c] = tmp;
01308         }
01309                                                                                 
01310     RETURN(0);
01311 }

00615 {
00616     int           c, nused;
00617 
00618 ENTRY("set_all_surf_vals");
00619 
00620     nused = sd->nused;
00621 
00622     if ( nused >= sd->nalloc )
00623     {
00624         fprintf(stderr,"** sasv: nused=%d >= nalloc=%d!\n", nused, sd->nalloc);
00625         RETURN(1);
00626     }
00627 
00628     if ( sd->nodes  )  sd->nodes[nused]  = node_ind;
00629     if ( sd->volind )  sd->volind[nused] = vind;
00630     if ( sd->i      )  sd->i[nused]      = i;
00631     if ( sd->j      )  sd->j[nused]      = j;
00632     if ( sd->k      )  sd->k[nused]      = k;
00633     if ( sd->nvals  )  sd->nvals[nused]  = sd->max_vals;
00634 
00635     for ( c = 0; c < sd->max_vals; c++ )
00636         sd->vals[c][nused] = fval;
00637 
00638     sd->nused++;
00639 
00640     RETURN(0);
00641 }

00540 {
00541     int           i, j, k, volind;
00542     int           c, max_index;
00543 ENTRY("set_surf_results");
00544 
00545     if ( sd->nused >= sd->nalloc )
00546     {
00547         fprintf(stderr,"** ssr: nused (%d) >= nalloc (%d)!\n",
00548                 sd->nused, sd->nalloc);
00549         RETURN(1);
00550     }
00551 
00552     i = r3res->i3arr[findex].ijk[0];
00553     j = r3res->i3arr[findex].ijk[1];
00554     k = r3res->i3arr[findex].ijk[2];
00555 
00556     /* now get 3D and 1D coordinates */
00557     volind = i + DSET_NX(p->gpar) * (j + DSET_NY(p->gpar) * k );
00558 
00559     /* set everything but the values */
00560     if (sd->nodes )  sd->nodes[sd->nused]  = node;
00561     if (sd->volind)  sd->volind[sd->nused] = volind;
00562     if (sd->i     )  sd->i[sd->nused]      = i;
00563     if (sd->j     )  sd->j[sd->nused]      = j;
00564     if (sd->k     )  sd->k[sd->nused]      = k;
00565     if (sd->nvals )  sd->nvals[sd->nused]  = r3res->ims.num;
00566 
00567     /* set max_index, and adjust in case max_vals has been restricted */
00568     max_index = p->over_steps ? r3res->ims.num : DSET_NVALS(p->gpar);
00569     if ( max_index > sd->max_vals ) max_index = sd->max_vals;
00570 
00571     if ( sopt->gp_index >= 0 )
00572         sd->vals[0][sd->nused] =
00573             v2s_apply_filter(r3res, sopt, sopt->gp_index, NULL);
00574     else
00575         for ( c = 0; c < max_index; c++ )
00576             sd->vals[c][sd->nused] = v2s_apply_filter(r3res, sopt, c, NULL);
00577 
00578     /* possibly fill line with default if by steps and short */
00579     if ( max_index < sd->max_vals )
00580         for ( c = max_index; c < sd->max_vals; c++ )
00581             sd->vals[c][sd->nused] = 0.0;
00582 
00583     /* rcr : should I nuke the MRI images, and just copy what is needed? */
00584     if ( node == sopt->dnode )
00585     {
00586         fprintf(stderr,
00587                 "--------------------------------------------------\n");
00588         if ( ! p->over_steps && sopt->gp_index >= 0 )
00589             fprintf(stderr,"+d dnode %d gets %f from gp_index %d\n",
00590                     node, sd->vals[0][sd->nused], sopt->gp_index);
00591         if ( sopt->debug > 1 )
00592             fprintf(stderr, "-d debug: node %d, findex %d, vol_index %d\n",
00593                     node, findex, volind );
00594         if ( sopt->use_norms )
00595         {
00596             float * fp = p->surf[0].norm[node].xyz;
00597             fprintf(stderr,"-d normal %f, %f, %f\n", fp[0], fp[1], fp[2]);
00598         }
00599         disp_mri_imarr( "+d raw data: ", &r3res->ims );
00600     }
00601 
00602     sd->nused++;
00603 
00604     RETURN(0);
00605 }

02059 {
02060     float   factor;
02061     byte  * tmask;
02062     int     nvox, type, c, size = 0;
02063 
02064 ENTRY("thd_mask_from_brick");
02065 
02066     if ( mask ) *mask = NULL;   /* to be sure */
02067 
02068     if ( !ISVALID_DSET(dset) || ! mask || volume < 0 )
02069         RETURN(-1);
02070 
02071     if ( volume >= DSET_NVALS(dset) )
02072     {
02073         fprintf(stderr,"** tmfb: sub-brick %d out-of-range\n", volume);
02074         RETURN(-1);
02075     }
02076 
02077     nvox = DSET_NVOX(dset);
02078     type = DSET_BRICK_TYPE(dset, volume);
02079 
02080     if ( type != MRI_byte && type != MRI_short &&
02081          type != MRI_int && type != MRI_float )
02082     {
02083         fprintf(stderr,"** tmfb: invalid dataset type %s, sorry...\n",
02084                 MRI_type_name[type]);
02085         RETURN(-1);
02086     }
02087 
02088     tmask = (byte *)calloc(nvox, sizeof(byte));
02089     if ( ! tmask )
02090     {
02091         fprintf(stderr,"** tmfb: failed to allocate mask of %d bytes\n", nvox);
02092         RETURN(-1);
02093     }
02094 
02095     factor = DSET_BRICK_FACTOR(dset, volume);
02096 
02097     /* cheat: adjust threshold, not data */
02098     if ( factor != 0.0 ) thresh /= factor;
02099 
02100     switch( DSET_BRICK_TYPE(dset, volume) )
02101     {
02102         case MRI_byte:
02103         {
02104             byte * dp  = DSET_ARRAY(dset, volume);
02105             byte   thr = BYTEIZE(thresh + 0.99999);  /* ceiling */
02106             for ( c = 0; c < nvox; c++ )
02107                 if ( dp[c] != 0 && ( dp[c] >= thr ) )
02108                 {
02109                     size++;
02110                     tmask[c] = 1;
02111                 }
02112         }
02113             break;
02114 
02115         case MRI_short:
02116         {
02117             short * dp  = DSET_ARRAY(dset, volume);
02118             short   thr = SHORTIZE(thresh + 0.99999);  /* ceiling */
02119             for ( c = 0; c < nvox; c++, dp++ )
02120                 if ( *dp != 0 && ( *dp >= thr || (absolute && *dp <= -thr) ) )
02121                 {
02122                     size++;
02123                     tmask[c] = 1;
02124                 }
02125         }
02126             break;
02127 
02128         case MRI_int:
02129         {
02130             int * dp  = DSET_ARRAY(dset, volume);
02131             int   thr = (int)(thresh + 0.99999);  /* ceiling */
02132             for ( c = 0; c < nvox; c++, dp++ )
02133                 if ( *dp != 0 && ( *dp >= thr || (absolute && *dp <= -thr) ) )
02134                 {
02135                     size++;
02136                     tmask[c] = 1;
02137                 }
02138         }
02139             break;
02140 
02141         case MRI_float:
02142         {
02143             float * dp = DSET_ARRAY(dset, volume);
02144             for ( c = 0; c < nvox; c++, dp++ )
02145                 if (*dp != 0 && (*dp >= thresh || (absolute && *dp <= -thresh)))
02146                 {
02147                     size++;
02148                     tmask[c] = 1;
02149                 }
02150         }
02151             break;
02152 
02153         default:                /* let's be sure */
02154         {
02155             fprintf(stderr,"** tmfb: invalid dataset type, sorry...\n");
02156             free(tmask);
02157         }
02158             break;
02159     }
02160 
02161     *mask = tmask;
02162 
02163     RETURN(size);
02164 }

00818 {
00819     THD_fvec3 f3_diff;
00820     float     dist, factor;
00821 
00822 ENTRY("v2s_adjust_endpts");
00823 
00824     if ( !sopt || !p1 || !pn )
00825     {
00826         fprintf(stderr,"** v2s_ae: invalid params (%p,%p,%p)\n", sopt, p1, pn);
00827         RETURN(-1);
00828     }
00829 
00830     /* first, get the difference, and distance */
00831     f3_diff.xyz[0] = pn->xyz[0] - p1->xyz[0];
00832     f3_diff.xyz[1] = pn->xyz[1] - p1->xyz[1];
00833     f3_diff.xyz[2] = pn->xyz[2] - p1->xyz[2];
00834 
00835     dist = dist_f3mm( p1, pn );
00836 
00837     if ( (sopt->f_p1_fr != 0.0) || (sopt->f_p1_mm != 0.0) )
00838     {
00839         if ( sopt->f_p1_fr != 0.0 )     /* what the heck, choose fr if both */
00840             factor = sopt->f_p1_fr;
00841         else
00842             factor = (dist == 0.0) ? 0.0 : sopt->f_p1_mm / dist;
00843 
00844         p1->xyz[0] += factor * f3_diff.xyz[0];
00845         p1->xyz[1] += factor * f3_diff.xyz[1];
00846         p1->xyz[2] += factor * f3_diff.xyz[2];
00847     }
00848 
00849     if ( (sopt->f_pn_fr != 0.0) || (sopt->f_pn_mm != 0.0) )
00850     {
00851         if ( sopt->f_pn_fr != 0.0 )
00852             factor = sopt->f_pn_fr;
00853         else
00854             factor = (dist == 0.0) ? 0.0 : sopt->f_pn_mm / dist;
00855 
00856         pn->xyz[0] += factor * f3_diff.xyz[0];
00857         pn->xyz[1] += factor * f3_diff.xyz[1];
00858         pn->xyz[2] += factor * f3_diff.xyz[2];
00859     }
00860 
00861     switch ( sopt->map )
00862     {
00863         default:
00864             fprintf(stderr,"** v2s_ae: mapping %d not ready\n", sopt->map );
00865             RETURN(-1);
00866 
00867         case E_SMAP_AVE:
00868         case E_SMAP_MAX:
00869         case E_SMAP_MAX_ABS:
00870         case E_SMAP_MIN:
00871         case E_SMAP_MASK:
00872         case E_SMAP_SEG_VALS:
00873         case E_SMAP_MEDIAN:
00874         case E_SMAP_MODE:
00875             break;
00876 
00877         case E_SMAP_MIDPT:
00878 
00879             /* set the first point to be the average of the two */
00880             p1->xyz[0] = (p1->xyz[0] + pn->xyz[0]) / 2.0;
00881             p1->xyz[1] = (p1->xyz[1] + pn->xyz[1]) / 2.0;
00882             p1->xyz[2] = (p1->xyz[2] + pn->xyz[2]) / 2.0;
00883 
00884             break;
00885     }
00886 
00887     RETURN(0);
00888 }

00899 {
00900     static float_list_t   flist = { 0, 0, NULL };  /* for sorting results */
00901     static int          * ind_list = NULL;         /* track index sources */
00902     double                tmp, comp = 0.0;
00903     float                 fval;
00904     int                   count, source;
00905     int                   brick_index = 0;
00906 
00907 ENTRY("v2s_apply_filter");
00908 
00909     if ( !rr || !sopt || index < 0 )
00910     {
00911         fprintf(stderr,"** v2s_cm2: invalid params (%p,%p,%d)\n",
00912                 rr, sopt, index);
00913         RETURN(0.0);
00914     }
00915     
00916     if ( rr->ims.num <= 0 )
00917         RETURN(0.0);
00918 
00919     /* if sorting is required for resutls, do it now */
00920     if ( v2s_map_needs_sort( sopt->map ) )
00921     {
00922         if ( float_list_alloc( &flist, &ind_list, rr->ims.num, 0 ) != 0 )
00923             RETURN(0.0);
00924         for ( count = 0; count < rr->ims.num; count++ )
00925         {
00926             flist.list[count] = MRI_FLOAT_PTR(rr->ims.imarr[count])[index];
00927             ind_list  [count] = count;   /* init index sources */
00928         }
00929         flist.nused = rr->ims.num;
00930         float_list_slow_sort( &flist, ind_list );
00931     }
00932 
00933     switch ( sopt->map )
00934     {
00935         default:
00936             if ( findex ) *findex = 0;
00937             RETURN(0.0);
00938 
00939         case E_SMAP_AVE:
00940             if ( findex ) *findex = 0;
00941             for ( count = 0; count < rr->ims.num; count++ )
00942                 comp += MRI_FLOAT_PTR(rr->ims.imarr[count])[index];
00943 
00944             comp = comp / rr->ims.num;
00945             break;
00946 
00947         case E_SMAP_MASK:
00948         case E_SMAP_MIDPT:
00949             if ( findex ) *findex = 0;
00950             /* we have only the one point */
00951             comp = MRI_FLOAT_PTR(rr->ims.imarr[0])[index];
00952             break;
00953 
00954         case E_SMAP_MAX:
00955             comp = MRI_FLOAT_PTR(rr->ims.imarr[0])[index];
00956             if ( findex ) *findex = 0;
00957 
00958             for ( count = 1; count < rr->ims.num; count++ )
00959             {
00960                 tmp = MRI_FLOAT_PTR(rr->ims.imarr[count])[index];
00961                 if ( tmp > comp )
00962                 {
00963                     if ( findex ) *findex = count;
00964                     comp = tmp;
00965                 }
00966             }
00967             break;
00968 
00969         case E_SMAP_MAX_ABS:
00970             comp = MRI_FLOAT_PTR(rr->ims.imarr[0])[index];
00971             if ( findex ) *findex = 0;
00972 
00973             for ( count = 1; count < rr->ims.num; count++ )
00974             {
00975                 tmp = MRI_FLOAT_PTR(rr->ims.imarr[count])[index];
00976                 if ( fabs(tmp) > fabs(comp) )
00977                 {
00978                     if ( findex ) *findex = count;
00979                     comp = tmp;
00980                 }
00981             }
00982             break;
00983 
00984         case E_SMAP_MIN:
00985             comp = MRI_FLOAT_PTR(rr->ims.imarr[0])[index];
00986             if ( findex ) *findex = 0;
00987 
00988             for ( count = 1; count < rr->ims.num; count++ )
00989             {
00990                 tmp = MRI_FLOAT_PTR(rr->ims.imarr[count])[index];
00991                 if ( tmp < comp )
00992                 {
00993                     if ( findex ) *findex = count;
00994                     comp = tmp;
00995                 }
00996             }
00997             break;
00998 
00999         case E_SMAP_SEG_VALS:
01000             /* if the user has specified a brick index, use it */
01001             if ( sopt->gp_index >= 0 ) brick_index = sopt->gp_index;
01002             if ( findex ) *findex = 0;
01003             comp = MRI_FLOAT_PTR(rr->ims.imarr[index])[brick_index];
01004             break;
01005 
01006         case E_SMAP_MEDIAN:
01007             count = flist.nused >> 1;
01008             if ( (flist.nused & 1) || (count == 0) )
01009             {
01010                 comp = flist.list[count];
01011                 if ( findex ) *findex = ind_list[count];
01012             }
01013             else
01014             {
01015                 comp = (flist.list[count-1] + flist.list[count]) / 2;
01016                 if ( findex ) *findex = ind_list[count-1]; /* take first */
01017             }
01018             break;
01019 
01020         case E_SMAP_MODE:
01021             float_list_comp_mode(&flist, &fval, &count, &source);
01022             comp = fval;
01023             if ( findex ) *findex = ind_list[source];
01024             break;
01025     }
01026 
01027     RETURN((float)comp);
01028 }

02174 {
02175 ENTRY("v2s_good_map_index");
02176 
02177     if ( map <= E_SMAP_INVALID || map >= E_SMAP_FINAL )
02178         RETURN(0);
02179 
02180     /* these have not (yet? do we care?) been implemented */
02181     if ( map == E_SMAP_MASK2 || map == E_SMAP_COUNT )
02182         RETURN(0);
02183 
02184     if ( from_afni && map == E_SMAP_SEG_VALS )
02185         RETURN(0);
02186 
02187     RETURN(1);
02188 }

01039 {
01040     if ( (map == E_SMAP_MEDIAN) || (map == E_SMAP_MODE) )
01041         return 1;
01042 
01043     return 0;
01044 }

02022 {
02023     v2s_map_nums map;
02024                                                                                 
02025 ENTRY("v2s_map_type");
02026                                                                                 
02027     if ( map_str == NULL )
02028     {
02029         fprintf( stderr, "** v2s_map_type: missing map_str parameter\n" );
02030         RETURN((int)E_SMAP_INVALID);
02031     }
02032                                                                                 
02033     if ( sizeof(gv2s_map_names) / sizeof(char *) != (int)E_SMAP_FINAL )
02034     {
02035         fprintf( stderr, "** error:  gv2s_map_names/v2s_map_num mis-match\n");
02036         RETURN((int)E_SMAP_INVALID);
02037     }
02038                                                                                 
02039     /* not ready for E_SMAP_COUNT yet (until someone wants it) */
02040     if ( !strcmp( map_str, gv2s_map_names[E_SMAP_COUNT] ) )
02041         RETURN((int)E_SMAP_INVALID);
02042                                                                                 
02043     for ( map = E_SMAP_INVALID; map < E_SMAP_FINAL; map++ )
02044         if ( !strcmp( map_str, gv2s_map_names[map] ) )
02045             RETURN((int)map);
02046 
02047     RETURN((int)E_SMAP_INVALID);
02048 }

02006 {
02007     if ( map == E_SMAP_SEG_VALS )
02008         return 1;
02009 
02010     return 0;
02011 }

02232 {
02233     FILE        * fp;
02234     int           c, c2;
02235 ENTRY("v2s_write_outfile_1D");
02236 
02237     fp = fopen( sopt->outfile_1D, "w" );
02238     if ( fp == NULL )
02239     {
02240         fprintf( stderr, "** failure to open '%s' for writing\n",
02241                  sopt->outfile_1D );
02242         RETURN(-1);
02243     }
02244 
02245     if ( ! sopt->no_head )
02246         print_header(fp, label, gv2s_map_names[sopt->map], sd);
02247 
02248     for ( c = 0; c < sd->nused; c++ )
02249     {
02250         /* keep old spacing */
02251         fputc(' ', fp);
02252         if ( sd->nodes  ) fprintf(fp, " %8d", sd->nodes[c]);
02253         if ( sd->volind ) fprintf(fp, "   %8d ", sd->volind[c]);
02254         if ( sd->i      ) fprintf(fp, "  %3d", sd->i[c]);
02255         if ( sd->j      ) fprintf(fp, "  %3d", sd->j[c]);
02256         if ( sd->k      ) fprintf(fp, "  %3d", sd->k[c]);
02257         if ( sd->nvals  ) fprintf(fp, "     %3d", sd->nvals[c]);
02258 
02259         for ( c2 = 0; c2 < sd->max_vals; c2++ )
02260             fprintf(fp, "  %10s", MV_format_fval(sd->vals[c2][c]));
02261         fputc('\n', fp);
02262     }
02263 
02264     fclose(fp);
02265 
02266     RETURN(0);
02267 }

02275                                                                                 {
02276     static char   v2s_name[] = "3dVol2Surf_dataset";
02277     NI_element  * nel = NULL;
02278     NI_stream     ns;
02279     char        * ni_name;
02280     int           c;
02281 
02282 ENTRY("v2s_write_outfile_niml");
02283 
02284     if ( !sopt->outfile_niml ) RETURN(0);
02285 
02286     nel = NI_new_data_element( v2s_name, sd->nused );
02287     if ( !nel )
02288     {
02289         fprintf(stderr,"** file NI_new_data_element, n = '%s', len = %d\n",
02290                 v2s_name, sd->nused);
02291         RETURN(1);
02292     }
02293 
02294     ni_name = (char *)calloc(strlen(sopt->outfile_niml)+6, sizeof(char));
02295     if ( !ni_name ) { fprintf(stderr,"** ni_name failed\n"); RETURN(1); }
02296     sprintf(ni_name, "file:%s", sopt->outfile_niml);
02297 
02298     ns = NI_stream_open(ni_name, "w");
02299 
02300     NI_add_column(nel,NI_INT,sd->nodes);
02301 
02302     for ( c = 0; c < sd->max_vals; c++ )
02303     {
02304         NI_add_column(nel, NI_FLOAT, sd->vals[c]);
02305         if ( free_vals ) { free(sd->vals[c]); sd->vals[c] = NULL; }
02306     }
02307     if ( free_vals ) { free(sd->vals); sd->vals = NULL; }
02308 
02309     if ( NI_write_element(ns, nel, NI_BINARY_MODE) < 0 )
02310     {
02311         fprintf(stderr,"** NI_write_element failed for: '%s'\n", ni_name);
02312         RETURN(1);
02313     }
02314 
02315     NI_free_element( nel );
02316     NI_stream_close( ns );
02317     free(ni_name);
02318 
02319     RETURN(0);
02320 }

01881 {
01882     int c;
01883 
01884 ENTRY("validate_v2s_inputs");
01885 
01886     if ( !sopt || !p )
01887     {
01888         fprintf(stderr,"** validate_v2s_inputs: bad params (%p,%p)\n", sopt, p);
01889         RETURN(-1);
01890     }
01891 
01892     /* validate surface options structure */
01893     if ( sopt->map <= E_SMAP_INVALID || sopt->map >= E_SMAP_FINAL )
01894     {
01895         fprintf(stderr,"** invalid mapping index %d\n", sopt->map);
01896         RETURN(1);
01897     }
01898     else if ( v2s_map_type(gv2s_map_names[sopt->map]) != sopt->map )
01899     {
01900         fprintf(stderr,"** mapping index failure for %d\n", sopt->map);
01901         RETURN(1);
01902     }
01903 
01904     if ( sopt->f_steps <= 0 || sopt->f_steps >= V2S_STEPS_TOOOOO_BIG )
01905     {
01906         fprintf(stderr,"** invalid f_steps = %d\n", sopt->f_steps);
01907         RETURN(1);
01908     }
01909 
01910     /* validate the contents of p, proper */
01911     if ( !p->gpar ) {fprintf(stderr,"** vv2si: no dset?\n"); RETURN(2);}
01912 
01913     if ( p->nvox != DSET_NVOX(p->gpar) )
01914     {
01915         fprintf(stderr,"** invalid voxel count (%d) for grid_parent\n",p->nvox);
01916         RETURN(2);
01917     }
01918 
01919     if ( sopt->gp_index >= DSET_NVALS(p->gpar) )
01920     {
01921         fprintf(stderr,"** gp_index (%d) > max grid_parent index (%d)\n",
01922                 sopt->gp_index, DSET_NVALS(p->gpar) - 1); 
01923         RETURN(2);
01924     }
01925 
01926     if ( p->nsurf < 1 || p->nsurf > 2 )  /* see V2S_MAX_SURFS */
01927     {
01928         fprintf(stderr,"** invalid: nsurf = %d must be in [%d,%d]\n",
01929                 p->nsurf, 1, 2);
01930         RETURN(2);
01931     }
01932 
01933     /* validate individual SUMA_surface structs */
01934     for ( c = 0; c < p->nsurf; c++ )
01935         if ( check_SUMA_surface(p->surf + c) )
01936             RETURN(3);
01937 
01938     /* now look for consistency */
01939     if ( p->nsurf > 1 )
01940     {
01941         if ( p->surf[0].num_ixyz != p->surf[1].num_ixyz )
01942         {
01943             fprintf(stderr,"** invalid surfaces, different # nodes (%d,%d)\n",
01944                     p->surf[0].num_ixyz, p->surf[1].num_ixyz);
01945             RETURN(4);
01946         }
01947     }
01948     else if ( sopt->use_norms && !p->surf[0].norm )
01949     {
01950         fprintf(stderr,"** missing surface normals for surface '%s'\n",
01951                 CHECK_EMPTY_STR(p->surf[0].label));
01952         RETURN(4);
01953     }
01954 
01955     RETURN(0);
01956 }

00268 {
00269     v2s_results * sd;
00270     int           rv;
00271 ENTRY("vol2surf");
00272 
00273     if ( sopt == NULL || p == NULL )
00274     {
00275         fprintf( stderr, "** smd_wo - bad params (%p,%p)\n", sopt, p );
00276         RETURN(NULL);
00277     }
00278 
00279     if ( validate_v2s_inputs(sopt, p) )
00280         RETURN(NULL);
00281 
00282     if ( sopt->map == E_SMAP_INVALID )
00283     {
00284         fprintf(stderr,"** v2s wo: invalid map %d\n", sopt->map);
00285         RETURN(NULL);
00286     }
00287 
00288     sd = alloc_output_mem( sopt, p );
00289     if ( !sd ) RETURN(NULL);
00290 
00291     if ( sopt->debug > 1 ) disp_v2s_param_t( "-d post alloc_output_mem : ", p );
00292 
00293     rv = dump_surf_3dt( sopt, p, sd );
00294 
00295     if ( compact_results(sd, sopt->debug) )
00296     {
00297         free_v2s_results(sd);           /* free whatever didn't get burned */
00298         RETURN(NULL);
00299     }
00300 
00301     if ( sopt->debug > 1 ) disp_v2s_results( "-d post surf creation : ", sd);
00302 
00303     RETURN(sd);
00304 }

 { "none", "mask", "midpoint", "mask2", "ave",
                            "count", "min", "max", "max_abs", "seg_vals",
                            "median", "mode" }