/* just make sure we have sufficient data for computations */
int check_dims(options_t * opts)
{
   int nt, nvox, nmask;

   ENTRY("check_dims");

   nt = DSET_NVALS(opts->inset);
   nvox = DSET_NVOX(opts->inset);
   if( opts->mask ) nmask = THD_countmask( nvox, opts->mask );
   else             nmask = nvox;

   /* make sure we have something to compute */
   if( nvox < 1 ) {
      ERROR_message("input dataset must have at least 1 voxel");
      RETURN(1);
   } else if( nmask < 1 ) {
      ERROR_message("input mask must have at least 1 voxel");
      RETURN(1);
   } else if( nt < 2 ) {
      ERROR_message("input dataset must have at least 2 time points");
      RETURN(1);
   }

   RETURN(0);
}

char *phelp(char *prog, TFORM targ, int verb) 
{
   char cmd[512], tout[128];
   char *help=NULL;
   
   ENTRY("phelp");
   
   if (!prog ) RETURN(help);
   
   if (!phelp_cmd(prog, targ, cmd, tout, verb)) {
      ERROR_message("Failed to get help command");
      RETURN(0);
   }

   if (system(cmd)) {
      if (0) {/* many programs finish help and set status afterwards. Naughty. */
         ERROR_message("Failed to get help for %s\nCommand: %s\n", prog, cmd);
         return 0;
      }
   }
   
   if (!(help = AFNI_suck_file(tout))) {
      if (verb) ERROR_message("File %s could not be read\n", tout);
      RETURN(help);
   }
                                 
   snprintf(cmd,500*sizeof(char),"\\rm -f %s", tout);
   system(cmd);
   
   help = sphelp(prog, &help, targ, verb);
   
   RETURN(help);
}

char *form_C_progopt_string_from_struct(PROG_OPTS po)
{
   char *sout=NULL, sbuf[128];
   int maxch=0, i, jj, N_opts=0;
   
   if (!po.program) return(NULL);
   
   maxch = strlen(po.program)+strlen(po.options)+100;
   if (!(sout = (char *)calloc((maxch+1), sizeof(char)))) {
      ERROR_message("Failed to allocate for %d chars!", maxch+1);
      return(NULL);
   }

   sout[0]='\0';
   strncat(sout,"{ \"", maxch-1);
   strncat(sout,po.program, maxch-strlen(sout)-1);
   strncat(sout,"\", \"", maxch-strlen(sout)-1);
   strncat(sout,po.options, maxch-strlen(sout)-1);
   sprintf(sbuf,"\", %d", N_opts); strncat(sout,sbuf, maxch-strlen(sout)-1);

   strncat(sout,"}", maxch-strlen(sout)-1);
   if (strlen(sout)>=maxch-1) {
      ERROR_message("Truncated complete string possible");
      free(sout); sout=NULL;
      return(sout);
   }

   return(sout);
   
}

/*! Find afni's bin directory if it exists.
    If not, NULL is returned.  If it exists, a pointer to the path is returned.
    Do free it with free()
------------------------------------------------------------------------------*/
char * THD_abindir (byte withslash) 
{
   char *afr = NULL, *af=NULL;
   int  nn = 0, N_afni=strlen("afni");
   THD_string_array *elist=NULL;
   
   if (!(elist = get_elist()) ||
       !(af = THD_find_executable("afni"))) {
      ERROR_message("Could not find afni, we're doomed daddy!");
      RETURN(NULL);
   }
   
   /* remove afni from the end to get the path */
   nn = strlen(af);
   if (strcmp(af+nn-N_afni,"afni")) {
      ERROR_message("This should not be (%s)!", af+nn-N_afni);
      RETURN(NULL);
   }
   
   afr = strdup(af);
   afr[strlen(af)-N_afni]='\0'; 
   
   /* remove slash */
   while ( (nn=strlen(afr)-1) && afr[nn] == '/') 
      afr[nn] = '\0';

   if (withslash) {
      nn=strlen(afr);
      afr[nn] = '/'; afr[nn+1]='\0';
   }
   return(afr);
}

void view_prog_help(char *prog)
{
   char *viewer=NULL, *hname=NULL;
   char *progname=NULL;
   
   if (!prog) return;
   if (!(progname = THD_find_executable(prog))) {
      ERROR_message("Could not find executable %s.\n",
                     prog);
      return;
   }
   if (!(viewer = GetAfniTextEditor())) {
      ERROR_message("No GUI editor defined, and guessing game failed.\n"
              "Set AFNI_GUI_EDITOR in your .afnirc for this option to work.\n"); 
      return;
   }
   
   hname = get_updated_help_file(0, 0, progname, -1);
   if (hname[0]=='\0') { /* failed, no help file ... */
      ERROR_message("No help file for %s\n", progname);
      return;
   }
   
   if (!(view_text_file(hname))) {
      ERROR_message("Failed to view %s\n", hname);
   }
   return;
}

/* 
   Get copy contents of sub-brick iv into an double array.
   if iv == -1, get the entire dset 
*/
double *THD_extract_to_double( int iv , THD_3dim_dataset *dset )
{
   MRI_IMAGE *im ;
   double *var=NULL, *vv=NULL;
   register int ii , nvox ;

   ENTRY("THD_extract_to_double") ;

   if (!dset) RETURN(var);
   if (iv >= 0) {
      if (!(im = THD_extract_double_brick(iv, dset))) RETURN(var);
      var = MRI_DOUBLE_PTR(im);mri_fix_data_pointer(NULL, im); 
                              mri_free(im);im=NULL;
   } else if (iv == -1) {
      if (!(var = (double *)calloc(DSET_NVOX(dset)*DSET_NVALS(dset),
                                  sizeof(double)))){
         ERROR_message("Failed to allocate");
         RETURN(NULL);
      }
      for (ii=0; ii<DSET_NVALS(dset); ++ii) {
         if (!(im = THD_extract_double_brick(ii, dset))) {
            ERROR_message("Failed toextract sb %d from dset", ii);
            if (var) free(var);
            RETURN(NULL);
         }
         vv = MRI_DOUBLE_PTR(im);
         memcpy(var+ii*DSET_NVOX(dset),vv, sizeof(double)*DSET_NVOX(dset)); 
         mri_free(im);im=NULL;
      }
   } else {
      ERROR_message("Bad value of %d\n", iv);
   }
   
   RETURN(var);
}

/*--------------- check for invalid dilations ---------------*/
int dilations_are_valid(int_list * D)
{
   int * ilist;
   int   ival, err = 0;

   ENTRY("dilations_are_valid");

   ilist = D->list;
   for( ival = 0; ival < D->num; ival++ )
      if( ilist[ival] == 0 ) {
         ERROR_message("dilation[%d] is zero", ival);
         err++;
      }
   for( ival = 1; ival < D->num; ival++ )
      if( ilist[ival-1]*ilist[ival] > 0){
         ERROR_message(
            "have sequential dilations of same sign (assuming mistake)\n"
            "   (d[%d] = %d, d[%d] = %d)\n",
            ival-1,ilist[ival-1],ival,ilist[ival]);
         err++;
      }
   if( err ) RETURN(0);

   RETURN(1);
}

void THD_vectim_indexed_to_dset( MRI_vectim *mrv, int nlist, int *ilist,
                                 THD_3dim_dataset *dset )
{
   int nvals , nvec ,  jj,kk ;
   float *tar , *var ;

ENTRY("THD_vectim_indexed_to_dset") ;

   if( mrv   == NULL || !ISVALID_DSET(dset) ||
       nlist <= 0    || ilist == NULL       || nlist > DSET_NVALS(dset)  ){
     ERROR_message("THD_vectim_indexed_to_dset: illegal inputs (nlist=%d)",nlist) ;
     EXRETURN ;
   }

   nvec  = mrv->nvec ;
   nvals = mrv->nvals ;

   for( kk=0 ; kk < nlist ; kk++ ){
     if( ilist[kk] < 0 || ilist[kk] >= nvals ){
       ERROR_message("THD_vectim_indexed_to_dset: illegal ilist[%d]=%d",kk,ilist[kk]) ;
       EXRETURN ;
     }
   }

   tar = (float *)malloc(sizeof(float)*nlist) ;

   for( kk=0 ; kk < nvec ; kk++ ){
     var = VECTIM_PTR(mrv,kk) ;
     for( jj=0 ; jj < nlist ; jj++ ) tar[jj] = var[ilist[jj]] ;
     THD_insert_series( mrv->ivec[kk] , dset ,
                        nlist , MRI_float , tar , 0 ) ;
   }

   free(tar) ; EXRETURN ;
}

char *form_C_progopt_string(char *prog, char **ws, int N_ws) 
{
   char *sout=NULL, sbuf[128];
   int maxch=0, i, jj, N_opts=0;
   NI_str_array *nisa=NULL;
   
   if (!prog || !ws) {
      return(NULL);
   }
   
   maxch = 256;
   for (i=0; i<N_ws; ++i) {
      if (ws[i]) {
         maxch+=strlen(ws[i])+10;
         if (strlen(ws[i]) > 127) {
            WARNING_message("Truncating atrocious option %s\n", ws[i]);
            ws[127] = '\0';
         }
      }
   }
   if (!(sout = (char *)calloc((maxch+1), sizeof(char)))) {
      ERROR_message("Failed to allocate for %d chars!", maxch+1);
      return(NULL);
   }
   sout[0]='\0';
   strncat(sout,"{ \"", maxch-1);
   strncat(sout,prog, maxch-strlen(sout)-1);
   strncat(sout,"\", \"", maxch-strlen(sout)-1);
   
   N_opts = 0;
   for (i=0; i<N_ws; ++i) {
      if (ws[i] && (nisa = NI_strict_decode_string_list(ws[i] ,"/"))) {
         for (jj=0; jj<nisa->num; ++jj) {
            if (ws[i][0]=='-' && nisa->str[jj][0] != '-') {
               snprintf(sbuf,127,"-%s; ", nisa->str[jj]);
            } else { 
               snprintf(sbuf,127,"%s; ", nisa->str[jj]);
            }
            ++N_opts;
            strncat(sout,sbuf, maxch-strlen(sout)-1);
            NI_free(nisa->str[jj]);
         }
         if (nisa->str) NI_free(nisa->str); 
         NI_free(nisa); nisa=NULL;
      }
   }
   sprintf(sbuf,"\", %d", N_opts); strncat(sout,sbuf, maxch-strlen(sout)-1);
   
   
   strncat(sout,"}", maxch-strlen(sout)-1);
   if (strlen(sout)>=maxch-1) {
      ERROR_message("Truncated complete string possible");
      free(sout); sout=NULL;
      return(sout);
   }

   return(sout);
}

MRI_IMAGE * mri_genARMA11( int nlen, int nvec, float ap, float lm, float sg )
{
   int kk,ii , do_rcmat ;
   double aa=ap, lam=lm , sig=sg ; int do_norm = (sg<=0.0f) ;
   double *rvec ;
   rcmat *rcm=NULL ;
   MRI_IMAGE *outim ;
   float     *outar , *vv ;
#if 0
   long seed=0 ;
   seed = (long)time(NULL)+(long)getpid() ;
   srand48(seed) ;
#endif

ENTRY("mri_genARMA11") ;

   if( nlen      <= 3    ){ ERROR_message("ARMA11 nlen < 4");    RETURN(NULL); }
   if( fabs(aa)  >= 0.98 ){ ERROR_message("ARMA11 a too big");   RETURN(NULL); }
   if( fabs(lam) >= 0.97 ){ ERROR_message("ARMA11 lam too big"); RETURN(NULL); }

   do_rcmat = ( lam != 0.0 ) ;

   /* setup */

   if( do_rcmat ){
     rcm = rcmat_arma11( nlen , NULL , aa , lam ) ;
     if( rcm == NULL ){
       ERROR_message("Can't setup ARMA11 matrix?!"); RETURN(NULL);
     }
     kk = rcmat_choleski( rcm ) ;
     if( kk > 0 ){
       ERROR_message("ARMA11 Choleski fails at row %d",kk); RETURN(NULL);
     }
   }

   /* simulate */

   outim = mri_new( nlen , nvec , MRI_float ) ; outar = MRI_FLOAT_PTR(outim) ;
   rvec  = (double *)malloc(sizeof(double)*nlen) ;

   for( kk=0 ; kk < nvec ; kk++ ){
     for( ii=0 ; ii < nlen ; ii++ ) rvec[ii] = zgaussian() ;
     if( do_rcmat ) rcmat_lowert_vecmul( rcm , rvec ) ;
     vv = outar + kk*nlen ;
     if( do_norm ){
       sig = 0.0 ;
       for( ii=0 ; ii < nlen ; ii++ ) sig += rvec[ii]*rvec[ii] ;
       sig = 1.0 / sqrt(sig) ;
     }
     if( sig != 1.0 ){ for( ii=0 ; ii < nlen ; ii++ ) vv[ii] = sig * rvec[ii]; }
   }

   free(rvec) ;
   if( do_rcmat ) rcmat_destroy(rcm) ;

   RETURN(outim) ;
}

bytevec * THD_create_mask_from_string( char *str )  /* Jul 2010 */
{
   bytevec *bvec=NULL ; int nstr ; char *buf=NULL ;

ENTRY("THD_create_mask") ;

   if( str == NULL || *str == '\0' ) RETURN(NULL) ;

   nstr = strlen(str) ;
   bvec = (bytevec *)malloc(sizeof(bytevec)) ;

   /* try to read it as a dataset */

   if( nstr < THD_MAX_NAME ){
     THD_3dim_dataset *dset = THD_open_one_dataset(str) ;
     if( dset != NULL ){
       bvec->nar = DSET_NVOX(dset) ;
       bvec->ar  = THD_makemask( dset , 0 , 1.0f,0.0f ) ;
       DSET_delete(dset) ;
       if( bvec->ar == NULL ){
         ERROR_message("Can't make mask from dataset '%s'",str) ;
         free(bvec) ; bvec = NULL ;
       }
       RETURN(bvec) ;
     }
   }

   /* if str is a filename, read that file;
      otherwise, use the string itself to find the mask */

   if( THD_is_file(str) ){
     buf = AFNI_suck_file(str) ;
     if( buf != NULL ) nstr = strlen(buf) ;
   } else {
     buf = str ;
   }

   /* try to read buf as a Base64 mask string */

   if( strrchr(buf,'=') != NULL ){
     int nvox ;
     bvec->ar = mask_from_b64string( buf , &nvox ) ;
     if( bvec->ar != NULL ){
       bvec->nar = nvox ;
     } else {
       ERROR_message("Can't make mask from string '%.16s' %s",buf,(nstr<=16)?" ":"...") ;
       free(bvec) ; bvec = NULL ;
     }
   } else {
     ERROR_message("Don't understand mask string '%.16s'",buf,(nstr<=16)?" ":"...") ;
     free(bvec) ; bvec = NULL ;
   }

   if( buf != str && buf != NULL ) free(buf) ;
   RETURN(bvec) ;
}

int list_afni_files(int type, int withpath, int withnum)
{
   int nprogs=0, ii=0;
   char *etr=NULL, s[12];
   THD_string_array *progs=NULL;
   
   switch (type) {
      case 0:
         if (!(progs = THD_get_all_afni_executables())) {
            ERROR_message(
               "Cannot get list of programs from your afni bin directory %s", 
               THD_abindir(1));
            RETURN(0);
         }
         break;
      case 1:
         if (!(progs = THD_get_all_afni_readmes())) {
            ERROR_message(
               "Cannot get list of readmes from your afni bin directory %s", 
               THD_abindir(1));
            RETURN(0);
         }
         break;
      case 2:
         if (!(progs = THD_get_all_afni_dsets())) {
            ERROR_message(
               "Cannot get list of dsets from your afni bin directory %s", 
               THD_abindir(1));
            RETURN(0);
         }
         break;
      default:
         ERROR_message("Whatchyoutalkinboutwillis?");
         RETURN(0);
         break;
   }
   
   for (ii=0; ii<progs->num ; ii++ ){
      if (withpath) etr = progs->ar[ii];
      else etr = THD_trailname( progs->ar[ii] , 0 ) ;
      if (withnum) {
         sprintf(s,"%d", ii);
         fprintf(stdout,"  %3s.   %s\n", s, etr);
      } else {
         fprintf(stdout,"%s\n", etr);
      }
   }
   nprogs = progs->num;
   
   DESTROY_SARR(progs);
   
   return(nprogs);
}

SEXP R_SUMA_HistString (SEXP SCallingFunc, SEXP Sarg, SEXP Shold) {
   char *CallingFunc=NULL;
   SEXP Rname = R_NilValue;
   char *fname = NULL, *hold=NULL, **arg=NULL, *res=NULL;
   int debug=0, nprot=0, narg=0, i= 0;
   
   if (!debug) debug = get_odebug();
   if (isNull(SCallingFunc) || isNull(Sarg)) {
      ERROR_message("Null input to R_SUMA_HistString");
      return(Rname);
   }
   /* get the executable name */
   PROTECT(SCallingFunc = AS_CHARACTER(SCallingFunc)); ++nprot;
   fname = R_alloc(strlen(CHAR(STRING_ELT(SCallingFunc,0)))+1, sizeof(char));
   strcpy(fname, CHAR(STRING_ELT(SCallingFunc,0)));
   if (debug) INFO_message("filename %s\n", fname);
  
   /* get the arg */
   PROTECT(Sarg = AS_CHARACTER(Sarg)); ++nprot;
   narg = (LENGTH(Sarg));
   arg = (char **)calloc(narg+1, sizeof(char *));
   if (fname) arg[0] = strdup(fname);
   else arg[0] = strdup("UnChevalSansNom");
   for (i=1; i<=narg; ++i) {
      arg[i] = (char *)calloc(strlen(CHAR(STRING_ELT(Sarg,i-1)))+1, 
                              sizeof(char));
      strcpy(arg[i], CHAR(STRING_ELT(Sarg,i-1)));
      if (debug) INFO_message("arg %d/%d %s\t", i, narg, arg[i]);
   }
   
   /* any old history ? */
   if (!isNull(Shold)) {
      PROTECT(Shold = AS_CHARACTER(Shold)); ++nprot;
      hold = R_alloc(strlen(CHAR(STRING_ELT(Shold,0)))+1, sizeof(char));
      strcpy(hold, CHAR(STRING_ELT(Shold,0)));
      if (debug) INFO_message("hold %s\n", hold);
   }
    
   if (( res = SUMA_HistString (fname, narg+1, arg, hold))) {
      PROTECT(Rname = allocVector(STRSXP, 1)); ++nprot;
      SET_STRING_ELT(Rname, 0, mkChar(res)); 
      if (debug) INFO_message("hist is %s\n", res);
      SUMA_free(res); res=NULL;
   } else {
      ERROR_message("Call to SUMA_HistString %s failed", fname);
   }
   
   for (i=0; i<=narg; ++i) { if (arg[i]) free(arg[i]); } free(arg); arg=NULL;
   
   UNPROTECT(nprot); 
   return(Rname);
}

int GenFeatureDist(SEG_OPTS *Opt) 
{
   
   ENTRY("GenFeatureDist");
   
   
   /* get the probability maps */
   if (!Opt->pset && Opt->DO_p) {
      if (!(Opt->pset = p_C_GIV_A(Opt))) {
         ERROR_message("Failed miserably");
         RETURN(0);
      }
   }
      
   /* Get the classes */
   if (!Opt->cset && Opt->crefix && Opt->DO_c) {
      if (!(SUMA_assign_classes_eng(Opt->pset, 
                           Opt->clss->str, Opt->clss->num, Opt->keys, 
                           Opt->cmask, &Opt->cset))) {
         ERROR_message("Failed aimlessly");
         RETURN(0);
      }
      EDIT_dset_items(Opt->cset, ADN_prefix, Opt->crefix, ADN_none);
      if( !THD_ok_overwrite() && THD_is_file( DSET_HEADNAME(Opt->cset) ) ){
      ERROR_exit("Output file %s already exists -- cannot continue!\n",
                  DSET_HEADNAME(Opt->cset) ) ;
      }
   }  
   
   /* group classes ? */
   if (Opt->group_classes) {
      THD_3dim_dataset *gcset=NULL;
      THD_3dim_dataset *gpset=NULL;
      if (!SUMA_Regroup_classes (Opt, 
                     Opt->clss->str, Opt->clss->num, Opt->keys,
                     Opt->group_classes->str,
                     Opt->group_classes->num,
                     Opt->group_keys, Opt->cmask,
                     Opt->pset, 
                     Opt->cset,
                     &gpset, 
                     &gcset) ) {
         ERROR_message("Failed to regroup");
         RETURN(0);
      }
      DSET_write(gpset);
      DSET_write(gcset);
   }
          
   RETURN(1);
}

SEXP R_SUMA_ParseModifyName(SEXP Sfname, SEXP Swhat, SEXP Sval, SEXP Scwd) 
{
   SEXP Rname = R_NilValue;
   char *fname = NULL, *what=NULL, *val=NULL, *cwd=NULL, *res=NULL;
   int debug=0, nprot=0;
   
   if (!debug) debug = get_odebug();
   if (isNull(Sfname) || isNull(Swhat) || isNull(Sval)) {
      ERROR_message("Null input to R_SUMA_ModifyName");
      return(Rname);
   }
   /* get the filename */
   PROTECT(Sfname = AS_CHARACTER(Sfname)); ++nprot;
   fname = R_alloc(strlen(CHAR(STRING_ELT(Sfname,0)))+1, sizeof(char));
   strcpy(fname, CHAR(STRING_ELT(Sfname,0)));
   if (debug) INFO_message("filename %s\n", fname);
  
   /* get the what */
   PROTECT(Swhat = AS_CHARACTER(Swhat)); ++nprot;
   what = R_alloc(strlen(CHAR(STRING_ELT(Swhat,0)))+1, sizeof(char));
   strcpy(what, CHAR(STRING_ELT(Swhat,0)));
   if (debug) INFO_message("what %s\n", what);

   /* get the val */
   PROTECT(Sval = AS_CHARACTER(Sval)); ++nprot;
   val = R_alloc(strlen(CHAR(STRING_ELT(Sval,0)))+1, sizeof(char));
   strcpy(val, CHAR(STRING_ELT(Sval,0)));
   if (debug) INFO_message("val %s\n", val);

   /* get the cwd */
   if (!isNull(Scwd)) {
      PROTECT(Scwd = AS_CHARACTER(Scwd)); ++nprot;
      cwd = R_alloc(strlen(CHAR(STRING_ELT(Scwd,0)))+1, sizeof(char));
      strcpy(cwd, CHAR(STRING_ELT(Scwd,0)));
      if (debug) INFO_message("cwd %s\n", cwd);
   } 
   
   if (debug) INFO_message("Modifying %s\n", fname);
   if ((res = SUMA_ModifyName(fname, what, val, cwd))) {
      PROTECT(Rname = allocVector(STRSXP, 1)); ++nprot;
      SET_STRING_ELT(Rname, 0, mkChar(res)); 
      SUMA_free(res); res=NULL;
   } else {
      ERROR_message("Call to SUMA_ModifyName %s %s %s failed", fname, what, val);
   }
   UNPROTECT(nprot); 
   
   return(Rname);
}

int fill_mask(options_t * opts)
{
   THD_3dim_dataset * mset;
   int nvox;

ENTRY("fill_mask");

   if( opts->automask ) {
      if( opts->verb ) INFO_message("creating automask...");

      opts->mask = THD_automask(opts->inset);
      if( ! opts->mask ) {
         ERROR_message("failed to apply -automask");
         RETURN(1);
      }

      RETURN(0);
   }

   if( opts->mask_name ) {
      if( opts->verb )
         INFO_message("reading mask dset from %s...", opts->mask_name);

      mset = THD_open_dataset( opts->mask_name );
      if( ! mset ) ERROR_exit("cannot open mask dset '%s'", opts->mask_name);
      nvox = DSET_NVOX(opts->inset);
      if( DSET_NVOX(mset) != nvox ) {
         ERROR_message("mask does not have the same voxel count as input");
         RETURN(1);
      }

      /* fill mask array and mask_nxyz, remove mask dset */
      DSET_load(mset); CHECK_LOAD_ERROR(mset);

      opts->mask = THD_makemask(mset, 0, 1, 0);
      DSET_delete(mset);

      if( ! opts->mask ) {
         ERROR_message("cannot make mask from '%s'", opts->mask_name);
         RETURN(1);
      }

      if( opts->verb > 1 )
         INFO_message("have mask with %d voxels", nvox);
   }

   RETURN(0);
}

/*
 * A hole is defined as a connected set of zero voxels that does
 * not reach an edge.
 *
 * The core functionality was added to libmri.a in THD_mask_fill_holes.
 */
int fill_holes(THD_3dim_dataset * dset, int verb)
{
   short * sptr;     /* to for filling holes */
   byte  * bmask;    /* computed result */
   int     nfilled;
   int     nx, ny, nz, nvox, index, fill=0;

   ENTRY("fill_holes");

   bmask = THD_makemask(dset, 0, 1, 0); /* copy input as byte mask */
   nx = DSET_NX(dset);  ny = DSET_NY(dset);  nz = DSET_NZ(dset);
   nvox = DSET_NVOX(dset);

   /* created filled mask */
   nfilled = THD_mask_fill_holes(nx,ny,nz, bmask, verb);
   if( nfilled < 0 ) { ERROR_message("failed to fill holes");  RETURN(1); }

   /* apply to short volume */
   sptr = DBLK_ARRAY(dset->dblk, 0);
   for( index = 0; index < nvox; index++ )
      if( !sptr[index] && bmask[index] ) { fill++;  sptr[index] = 1; }

   if(verb>2) INFO_message("final check: fill=%d, nfilled=%d", fill, nfilled);

   RETURN(0);
}

void AFNI_coord_filer_setup( Three_D_View *im3d )
{
   char ename[32] , *eval ; int ic ;

ENTRY("AFNI_coord_filer_setup") ;

   if( !IM3D_OPEN(im3d) ) EXRETURN ;
   ic = AFNI_controller_index(im3d) ;
   if( ic < 0 || ic >= MAX_CONTROLLERS || fpc[ic] != NULL ) EXRETURN ;

   sprintf(ename,"AFNI_FILE_COORDS_%c",abet[ic]) ;
   eval = my_getenv(ename) ;
   if( eval == NULL || *eval == '\0' ){ fpc[ic] = NULL ; EXRETURN ; }

   if( strcmp(eval,"-") == 0 || strncmp(eval,"stdout",6) == 0 )
     fpc[ic] = stdout ;
   else {
     fpc[ic] = fopen( eval , "w" ) ;
     if( fpc[ic] == NULL ){
       ERROR_message("Unable to open file %s from %s",eval,ename) ;
       EXRETURN ;
     }
   }

   AFNI_receive_init( im3d , RECEIVE_VIEWPOINT_MASK ,
                             AFNI_filer_viewpoint_CB ,
                             im3d , "AFNI_filer_viewpoint_CB" ) ;

   INFO_message("Logging [%c] viewpoint changes to '%s'",abet[ic],eval) ;
   EXRETURN ;
}

/*
 * for each input dataset name
 *    open (check dims, etc.)
 *    dilate (zeropad, make binary, dilate, unpad, apply)
 *    fill list of bytemask datasets
 *
 * also, count total volumes
 */
int process_input_dsets(param_t * params)
{
   THD_3dim_dataset * dset, * dfirst=NULL;
   int                iset, nxyz;

   ENTRY("process_input_dsets");

   if( !params ) ERROR_exit("NULL inputs to PID");

   if( params->ndsets <= 0 ) {
      ERROR_message("process_input_dsets: no input datasets");
      RETURN(1);
   }

   /* allocate space for dsets array */
   params->dsets = (THD_3dim_dataset **)malloc(params->ndsets*
                                               sizeof(THD_3dim_dataset*));
   if( !params->dsets ) ERROR_exit("failed to allocate dset pointers");

   if( params->verb ) INFO_message("processing %d input datasets...",
                                   params->ndsets);
   
   /* warn user of dilations */
   if(params->verb && params->ndsets) {
      int pad = needed_padding(&params->IND);
      INFO_message("padding all datasets by %d (for dilations)", pad);
   }

   /* process the datasets */
   nxyz = 0;
   for( iset=0; iset < params->ndsets; iset++ ) {
      /* open and verify dataset */
      dset = THD_open_dataset(params->inputs[iset]);
      if( !dset ) ERROR_exit("failed to open mask dataset '%s'",
                             params->inputs[iset]);
      DSET_load(dset);  CHECK_LOAD_ERROR(dset);

      if( params->verb>1 ) INFO_message("loaded dset %s, with %d volumes",
                                        DSET_PREFIX(dset), DSET_NVALS(dset));

      if( nxyz == 0 ) { /* make an empty copy of the first dataset */
         nxyz = DSET_NVOX(dset);
         dfirst = EDIT_empty_copy(dset);
      }

      /* check for consistency in voxels and grid */
      if( DSET_NVOX(dset) != nxyz ) ERROR_exit("nvoxel mis-match");
      if( ! EQUIV_GRIDS(dset, dfirst) )
         WARNING_message("grid from dset %s does not match that of dset %s",
                         DSET_PREFIX(dset), DSET_PREFIX(dfirst));

      /* apply dilations to all volumes, returning bytemask datasets */
      params->dsets[iset] = apply_dilations(dset, &params->IND,1,params->verb);
      if( ! params->dsets[iset] ) RETURN(1);
   } 

   DSET_delete(dfirst); /* and nuke */

   RETURN(0);
}

float THD_dset_max(THD_3dim_dataset *dset, int scl) {
   float max,min;
   if (!THD_dset_minmax(dset, scl,&min, &max)) {
      ERROR_message("Could not get dset min max");
   }
   return(max);
}