Merge branch 'dev' of github.com:freesurfer/freesurfer into dev

Author: Douglas Greve

distribution/FreeSurferColorLUT.txt

@@ -277,6 +277,15 @@
 265 Left-Eyeball                          60  60  60  0
 266 Right-Eyeball                         60  60  60  0
 
+270 ctl-lh-infragranular                  205 70  78  0
+271 ctl-lh-layer1                         210 80  78  0
+272 ctl-lh-layer2                         215 90  78  0
+273 ctl-lh-layer3                         220 100  78  0
+274 ctl-lh-layer4                         225 110  78  0
+275 ctl-lh-layer5                         230 120  78  0
+276 ctl-lh-layer6                         235 135  78  0
+277 ctx-lh-supragranular                  240 140  78  0
+
 # lymph node and vascular labels
 331 Aorta                                   255 0   0   0
 332 Left-Common-IliacA                      255 80  0   0
@@ -308,6 +317,15 @@
 358 Pos-Lymph                               20  130 180 0
 359 Neg-Lymph                               20  180 130 0
 
+370 ctl-rh-infragranular                  205 70  130  0
+371 ctl-rh-layer1                         210 80  130  0
+372 ctl-rh-layer2                         215 90  130  0
+373 ctl-rh-layer3                         220 100  130  0
+374 ctl-rh-layer4                         225 110  130  0
+375 ctl-rh-layer5                         230 120  130  0
+376 ctl-rh-layer6                         235 135  130  0
+377 ctx-rh-supragranular                  240 140  130  0
+
 400 V1                                      206 62  78  0
 401 V2                                      121 18  134 0
 402 BA44                                    199 58  250 0

mri_compute_volume_fractions/mri_compute_layer_fractions.cpp

@@ -168,7 +168,7 @@ main(int argc, char *argv[])
 
   MRIreInitCache(mri_layers) ; 
 
-  if (FS_names && nlayers != 1)
+  if (FS_names && nlayers > 2)
     ErrorExit(ERROR_UNSUPPORTED, "%s: if specifying FS_names must use -nlayers 1", Progname) ;
   printf("reading laminar surfaces from %s.?\n", LAMINAR_NAME) ;
   for (i = 0 ; i <= nlayers ; i++) {
@@ -186,8 +186,15 @@ main(int argc, char *argv[])
 	}
       }
     } else {
-      int req = snprintf(fname, STRLEN,
-			 "%s/%s/surf/%s.%s.%d", sdir, subject,hemi,LAMINAR_NAME,i) ;
+      int req ; 
+      if (FS_names && i == 0)
+	req = snprintf(fname, STRLEN, "%s/%s/surf/%s.white", sdir, subject,hemi) ;
+      else if (FS_names && i == 2)
+	req = snprintf(fname, STRLEN, "%s/%s/surf/%s.pial", sdir, subject,hemi) ;
+      else if (FS_names)
+	req = snprintf(fname, STRLEN, "%s/%s/surf/%s.%s", sdir, subject,hemi,LAMINAR_NAME) ;
+      else
+	req = snprintf(fname, STRLEN, "%s/%s/surf/%s.%s.%d", sdir, subject,hemi,LAMINAR_NAME,i) ;
       if( req >= STRLEN ) {
         std::cerr << __FUNCTION__ << ": Truncation on line " << __LINE__ << std::endl;
       }

mris_convert/mris_convert.cpp

@@ -465,13 +465,12 @@ main(int argc, char *argv[])
       mris->vg.valid = 0;
     }
 
-    // ??? should we call MRIcopyMRIS() instead ???
     if(MRISfileNameType(out_fname) == MRIS_VOLUME_FILE) {
       printf("Saving surface xyz %s as a volume format\n",out_fname);
       MRI *vol = MRIallocSequence(mris->nvertices, 1, 1, MRI_FLOAT, 3);
-      MRIScopyMRI(mris,vol,0,"x");
-      MRIScopyMRI(mris,vol,1,"y");
-      MRIScopyMRI(mris,vol,2,"z");
+      MRIcopyMRIS(vol,mris,0,"x");
+      MRIcopyMRIS(vol,mris,1,"y");
+      MRIcopyMRIS(vol,mris,2,"z");
       MRIwrite(vol,out_fname);
       MRIfree(&vol);
     }

mris_expand/mris_expand.cpp

@@ -36,6 +36,7 @@
 #include "const.h"
 #include "timer.h"
 #include "version.h"
+#include "randomfields.h"
 
 int main(int argc, char *argv[]) ;
 static int get_option(int argc, char *argv[]) ;
@@ -47,7 +48,12 @@ static int use_thickness = 0 ;
 static int nsurfaces = 1 ;
 static const char *thickness_name = "thickness" ;
 static const char *pial_name = "pial" ;
-static const char *tmap_fname = NULL ;
+static char *tmap_fname = NULL ;
+static char *tmap_write_fname = NULL ;
+static float tmap_std = 0.0 ;
+static float tmap_min = .25 ;
+static float tmap_max = .75 ;
+static int tmap_avgs = 0 ;
 static int nbrs = 2 ;
 
 static char *orig_name = NULL ;
@@ -173,12 +179,50 @@ main(int argc, char *argv[])
     parms.mri_dtrans = MRIread(tmap_fname);
     if (parms.mri_dtrans == NULL)
       ErrorExit(ERROR_NOFILE,
-                "%s: could not tmap vertex percentages from %s\n", Progname, tmap_fname) ;
+                "%s: could not read tmap vertex percentages from %s\n", Progname, tmap_fname) ;
     if (parms.mri_dtrans->width != mris->nvertices)
       ErrorExit(ERROR_NOFILE,
-                "%s: could not tmap width %d != mris->nvertices %d in %s\n", 
+                "%s: tmap width %d != mris->nvertices %d in %s\n", 
 		Progname, parms.mri_dtrans->width, mris->nvertices, tmap_fname) ;
-    printf("setting vno distance to %f\n", MRIgetVoxVal(parms.mri_dtrans, 0, 0, 0, 0));
+  }
+  else if (!FZERO(tmap_std))    // create a random map of distances
+  {
+    RFS    *rfs = RFspecInit(0, NULL);
+    int    vno ;
+
+    parms.mri_dtrans = MRIalloc(mris->nvertices, 1, 1, MRI_FLOAT) ;
+    if (parms.mri_dtrans == NULL)
+      ErrorExit(ERROR_NOFILE,
+                "%s: could not allocate tmap vertex percentages from %s\n", Progname, tmap_fname) ;
+
+    printf("creating random tmap distances\n") ;
+    rfs->name = strcpyalloc("gaussian");
+    rfs->params[0] = mm_out ;         // mean is halfway through the ribbon
+    rfs->params[1] = tmap_std;     // std
+
+    for (vno = 0 ; vno < mris->nvertices ; vno++)
+      mris->vertices[vno].val = RFdrawVal(rfs);
+
+    MRISaverageVals(mris, tmap_avgs) ;
+    for (vno = 0 ; vno < mris->nvertices ; vno++)
+    {
+      float dist = mris->vertices[vno].val ;
+
+      dist = MIN(MAX(dist, tmap_min), tmap_max) ;
+      if (mris->vertices[vno].ripflag)
+	dist = 0 ;
+
+      if (vno == Gdiag_no)
+	printf("vno %d: val %f, ripflag %d\n", vno, dist, mris->vertices[vno].ripflag) ;
+
+      mris->vertices[vno].val = dist ;
+      MRIsetVoxVal(parms.mri_dtrans, vno, 0, 0, 0, dist) ;
+    }
+    if (tmap_write_fname != NULL)
+    {
+      printf("writing random tmap to %s\n", tmap_write_fname) ;
+      MRISwriteValues(mris, tmap_write_fname) ;
+    }
   }
   else
     parms.mri_dtrans = NULL ;
@@ -329,6 +373,12 @@ get_option(int argc, char *argv[])
     nargs = 1 ;
     fprintf(stderr, "l_surf_repulse = %2.3f\n", parms.l_surf_repulse) ;
   }
+  else if (!stricmp(option, "wd"))
+  {
+    tmap_write_fname = argv[2] ;
+    nargs = 1 ;
+    fprintf(stderr, "writing random tmap to %s\n", tmap_write_fname) ;
+  }
   else if (!stricmp(option, "pial"))
   {
     pial_name = argv[2] ;
@@ -338,9 +388,22 @@ get_option(int argc, char *argv[])
   else if (!stricmp(option, "tmap"))
   {
     use_thickness = -1 ;
-    tmap_fname = argv[2] ;
-    printf("reading thickness target percent map from %s\n", tmap_fname) ;
-    nargs = 1 ;
+    if (!stricmp(argv[2], "random"))
+    {
+      tmap_std = atof(argv[3]);
+      tmap_min = atof(argv[4]);
+      tmap_max = atof(argv[5]);
+      tmap_avgs = atoi(argv[6]);
+      printf("creating random tmap in [%2.2f, %2.2f] with std %2.2f and %d averages\n",
+	     tmap_min, tmap_max, tmap_std, tmap_avgs) ;
+      nargs = 5 ;
+    }
+    else
+    {
+      tmap_fname = argv[2] ;
+      printf("reading thickness target percent map from %s\n", tmap_fname) ;
+      nargs = 1 ;
+    }
   }
   else switch (toupper(*option))
     {
@@ -414,5 +477,10 @@ usage_exit(int code)
   printf("  Example: mris_expand -thickness lh.white 0.5 lh.graymid\n");
   printf("  Example: mris_expand -label labelfile lh.white 0.5 lh.graymid\n");
   printf("  Example: mris_expand -tmap thickness_pct_target.mgz lh.white 0.5 lh.graymid\n");
+  printf("     use a prespecified map of percent thickness to compute the target locations for expansion\n");
+  printf("  Example: mris_expand -tmap random 2 .25 .75 100 -wd tmap.mgz lh.white 0.5 lh.graymid\n");
+  printf("     creates a random target distance map with gaussian sampling (mean=.5, std=2) and cropping to .25/.75\n");
+  printf("     and spatial averaging 100 times (cropping is average averaging). The map will be written to \n");
+  printf("     file tmap.mgz\n");
   exit(code) ;
 }

utils/gcsa.cpp

@@ -2248,9 +2248,7 @@ int GCSArelabelWithAseg(GCSA *gcsa, MRI_SURFACE *mris, MRI *mri_aseg)
 {
   int old_index, vno, vno_classifier, vno_prior, label, index, changed, cc_annotation;
   VERTEX *v, *v_classifier, *v_prior;
-  GCSA_NODE *gcsan;
-  CP_NODE *cpn;
-  double v_inputs[100], p;
+  double v_inputs[100], p[3];
   double x, y, z;
 
   for (changed = vno = 0; vno < mris->nvertices; vno++) {
@@ -2267,11 +2265,11 @@ int GCSArelabelWithAseg(GCSA *gcsa, MRI_SURFACE *mris, MRI *mri_aseg)
     v_classifier = GCSAsourceToClassifierVertex(gcsa, v_prior);
     vno_classifier = v_classifier - gcsa->mris_classifiers->vertices;
     if (vno_classifier == Gdiag_no) DiagBreak();
-    gcsan = &gcsa->gc_nodes[vno_classifier];
+    GCSA_NODE *gcsan = &gcsa->gc_nodes[vno_classifier];
 
     CTABfindAnnotation(mris->ct, v->annotation, &old_index);
 
-    cpn = &gcsa->cp_nodes[vno_prior];
+    CP_NODE *cpn = &gcsa->cp_nodes[vno_prior];
     if (IS_CC(label)) {
       CTABfindName(mris->ct, "corpuscallosum", &index);
       if (index < 0) CTABfindName(mris->ct, "Medial_wall", &index);
@@ -2291,7 +2289,7 @@ int GCSArelabelWithAseg(GCSA *gcsa, MRI_SURFACE *mris, MRI *mri_aseg)
     else if (old_index >= 0 && mris->ct && !stricmp(mris->ct->entries[old_index]->name, "corpuscallosum")) {
       // find 2nd most likely label that isn't callosum
       CTABannotationAtIndex(mris->ct, old_index, &cc_annotation);
-      label = GCSANclassify(gcsan, cpn, v_inputs, gcsa->ninputs, &p, &cc_annotation, 1, vno);
+      label = GCSANclassify(gcsan, cpn, v_inputs, gcsa->ninputs, p, &cc_annotation, 1, vno);
       if (label != v->annotation) {
         changed++;
         v->annotation = label;
@@ -2305,9 +2303,7 @@ int GCSAreclassifyMarked(GCSA *gcsa, MRI_SURFACE *mris, int mark, int *exclude_l
 {
   int old_index, vno, vno_classifier, vno_prior, label, index, changed, num, n;
   VERTEX *v, *v_classifier, *v_prior, *vn;
-  GCSA_NODE *gcsan;
-  CP_NODE *cpn;
-  double v_inputs[100], p;
+  double v_inputs[100], p[3];
 
   for (changed = vno = 0; vno < mris->nvertices; vno++) {
     v = &mris->vertices[vno];
@@ -2321,12 +2317,12 @@ int GCSAreclassifyMarked(GCSA *gcsa, MRI_SURFACE *mris, int mark, int *exclude_l
     v_classifier = GCSAsourceToClassifierVertex(gcsa, v_prior);
     vno_classifier = v_classifier - gcsa->mris_classifiers->vertices;
     if (vno_classifier == Gdiag_no) DiagBreak();
-    gcsan = &gcsa->gc_nodes[vno_classifier];
+    GCSA_NODE *gcsan = &gcsa->gc_nodes[vno_classifier];
 
     CTABfindAnnotation(mris->ct, v->annotation, &old_index);
 
-    cpn = &gcsa->cp_nodes[vno_prior];
-    label = GCSANclassify(gcsan, cpn, v_inputs, gcsa->ninputs, &p, exclude_list, nexcluded, vno);
+    CP_NODE *cpn = &gcsa->cp_nodes[vno_prior];
+    label = GCSANclassify(gcsan, cpn, v_inputs, gcsa->ninputs, p, exclude_list, nexcluded, vno);
     if (label >= 0 && label != v->annotation) {
       changed++;
       v->annotation = label;
@@ -2337,9 +2333,13 @@ int GCSAreclassifyMarked(GCSA *gcsa, MRI_SURFACE *mris, int mark, int *exclude_l
      that had nonzero priors and reclassify them based on nbrs. */
   do {
     for (num = vno = 0; vno < mris->nvertices; vno++) {
-      VERTEX_TOPOLOGY const * const vt = &mris->vertices_topology[vno];
       v = &mris->vertices[vno];
       if (v->ripflag || v->marked != mark) continue;
+      VERTEX_TOPOLOGY const * const vt = &mris->vertices_topology[vno];
+      /* moved to before assigning VERTEX_TOPOLOGY const * const vt
+       * v = &mris->vertices[vno];
+       * if (v->ripflag || v->marked != mark) continue;
+       */
       if (vno == Gdiag_no) DiagBreak();
       for (n = 0; n < vt->vnum; n++) {
         vn = &mris->vertices[vt->v[n]];