[gnuastro-commits] master 3e76e1f 010/113: MakeProfiles --kernel builds 3D kernels also

[Mohammad Akhlaghi]

branch: master
commit 3e76e1f86984197df1830c50f00dd59ec8b47227
Author: Mohammad Akhlaghi <akhlaghi@gnu.org>
Commit: Mohammad Akhlaghi <akhlaghi@gnu.org>

    MakeProfiles --kernel builds 3D kernels also
    
    The `--kernel' option can now build a 3D kernel as well as a 2D kernel. A
    test has also been created to check this feature in `make check'.
---
 NEWS                                    |   9 +-
 bin/mkprof/main.h                       |   5 +
 bin/mkprof/mkprof.c                     | 114 +++++++-----
 bin/mkprof/oneprofile.c                 |  19 +-
 bin/mkprof/ui.c                         | 314 +++++++++++++++++++++++---------
 doc/gnuastro.texi                       |  62 +++++--
 tests/Makefile.am                       |   8 +-
 tests/mknoise/addnoise-3d.sh            |   2 +-
 tests/mkprof/{3dcat.sh => 3d-cat.sh}    |   2 +-
 tests/mkprof/{3dcat.txt => 3d-cat.txt}  |   0
 tests/mkprof/{3dcat.sh => 3d-kernel.sh} |   7 +-
 11 files changed, 367 insertions(+), 175 deletions(-)

diff --git a/NEWS b/NEWS
index dd65fa6..64cff8e 100644
--- a/NEWS
+++ b/NEWS
@@ -18,13 +18,12 @@ GNU Astronomy Utilities NEWS                          -*- 
outline -*-
   Gnuastro. The `--config' option can be used to specify default values for
   3D profiles.
 
-  MakeProfiles: the new `--kernel' option can make a kernel image without
-  the need to define a catalog. With this option, a catalog (or
-  accompanying background image) must not be given.
+  MakeProfiles: the new `--kernel' option can make a kernel image or cube
+  without the need to define a catalog.
 
   MakeProfiles: the new `--pc', `--cunit' and `--ctype' options can be used
-  to specify the PC matrix, CUNIT and CTYPE world coordinate system
-  keywords of the output FITS file.
+  along with `--crpix', `--crval' and `--cdelt' to fully specify the world
+  coordinate system (WCS) keywords of the output FITS file.
 
   MakeProfiles: can build a new raidal profile called `distance' (code
   7). The value this profile gives to each pixel is equal to its distance
diff --git a/bin/mkprof/main.h b/bin/mkprof/main.h
index 9d4e989..64b7ccb 100644
--- a/bin/mkprof/main.h
+++ b/bin/mkprof/main.h
@@ -185,6 +185,11 @@ struct mkprofparams
   char           *wcsheader;  /* The WCS header information for main img. */
   int            wcsnkeyrec;  /* The number of keywords in the WCS header.*/
   char       *mergedimgname;  /* Name of merged image.                    */
+  int                  nwcs;  /* for WCSLIB: no. coord. representations.  */
+  struct wcsprm        *wcs;  /* WCS information for this dataset.        */
+  size_t               ndim;  /* Number of dimensions (for `nomerged').   */
+                              /* We can't put it in `out' because it is   */
+                              /* meaning ful there.                       */
 };
 
 #endif
diff --git a/bin/mkprof/mkprof.c b/bin/mkprof/mkprof.c
index dbdf4ce..b1332e2 100644
--- a/bin/mkprof/mkprof.c
+++ b/bin/mkprof/mkprof.c
@@ -79,11 +79,17 @@ builtqueue_addempty(struct builtqueue **bq)
     error(EXIT_FAILURE, 0, "%s: allocating %zu bytes for `tbq'",
           __func__, sizeof *tbq);
 
-  /* Initialize some of the values. */
-  tbq->image=NULL;
-  tbq->numaccu=0;
-  tbq->accufrac=0.0f;
-  tbq->indivcreated=0;
+  /* Initialize the values (same order as in structure definition). */
+  tbq->id           = GAL_BLANK_SIZE_T;
+  tbq->ispsf        = 0;
+  tbq->overlaps     = 0;
+  tbq->image        = NULL;
+  tbq->overlap_i    = NULL;
+  tbq->overlap_m    = NULL;
+  tbq->func         = PROFILE_MAXIMUM_CODE;
+  tbq->indivcreated = 0;
+  tbq->numaccu      = 0;
+  tbq->accufrac     = 0.0f;
 
   /* Set its next element to the input bq and re-set the input bq. */
   tbq->next=*bq;
@@ -120,7 +126,7 @@ saveindividual(struct mkonthread *mkp)
 
   double *crpix;
   long os=p->oversample;
-  size_t i, ndim=p->out->ndim;
+  size_t i, ndim=p->ndim;
   struct builtqueue *ibq=mkp->ibq;
   char *filename, *jobname, *outdir=p->outdir;
 
@@ -202,7 +208,7 @@ mkprof_build_single(struct mkonthread *mkp, long *fpixel_i, 
long *lpixel_i,
 
   void *ptr;
   int needs_crop=0;
-  size_t i, ind, fits_i, ndim=p->out->ndim;
+  size_t i, ind, fits_i, ndim=p->ndim;
   size_t start_indiv[3], start_mrg[3], dsize[3], os=p->oversample;
 
   /* Make a copy of the main random number generator to use for this
@@ -228,7 +234,7 @@ mkprof_build_single(struct mkonthread *mkp, long *fpixel_i, 
long *lpixel_i,
   /* If we want a merged image, then a tile needs to be defined over the
      individual profile array and the output merged array to define the
      overlapping region. */
-  if(p->out->array)
+  if(p->out)
     {
       /* Note that `fpixel_o' and `lpixel_o' were in the un-oversampled
          image, they are also in the FITS coordinates. */
@@ -367,7 +373,7 @@ mkprof_build(void *inparam)
   struct mkonthread *mkp=(struct mkonthread *)inparam;
   struct mkprofparams *p=mkp->p;
 
-  size_t i, id, ndim=p->out->ndim;
+  size_t i, id, ndim=p->ndim;
   struct builtqueue *ibq, *fbq=NULL;
   double center[3], semiaxes[3], euler_deg[3];
   long fpixel_i[3], lpixel_i[3], fpixel_o[3], lpixel_o[3];
@@ -396,33 +402,43 @@ mkprof_build(void *inparam)
       if( p->f[id] == PROFILE_POINT )
         mkp->width[0]=mkp->width[1]=1;
       else
-        {
-          if(p->out->ndim==2)
+        switch(ndim)
+          {
+          case 2:
             gal_box_bound_ellipse(mkp->truncr, mkp->q[0]*mkp->truncr,
                                   p->p1[id], mkp->width);
-          else
-            {
-              euler_deg[0] = p->p1[id];
-              euler_deg[1] = p->p2[id];
-              euler_deg[2] = p->p3[id];
-              semiaxes[0]  = mkp->truncr;
-              semiaxes[1]  = mkp->truncr * mkp->q[0];
-              semiaxes[2]  = mkp->truncr * mkp->q[1];
-              gal_box_bound_ellipsoid(semiaxes, euler_deg, mkp->width);
-            }
-        }
+            break;
+
+          case 3:
+            euler_deg[0] = p->p1[id];
+            euler_deg[1] = p->p2[id];
+            euler_deg[2] = p->p3[id];
+            semiaxes[0]  = mkp->truncr;
+            semiaxes[1]  = mkp->truncr * mkp->q[0];
+            semiaxes[2]  = mkp->truncr * mkp->q[1];
+            gal_box_bound_ellipsoid(semiaxes, euler_deg, mkp->width);
+            break;
+
+          default:
+            error(EXIT_FAILURE, 0, "%s: a bug! please contact us at %s to "
+                  "address the issue. %zu is not recognized for `ndim'",
+                  __func__, PACKAGE_BUGREPORT, ndim);
+          }
 
 
       /* Get the overlapping pixels using the starting points (NOT
          oversampled). */
-      center[0]=p->x[id];
-      center[1]=p->y[id];
-      if(ndim==3) center[2]=p->z[id];
-      gal_box_border_from_center(center, ndim, mkp->width, fpixel_i,
-                                 lpixel_i);
-      memcpy(mkp->fpixel_i, fpixel_i, ndim*sizeof *fpixel_i);
-      ibq->overlaps = gal_box_overlap(mkp->onaxes, fpixel_i, lpixel_i,
-                                      fpixel_o, lpixel_o, ndim);
+      if(p->out)
+        {
+          center[0]=p->x[id];
+          center[1]=p->y[id];
+          if(ndim==3) center[2]=p->z[id];
+          gal_box_border_from_center(center, ndim, mkp->width, fpixel_i,
+                                     lpixel_i);
+          memcpy(mkp->fpixel_i, fpixel_i, ndim*sizeof *fpixel_i);
+          ibq->overlaps = gal_box_overlap(mkp->onaxes, fpixel_i, lpixel_i,
+                                          fpixel_o, lpixel_o, ndim);
+        }
 
 
       /* Build the profile if necessary. */
@@ -504,7 +520,7 @@ mkprof_write(struct mkprofparams *p)
          both the individual array and the final merged array, here we will
          use those to put the required profile pixels into the final
          array. */
-      if(ibq->overlaps && out->array)
+      if(ibq->overlaps && out)
         GAL_TILE_PO_OISET(float,float,ibq->overlap_i,ibq->overlap_m,1,0, {
             *o  = p->replace ? ( *i==0.0f ? *o : *i ) :  (*i + *o);
             sum += *i;
@@ -562,14 +578,18 @@ mkprof_write(struct mkprofparams *p)
 
   /* Write the final array to the output FITS image if a merged image is to
      be created. */
-  if(out->array)
+  if(out)
     {
       /* Get the current time for verbose output. */
       if(!p->cp.quiet) gettimeofday(&t1, NULL);
 
-      /* Write the final image into a FITS file with the requested type. */
+      /* Write the final image into a FITS file with the requested
+         type. Until now, we were using `p->wcs' for the WCS, but from now
+         on, will put it in `out' to also free it while freeing `out'. */
+      out->wcs=p->wcs;
       gal_fits_img_write_to_type(out, p->mergedimgname, NULL,
                                  PROGRAM_STRING, p->cp.type);
+      p->wcs=NULL;
 
       /* Clean up */
       gal_data_free(out);
@@ -582,11 +602,6 @@ mkprof_write(struct mkprofparams *p)
           free(jobname);
         }
     }
-
-  /* Even with no merged image, there might still be pointers in `out' that
-     need to be freed. */
-  else
-    gal_data_free(p->out);
 }
 
 
@@ -621,9 +636,9 @@ mkprof(struct mkprofparams *p)
   pthread_barrier_t b;
   struct mkonthread *mkp;
   gal_list_str_t *comments=NULL;
-  long *onaxes, os=p->oversample;
   size_t i, fi, *indexs, thrdcols;
-  size_t nb, ndim=p->out->ndim, nt=p->cp.numthreads;
+  long *onaxes=NULL, os=p->oversample;
+  size_t nb, ndim=p->ndim, nt=p->cp.numthreads;
 
   /* Allocate the arrays to keep the thread and parameters for each
      thread. Note that we only want nt-1 threads to do the
@@ -641,13 +656,18 @@ mkprof(struct mkprofparams *p)
   gal_threads_dist_in_threads(p->num, nt, &indexs, &thrdcols);
 
 
-  /* onaxes are sides of the image without over-sampling in FITS order. */
-  onaxes=gal_data_malloc_array(GAL_TYPE_LONG, ndim, __func__, "onaxes");
-  for(fi=0; fi < ndim; ++fi)
+  /* `onaxes' are size of the merged output image without over-sampling or
+     shifting in FITS order. When no output merged image is needed, we can
+     ignore it. */
+  if(p->out)
     {
-      i=ndim-fi-1;
-      onaxes[fi] = ( ( p->dsize[i] - 2 * p->shift[i] ) / os
-                     + 2 * p->shift[i]/os );
+      onaxes=gal_data_malloc_array(GAL_TYPE_LONG, ndim, __func__, "onaxes");
+      for(fi=0; fi < ndim; ++fi)
+        {
+          i=ndim-fi-1;
+          onaxes[fi] = ( ( p->dsize[i] - 2 * p->shift[i] ) / os
+                         + 2 * p->shift[i]/os );
+        }
     }
 
 
@@ -724,5 +744,5 @@ mkprof(struct mkprofparams *p)
   /* Clean up. */
   free(mkp);
   free(indexs);
-  free(onaxes);
+  if(onaxes) free(onaxes);
 }
diff --git a/bin/mkprof/oneprofile.c b/bin/mkprof/oneprofile.c
index 60698b8..8e20236 100644
--- a/bin/mkprof/oneprofile.c
+++ b/bin/mkprof/oneprofile.c
@@ -63,10 +63,10 @@ oneprofile_center_oversampled(struct mkonthread *mkp)
 
   double *dim;
   long os=p->oversample;
+  size_t i, id=mkp->ibq->id;
   double val, pixfrac, intpart;
-  size_t i, id=mkp->ibq->id, ndim=p->out->ndim;
 
-  for(i=0;i<ndim;++i)
+  for(i=0;i<p->ndim;++i)
     {
       dim = i==0 ? p->x : (i==1 ? p->y : p->z);
 
@@ -116,7 +116,7 @@ oneprofile_r_el(struct mkonthread *mkp)
   double c3=mkp->c[2],   s3=mkp->s[2];
   double x=mkp->coord[0], y=mkp->coord[1], z=mkp->coord[2];
 
-  switch(mkp->p->out->ndim)
+  switch(mkp->p->ndim)
     {
     case 2:
       /* The parenthesis aren't necessary, but help in readability and
@@ -136,8 +136,7 @@ oneprofile_r_el(struct mkonthread *mkp)
     default:
       error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to fix "
             "the problem. The value %zu is not recognized for "
-            "`mkp->p->out->ndim'", __func__, PACKAGE_BUGREPORT,
-            mkp->p->out->ndim);
+            "`mkp->p->ndim'", __func__, PACKAGE_BUGREPORT, mkp->p->ndim);
     }
 }
 
@@ -197,7 +196,7 @@ float
 oneprofile_randompoints(struct mkonthread *mkp)
 {
   double range[3], sum=0.0f;
-  size_t i, j, numrandom=mkp->p->numrandom, ndim=mkp->p->out->ndim;
+  size_t i, j, numrandom=mkp->p->numrandom, ndim=mkp->p->ndim;
 
   /* Set the range in each dimension. */
   for(i=0;i<ndim;++i)
@@ -322,7 +321,7 @@ oneprofile_center_pix_index(struct mkonthread *mkp)
 {
   double pixfrac, intpart;
   size_t *dsize=mkp->ibq->image->dsize;
-  size_t i, coord[3], ndim=mkp->ibq->image->ndim;
+  size_t i, coord[3], ndim=mkp->p->ndim;
 
   /* Find the coordinates of the center point. Note `mkp->center' is in
      FITS coordinates, while coord must be in C coordinates (to be used in
@@ -523,7 +522,7 @@ oneprof_set_prof_params(struct mkonthread *mkp)
 
   double sigma;
   int tp=p->tunitinp;
-  size_t id=mkp->ibq->id, ndim=p->out->ndim;
+  size_t id=mkp->ibq->id, ndim=p->ndim;
 
   /* Fill the most basic dimension and profile agnostic parameters. */
   mkp->brightness = pow( 10, (p->zeropoint - p->m[id]) / 2.5f );
@@ -563,7 +562,7 @@ oneprof_set_prof_params(struct mkonthread *mkp)
     default:
       error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to "
             "address the problem. The value `%zu' is not recognized for "
-            "`p->out->ndim'", __func__, PACKAGE_BUGREPORT, p->out->ndim);
+            "`ndim'", __func__, PACKAGE_BUGREPORT, ndim);
     }
 
 
@@ -705,7 +704,7 @@ oneprofile_make(struct mkonthread *mkp)
 
   double sum;
   float *f, *ff;
-  size_t i, dsize[3], ndim=p->out->ndim;
+  size_t i, dsize[3], ndim=p->ndim;
 
 
   /* Find the profile center in the over-sampled image in C
diff --git a/bin/mkprof/ui.c b/bin/mkprof/ui.c
index 2b9d2ec..6138fe7 100644
--- a/bin/mkprof/ui.c
+++ b/bin/mkprof/ui.c
@@ -289,8 +289,8 @@ ui_parse_kernel(struct argp_option *option, char *arg,
   long profcode;
   double *darray;
   gal_data_t *kernel;
-  size_t i, nc, numneeded=0;
-  char *c, *profile, *tailptr;
+  size_t i, nc, need=0;
+  char *c, *dstr, *profile, *tailptr;
   char *str, sstr[GAL_OPTIONS_STATIC_MEM_FOR_VALUES];
 
   /* We want to print the stored values. */
@@ -302,7 +302,16 @@ ui_parse_kernel(struct argp_option *option, char *arg,
 
       /* First write the profile function code into the output string. */
       nc=0;
-      nc += sprintf(sstr+nc, "%s,", ui_profile_name_write(kernel->status));
+      profile=ui_profile_name_write(kernel->status);
+      switch(kernel->flag)
+        {
+        case 2: nc += sprintf(sstr+nc, "%s,",    profile); break;
+        case 3: nc += sprintf(sstr+nc, "%s-3d,", profile); break;
+        default:
+          error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to "
+                "fix the problem. %u is not a recognized kernel "
+                "dimensionality", __func__, PACKAGE_BUGREPORT, kernel->flag);
+        }
 
       /* Write the values into a string. */
       for(i=0;i<kernel->size;++i)
@@ -329,13 +338,51 @@ ui_parse_kernel(struct argp_option *option, char *arg,
          the rest.*/
       c=arg;while(*c!='\0' && *c!=',') ++c;
       profile=arg;
-      arg = (*c=='\0') ? NULL : c+1;  /* `point' doesn't need any numbers. */
-      *c='\0';
+      arg = (*c=='\0') ? NULL : c+1;  /* the `point' profile doesn't need */
+      *c='\0';                        /* any numbers.                     */
+
 
       /* Read the parameters. */
       kernel=gal_options_parse_list_of_numbers(arg, filename, lineno);
       *(gal_data_t **)(option->value) = kernel;
 
+
+      /* All parameters must be positive. */
+      darray=kernel->array;
+      for(i=0;i<kernel->size;++i)
+        if(darray[i]<=0)
+          error(EXIT_FAILURE, 0, "value number %zu (%g) in the given list "
+                "of kernel parameters (`%s') is not acceptable. All "
+                "parameters to the `--kernel' option must be non-zero and "
+                "positive", i+1, darray[i], arg);
+
+
+      /* See if a 2D kernel is requested or a 3D kernel and keep the value
+         in `kernel->flag'. If no dimensionality is defined, then by
+         default, we'll assume it is 2D.*/
+      c=profile;while(*c!='\0' && *c!='-') ++c;
+      if(*c=='\0')
+        kernel->flag=2;
+      else
+        {
+          *c='\0';
+          dstr=c+1;
+          if( (dstr[1]!='d' && dstr[1]!='D') || dstr[2]!='\0')
+            error(EXIT_FAILURE, 0, "bad formatting in `--kernel' "
+                  "dimensionality. The dimensionality suffix must be either "
+                  "2d, 3d (not case sensitive). You have given `%s'", dstr);
+          switch(dstr[0])
+            {
+            case '2': kernel->flag=2; break;
+            case '3': kernel->flag=3; break;
+            default:
+              error(EXIT_FAILURE, 0, "only 2 or 3 dimensional kernels can "
+                    "currently be built, you have asked for a %c "
+                    "dimensional kernel", dstr[0]);
+            }
+        }
+
+
       /* Write the profile type code into `kernel->status'. If it starts
          with a digit, then the user might have given the code of the
          profile directly. In that case, parse the number. Otherwise,
@@ -356,15 +403,17 @@ ui_parse_kernel(struct argp_option *option, char *arg,
       else
         kernel->status=ui_profile_name_read(profile, 0);
 
+
       /* Make sure the number of parameters conforms with the profile. */
       switch(kernel->status)
         {
-        case PROFILE_SERSIC:        numneeded=3;     break;
-        case PROFILE_MOFFAT:        numneeded=3;     break;
-        case PROFILE_GAUSSIAN:      numneeded=2;     break;
-        case PROFILE_POINT:         numneeded=0;     break;
-        case PROFILE_FLAT:          numneeded=1;     break;
-        case PROFILE_CIRCUMFERENCE: numneeded=1;     break;
+        case PROFILE_SERSIC:        need = kernel->flag==2 ? 3 : 4;  break;
+        case PROFILE_MOFFAT:        need = kernel->flag==2 ? 3 : 4;  break;
+        case PROFILE_GAUSSIAN:      need = kernel->flag==2 ? 2 : 3;  break;
+        case PROFILE_POINT:         need = 0;                        break;
+        case PROFILE_FLAT:          need = kernel->flag==2 ? 1 : 2;  break;
+        case PROFILE_CIRCUMFERENCE: need = kernel->flag==2 ? 1 : 2;  break;
+        case PROFILE_DISTANCE:      need = kernel->flag==2 ? 1 : 2;  break;
         default:
           error_at_line(EXIT_FAILURE, 0, filename, lineno, "%s: a bug! "
                         "Please contact us at %s to correct the issue. "
@@ -372,13 +421,31 @@ ui_parse_kernel(struct argp_option *option, char *arg,
                         PACKAGE_BUGREPORT, kernel->status);
         }
 
+
       /* Make sure the number of parameters given are the same number that
          are needed. */
-      if( kernel->size != numneeded )
-        error_at_line(EXIT_FAILURE, 0, filename, lineno, "as a kernel, a "
-                      "`%s' profile needs %zu parameters, but %zu is given",
-                      ui_profile_name_write(kernel->status), numneeded,
-                      kernel->size);
+      if( kernel->size != need )
+        error_at_line(EXIT_FAILURE, 0, filename, lineno, "as a %uD kernel, "
+                      "a `%s' profile needs %zu parameters, but only %zu "
+                      "parameter%s given to `--kernel'", kernel->flag,
+                      ui_profile_name_write(kernel->status), need,
+                      kernel->size, kernel->size>1?"s are":" is");
+
+
+      /* If we want a 3D kernel, the last value (axis ratio) must be
+         positive and smaller and equal to 1. Note that we have already
+         checked for everything to be positive, so it can't be negative at
+         this point.
+
+         IMPORTANT NOTE: a `point' profile can have a kernel parameter list
+         of zero elements. So we need to make sure there are actually
+         numbers in the list also. */
+      if(kernel->flag==3 && kernel->size && darray[kernel->size-1]>1.0f)
+        error(EXIT_FAILURE, 0, "%g (last value in the list of kernel "
+              "parameters `%s') is interpretted as the kernel axis ratio "
+              "along the third dimension. Hence it must not be larger than 1",
+              darray[kernel->size-1], arg);
+
 
       /* Our job is done, return NULL. */
       return NULL;
@@ -633,12 +700,12 @@ ui_read_cols_2d(struct mkprofparams *p)
   int checkblank;
   char *colname=NULL, **strarr;
   gal_list_str_t *colstrs=NULL, *ccol;
+  size_t i, counter=0, coordcounter=0;
   gal_data_t *cols, *tmp, *corrtype=NULL;
-  size_t i, counter=0, coordcounter=0, ndim=p->out->ndim;
 
   /* The coordinate columns are a linked list of strings. */
   ccol=p->ccol;
-  for(i=0;i<ndim;++i)
+  for(i=0; i<p->ndim; ++i)
     {
       gal_list_str_add(&colstrs, ccol->v, 0);
       ccol=ccol->next;
@@ -693,6 +760,7 @@ ui_read_cols_2d(struct mkprofparams *p)
             }
           break;
 
+
         case 3:
           if(tmp->type==GAL_TYPE_STRING)
             {
@@ -727,12 +795,13 @@ ui_read_cols_2d(struct mkprofparams *p)
             }
           break;
 
+
         case 4:
           colname="radius (`rcol')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
           p->r=corrtype->array;
 
-          /* Check if there is no zero-radius profile. */
+          /* Check if there is no negative or zero-radius profile. */
           for(i=0;i<p->num;++i)
             if(p->r[i]<=0.0f)
               error(EXIT_FAILURE, 0, "%s: row %zu, the radius value %g is "
@@ -740,24 +809,35 @@ ui_read_cols_2d(struct mkprofparams *p)
                     i+1, p->r[i]);
           break;
 
+
         case 5:
           colname="index (`ncol')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
           p->n=corrtype->array;
           break;
 
+
         case 6:
           colname="position angle (`pcol')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
           p->p1=corrtype->array;
           break;
 
+
         case 7:
           colname="axis ratio (`qcol')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
           p->q1=corrtype->array;
+
+          /* Check if there is no negative or >1.0f axis ratio. */
+          for(i=0;i<p->num;++i)
+            if(p->q1[i]<=0.0f || p->q1[i]>1.0f)
+              error(EXIT_FAILURE, 0, "%s: row %zu, the axis ratio value %g "
+                    "is not acceptable. It has to be >0 and <=1", p->catname,
+                    i+1, p->q1[i]);
           break;
 
+
         case 8:
           colname="magnitude (`mcol')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
@@ -765,12 +845,21 @@ ui_read_cols_2d(struct mkprofparams *p)
           checkblank=0;          /* Magnitude can be NaN: to mask regions. */
           break;
 
+
         case 9:
           colname="truncation (`tcol')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
           p->t=corrtype->array;
+
+          /* Check if there is no negative or zero truncation radius. */
+          for(i=0;i<p->num;++i)
+            if(p->t[i]<=0.0f)
+              error(EXIT_FAILURE, 0, "%s: row %zu, the truncation radius "
+                    "value %g is not acceptable. It has to be larger than 0",
+                    p->catname, i+1, p->t[i]);
           break;
 
+
         /* If the index isn't recognized, then it is larger, showing that
            there was more than one match for the given criteria */
         default:
@@ -812,8 +901,8 @@ ui_read_cols_3d(struct mkprofparams *p)
   int checkblank;
   char *colname=NULL, **strarr;
   gal_list_str_t *colstrs=NULL, *ccol;
+  size_t i, counter=0, coordcounter=0;
   gal_data_t *cols, *tmp, *corrtype=NULL;
-  size_t i, counter=0, coordcounter=0, ndim=p->out->ndim;
 
   /* The 3D-specific columns are not mandatory in `args.h', so we need to
      check here if they are given or not before starting to read them. */
@@ -824,7 +913,7 @@ ui_read_cols_3d(struct mkprofparams *p)
 
   /* The coordinate columns are a linked list of strings. */
   ccol=p->ccol;
-  for(i=0;i<ndim;++i)
+  for(i=0; i<p->ndim; ++i)
     {
       gal_list_str_add(&colstrs, ccol->v, 0);
       ccol=ccol->next;
@@ -954,12 +1043,26 @@ ui_read_cols_3d(struct mkprofparams *p)
           colname="axis ratio 1 (`qcol')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
           p->q1=corrtype->array;
+
+          /* Check if there is no negative or >1.0f axis ratio. */
+          for(i=0;i<p->num;++i)
+            if(p->q1[i]<=0.0f || p->q1[i]>1.0f)
+              error(EXIT_FAILURE, 0, "%s: row %zu, the first axis ratio "
+                    "(`--qcol') value %g is not acceptable. It has to be "
+                    ">0 and <=1", p->catname, i+1, p->q1[i]);
           break;
 
         case 11:
           colname="axis ratio 2 (`q2col')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
           p->q2=corrtype->array;
+
+          /* Check if there is no negative or >1.0f axis ratio. */
+          for(i=0;i<p->num;++i)
+            if(p->q2[i]<=0.0f || p->q2[i]>1.0f)
+              error(EXIT_FAILURE, 0, "%s: row %zu, the second axis ratio "
+                    "(`--q2col') value %g is not acceptable. It has to be "
+                    ">0 and <=1", p->catname, i+1, p->q2[i]);
           break;
 
         case 12:
@@ -973,6 +1076,13 @@ ui_read_cols_3d(struct mkprofparams *p)
           colname="truncation (`tcol')";
           corrtype=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT32);
           p->t=corrtype->array;
+
+          /* Check if there is no negative or zero truncation radius. */
+          for(i=0;i<p->num;++i)
+            if(p->t[i]<=0.0f)
+              error(EXIT_FAILURE, 0, "%s: row %zu, the truncation radius "
+                    "value %g is not acceptable. It has to be larger than 0",
+                    p->catname, i+1, p->t[i]);
           break;
 
         /* If the index isn't recognized, then it is larger, showing that
@@ -1017,7 +1127,6 @@ ui_prepare_columns(struct mkprofparams *p)
 {
   double *karr;
   float r, n, t;
-  size_t ndim=p->out->ndim;
 
   /* If the kernel option was called, then we need to build a series of
      single element columns to create an internal catalog. */
@@ -1027,23 +1136,33 @@ ui_prepare_columns(struct mkprofparams *p)
       p->num=1;
 
       /* Allocate the necessary columns. */
-      p->x  = gal_data_malloc_array(GAL_TYPE_FLOAT64, 1, __func__, "p->x");
-      p->y  = gal_data_malloc_array(GAL_TYPE_FLOAT64, 1, __func__, "p->y");
-      p->f  = gal_data_malloc_array(GAL_TYPE_UINT8,   1, __func__, "p->f");
-      p->r  = gal_data_malloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->r");
-      p->n  = gal_data_malloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->n");
-      p->p1 = gal_data_malloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->p");
-      p->q1 = gal_data_malloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->q");
-      p->m  = gal_data_malloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->m");
-      p->t  = gal_data_malloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->t");
-
-      /* Set the values that need special consideration. */
+      p->x  = gal_data_calloc_array(GAL_TYPE_FLOAT64, 1, __func__, "p->x");
+      p->y  = gal_data_calloc_array(GAL_TYPE_FLOAT64, 1, __func__, "p->y");
+      p->f  = gal_data_calloc_array(GAL_TYPE_UINT8,   1, __func__, "p->f");
+      p->r  = gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->r");
+      p->n  = gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->n");
+      p->p1 = gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->p1");
+      p->q1 = gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->q1");
+      p->m  = gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->m");
+      p->t  = gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->t");
+      if(p->ndim==3)
+        {
+          p->z =gal_data_calloc_array(GAL_TYPE_FLOAT64, 1, __func__, "p->z");
+          p->p2=gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->p2");
+          p->p3=gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->p3");
+          p->q2=gal_data_calloc_array(GAL_TYPE_FLOAT32, 1, __func__, "p->q2");
+        }
+
+      /* For profiles that need a different number of input values. Note
+         that when a profile doesn't need a value, it will be ignored. */
+      karr=p->kernel->array;
       if(p->kernel->size)
         {
-          karr=p->kernel->array;
           r = karr[0];
           n = p->kernel->size==2 ? 0.0f : karr[1];
-          t = p->kernel->size==1 ? 1.0f : karr[ p->kernel->size - 1 ];
+          t = ( p->ndim==2
+                ? p->kernel->size==1 ? 1.0f : karr[ p->kernel->size - 1 ]
+                : p->kernel->size==1 ? 1.0f : karr[ p->kernel->size - 2 ] );
         }
       else r=n=t=0.0f;
 
@@ -1057,6 +1176,13 @@ ui_prepare_columns(struct mkprofparams *p)
       p->q1[0] = 1.0f;
       p->m[0]  = 0.0f;
       p->t[0]  = t;
+      if(p->ndim==3)
+        {
+          p->z[0]  = 0.0f;
+          p->p2[0] = 0.0f;
+          p->p3[0] = 0.0f;
+          p->q2[0] = p->kernel->size ? karr[ p->kernel->size - 1 ] : 0.0f;
+        }
     }
   else
     {
@@ -1064,20 +1190,20 @@ ui_prepare_columns(struct mkprofparams *p)
          dimensions in outputs are the same. There is no problem if it is
          more than `ndim'. In that case, the last values (possibly in
          configuration files) will be ignored.*/
-      if( gal_list_str_number(p->ccol) < ndim )
+      if( gal_list_str_number(p->ccol) < p->ndim )
         error(EXIT_FAILURE, 0, "%zu coordinate columns (calls to "
               "`--coordcol') given but output has %zu dimensions",
-              gal_list_str_number(p->ccol), p->out->ndim);
+              gal_list_str_number(p->ccol), p->ndim);
 
       /* Call the respective function. */
-      switch(ndim)
+      switch(p->ndim)
         {
         case 2: ui_read_cols_2d(p);   break;
         case 3: ui_read_cols_3d(p);   break;
         default:
           error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to "
-                "resolve the issue. %zu not recognized for `p->out->ndim'",
-                __func__, PACKAGE_BUGREPORT, ndim);
+                "resolve the issue. %zu not recognized for `p->ndim'",
+                __func__, PACKAGE_BUGREPORT, p->ndim);
         }
     }
 }
@@ -1093,7 +1219,7 @@ ui_prepare_columns(struct mkprofparams *p)
 static int
 ui_wcs_sanity_check(struct mkprofparams *p)
 {
-  size_t ndim=p->out->ndim;
+  size_t ndim=p->ndim;
 
   if(p->crpix)
     {
@@ -1161,14 +1287,10 @@ ui_prepare_wcs(struct mkprofparams *p)
   int status;
   struct wcsprm *wcs;
   char **cunit, **ctype;
-  size_t i, ndim=p->out->ndim;
+  size_t i, ndim=p->ndim;
   double *crpix, *crval, *cdelt, *pc;
 
 
-  /* If the WCS structure is already set, then return. */
-  if(p->out->wcs) return;
-
-
   /* Check and initialize the WCS information. If any of the necessary WCS
      parameters are missing, then don't build any WCS. */
   if( ui_wcs_sanity_check(p) ) return;
@@ -1182,7 +1304,7 @@ ui_prepare_wcs(struct mkprofparams *p)
 
   /* Allocate the memory necessary for the wcsprm structure. */
   errno=0;
-  wcs=p->out->wcs=malloc(sizeof *wcs);
+  wcs=p->wcs=malloc(sizeof *wcs);
   if(wcs==NULL)
     error(EXIT_FAILURE, errno, "%zu for wcs in preparewcs", sizeof *wcs);
 
@@ -1227,10 +1349,11 @@ static void
 ui_prepare_canvas(struct mkprofparams *p)
 {
   float *f, *ff;
+  gal_data_t *keysll;
   long width[2]={1,1};
   int status=0, setshift=0;
-  size_t i, ndim, nshift=0, *dsize=NULL;
   double truncr, semiaxes[3], euler_deg[3];
+  size_t i, nshift=0, *dsize=NULL, ndim_counter;
 
   /* If a background image is specified, then use that as the output
      image to build the profiles over. */
@@ -1248,26 +1371,33 @@ ui_prepare_canvas(struct mkprofparams *p)
 
       /* Read in the background image and its coordinates, note that when
          no merged image is desired, we just need the WCS information of
-         the background image. So `ndim==0' and what `dsize' points to is
-         irrelevant. */
+         the background image and the number of its dimensions. So
+         `ndim==0' and what `dsize' points to is irrelevant. */
       if(p->nomerged)
-        p->out=gal_data_alloc(NULL, GAL_TYPE_FLOAT32, 0, dsize, NULL,
-                              1, p->cp.minmapsize, NULL, NULL, NULL);
+        {
+          /* Get the number of the background image's dimensions. */
+          keysll=gal_data_array_calloc(1);
+          keysll->name="NAXIS";   keysll->type=GAL_TYPE_SIZE_T;
+          gal_fits_key_read(p->backname, p->backhdu, keysll, 0, 0);
+          p->ndim = *(size_t *)(keysll->array);
+          keysll->name=NULL;
+          gal_data_array_free(keysll, 1, 1);
+        }
       else
         {
           /* Read the image. */
           p->out=gal_fits_img_read_to_type(p->backname, p->backhdu,
                                            GAL_TYPE_FLOAT32,
                                            p->cp.minmapsize);
-          ndim=p->out->ndim;
+          p->ndim=p->out->ndim;
 
           /* If p->dsize was given as an option, free it. */
           if( p->dsize ) free(p->dsize);
 
           /* Write the size of the background image into `dsize'. */
-          p->dsize=gal_data_malloc_array(GAL_TYPE_SIZE_T, p->out->ndim,
-                                         __func__, "p->dsize");
-          for(i=0;i<p->out->ndim;++i) p->dsize[i] = p->out->dsize[i];
+          p->dsize=gal_data_malloc_array(GAL_TYPE_SIZE_T, p->ndim, __func__,
+                                         "p->dsize");
+          for(i=0;i<p->ndim;++i) p->dsize[i] = p->out->dsize[i];
 
           /* Set all pixels to zero if the user wanted a clear canvas. */
           if(p->clearcanvas)
@@ -1278,17 +1408,18 @@ ui_prepare_canvas(struct mkprofparams *p)
          there is no shifts. */
       p->oversample=1;
       if(p->shift) free(p->shift);
-      p->shift=gal_data_calloc_array(GAL_TYPE_SIZE_T, p->out->ndim,
-                                     __func__, "p->shift (1)");
+      p->shift=gal_data_calloc_array(GAL_TYPE_SIZE_T, p->ndim, __func__,
+                                     "p->shift (1)");
 
       /* Read the WCS structure of the background image. */
-      p->out->wcs=gal_wcs_read(p->backname, p->backhdu, 0, 0, &p->out->nwcs);
-
+      p->wcs=gal_wcs_read(p->backname, p->backhdu, 0, 0, &p->nwcs);
     }
   else
     {
       /* Get the number of dimensions. */
-      ndim=0; for(i=0;p->dsize[i]!=GAL_BLANK_SIZE_T;++i) ++ndim;
+      ndim_counter=0;
+      for(i=0;p->dsize[i]!=GAL_BLANK_SIZE_T;++i) ++ndim_counter;
+      p->ndim=ndim_counter;
 
       /* If any of the shift elements are zero, the others should be too!*/
       if(p->shift && p->shift[0] && p->shift[1])
@@ -1301,10 +1432,10 @@ ui_prepare_canvas(struct mkprofparams *p)
             }
 
           /* Make sure it has the same number of elements as naxis. */
-          if(ndim!=nshift)
+          if(p->ndim!=nshift)
             error(EXIT_FAILURE, 0, "%zu and %zu elements given to `--ndim' "
                   "and `--shift' respectively. These two numbers must be the "
-                  "same", ndim, nshift);
+                  "same", p->ndim, nshift);
         }
       else
         {
@@ -1327,7 +1458,7 @@ ui_prepare_canvas(struct mkprofparams *p)
                        half of it for the shift. */
                     setshift=1;
                     truncr = p->tunitinp ? p->t[i] : p->t[i] * p->r[i]/2;
-                    if(ndim==2)
+                    if(p->ndim==2)
                       gal_box_bound_ellipse(truncr, p->q1[i]*truncr,
                                             p->p1[i], width);
                     else
@@ -1347,7 +1478,7 @@ ui_prepare_canvas(struct mkprofparams *p)
                  shifts (from zero). So, we'll just free it and reset
                  it. */
               if(p->shift) free(p->shift);
-              p->shift=gal_data_calloc_array(GAL_TYPE_SIZE_T, p->out->ndim,
+              p->shift=gal_data_calloc_array(GAL_TYPE_SIZE_T, p->ndim,
                                              __func__, "p->shift (2)");
               if(setshift)
                 {
@@ -1359,43 +1490,40 @@ ui_prepare_canvas(struct mkprofparams *p)
 
       /* If shift has not been set until now, set it. */
       if(p->shift==NULL)
-        p->shift=gal_data_calloc_array(GAL_TYPE_SIZE_T, ndim, __func__,
+        p->shift=gal_data_calloc_array(GAL_TYPE_SIZE_T, p->ndim, __func__,
                                        "p->shift (3)");
 
       /* Prepare the sizes of the final merged image (if it is to be
          made). Note that even if we don't want a merged image, we still
          need its WCS structure. */
-      if(p->nomerged)
-        ndim=0;
-      else
+      if(p->nomerged==0)
         {
-          ndim=0;
+          ndim_counter=0;
           for(i=0;p->dsize[i]!=GAL_BLANK_SIZE_T;++i)
             {
               /* Count the number of dimensions. */
-              ++ndim;
+              ++ndim_counter;
 
               /* Correct dsize. */
               p->dsize[i] = (p->dsize[i]*p->oversample) + (2*p->shift[i]);
             }
           dsize = p->dsize;
-        }
 
-      /* Make the output structure. */
-      p->out=gal_data_alloc(NULL, GAL_TYPE_FLOAT32, ndim, dsize, NULL, 1,
-                            p->cp.minmapsize, NULL, NULL, NULL);
+          /* Make the output structure. */
+          p->out=gal_data_alloc(NULL, GAL_TYPE_FLOAT32, ndim_counter, dsize,
+                                NULL, 1, p->cp.minmapsize, NULL, NULL, NULL);
+        }
     }
 
 
   /* Make the WCS structure of the output data structure if it has not
      been set yet. */
-  ui_prepare_wcs(p);
+  if(p->wcs==NULL)
+    ui_prepare_wcs(p);
 
 
-  /* Set the name, comments and units of the output structure/file. Note
-     that when no merged image is to be created, the array in `p->out' will
-     be NULL.*/
-  if(p->out->array)
+  /* Set the name, comments and units of the final merged output. */
+  if(p->out)
     {
       if(p->out->name) free(p->out->name);
       gal_checkset_allocate_copy("Mock profiles", &p->out->name);
@@ -1409,9 +1537,9 @@ ui_prepare_canvas(struct mkprofparams *p)
      string to speed up the process: if we don't do it here, this process
      will be necessary on every individual profile's output. So it is much
      more efficient done once here. */
-  if(p->individual && p->out->wcs)
+  if(p->individual && p->wcs)
     {
-      status=wcshdo(WCSHDO_safe, p->out->wcs, &p->wcsnkeyrec, &p->wcsheader);
+      status=wcshdo(WCSHDO_safe, p->wcs, &p->wcsnkeyrec, &p->wcsheader);
       if(status)
         error(EXIT_FAILURE, 0, "wcshdo error %d: %s", status,
               wcs_errmsg[status]);
@@ -1425,10 +1553,10 @@ ui_prepare_canvas(struct mkprofparams *p)
 static void
 ui_finalize_coordinates(struct mkprofparams *p)
 {
+  size_t i, ndim=p->ndim;
   double *x=NULL, *y=NULL;
   uint8_t os=p->oversample;
-  size_t i, ndim=p->out->ndim;
-  double *cdelt=p->out->wcs->cdelt, *crpix=p->out->wcs->crpix;
+  double *cdelt=p->wcs->cdelt, *crpix=p->wcs->crpix;
 
   /* When the user specified RA and Dec columns, the respective values
      where stored in the `p->x' and `p->y' arrays. So before proceeding, we
@@ -1438,7 +1566,7 @@ ui_finalize_coordinates(struct mkprofparams *p)
       /* Note that we read the RA and Dec columns into the `p->x' and `p->y'
          arrays temporarily before. Here, we will convert them, free the old
          ones and replace them with the proper X and Y values. */
-      gal_wcs_world_to_img(p->out->wcs, p->x, p->y, &x, &y, p->num);
+      gal_wcs_world_to_img(p->wcs, p->x, p->y, &x, &y, p->num);
 
       /* If any conversions created a WCSLIB error, both the outputs will be
          set to NaN. */
@@ -1459,7 +1587,7 @@ ui_finalize_coordinates(struct mkprofparams *p)
      background image, oversample is set to 1. This is done here because
      the conversion of WCS to pixel coordinates needs to be done with the
      non-over-sampled image.*/
-  for(i=0;i<p->out->ndim;++i)
+  for(i=0;i<p->ndim;++i)
     {
       /* Oversampling has already been applied in `p->shift'. Also note
          that shift is in the C dimension ordring, while crpix is in FITS
@@ -1532,19 +1660,25 @@ ui_preparations(struct mkprofparams *p)
      over-written: */
   if(p->kernel)
     {
+      /* Set the necessary constants. */
       p->nomerged=1;
       p->psfinimg=0;
       p->individual=1;
-    }
+      p->ndim=p->kernel->flag;
 
-  /* Prepare the output canvas. */
-  ui_prepare_canvas(p);
+      /* Set the shift array. */
+      p->shift=gal_data_calloc_array(GAL_TYPE_SIZE_T, p->ndim,
+                                     __func__, "p->shift");
+    }
+  else
+    ui_prepare_canvas(p);
 
   /* Read in all the columns. */
   ui_prepare_columns(p);
 
   /* Read the (possible) RA/Dec inputs into X and Y for the builder.*/
-  ui_finalize_coordinates(p);
+  if(p->wcs)
+    ui_finalize_coordinates(p);
 
   /* Allocate the random number generator: */
   gsl_rng_env_setup();
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 21520f4..194eb85 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -15194,12 +15194,13 @@ Only build one kernel profile with the parameters 
given as the values to
 this option. The different values must be separated by a comma
 (@key{,}). The first value identifies the radial function of the profile,
 either through a string or through a number (see description of
-@option{--fcol} below). Each radial profile needs a different total number
-of parameters: S@'ersic and Moffat functions need 3 parameters: radial,
-S@'ersic index or Moffat @mymath{\beta}, and truncation radius. The
-Gaussian function needs two parameters: radial and truncation radius. The
-point function doesn't need any parameters and flat and circumference
-profiles just need one parameter (truncation radius).
+@option{--fcol} in @ref{MakeProfiles catalog}). Each radial profile needs a
+different total number of parameters: S@'ersic and Moffat functions need 3
+parameters: radial, S@'ersic index or Moffat @mymath{\beta}, and truncation
+radius. The Gaussian function needs two parameters: radial and truncation
+radius. The point function doesn't need any parameters and flat,
+circumference and distance profiles just need one parameter (truncation
+radius).
 
 The PSF or kernel is a unique (and highly constrained) type of profile: the
 sum of its pixels must be one, its center must be the center of the central
@@ -15207,24 +15208,54 @@ pixel (in an image with an odd number of pixels on 
each side), and commonly
 it is circular, so its axis ratio and position angle are one and zero
 respectively. Kernels are commonly necessary for various data analysis and
 data manipulation steps (for example see @ref{Convolve}, and
-@ref{NoiseChisel}. Because of this it is inconvenient to define a catalog
+@ref{NoiseChisel}). Because of this it is inconvenient to define a catalog
 with one row and many zero valued columns (for all the non-necessary
-parameters). Hence, with this option, it is possible to create a kernel
-with MakeProfiles without the need to create a catalog. Here are some
-examples:
+parameters). With this option, it is possible to easily create a kernel
+with MakeProfiles without a catalog. Here are some examples:
 
 @table @option
 @item --kernel=moffat,3,2.8,5
-A Moffat kernel with FWHM of 3 pixels, @mymath{\beta=2.8} which is
+A circular Moffat kernel with FWHM of 3 pixels, @mymath{\beta=2.8} which is
 truncated at 5 times the FWHM.
 
 @item --kernel=gaussian,2,3
-A Gaussian kernel with FWHM of 2 pixels and truncated at 3 times the FWHM.
+A circular Gaussian kernel with FWHM of 2 pixels and truncated at 3 times
+the FWHM.
+@end table
+
+This option may also be used to create a 3D kernel. To do that, to small
+modifications are necessary: add a @code{-3d} (or @code{-3D}) to the
+profile name (for example @code{moffat-3d}) and add a fraction (value
+between 0 and 1) to the end of the parameters for all profiles except
+@code{point}. The fractional value is the axis ratio in the third
+dimension. The main reason behind the fractional value is that commonly in
+3D astronomical datasets the third dimension doesn't have the same nature
+as the first and second.
+
+For example in IFU datacubes, the first and second dimensions are angular
+positions (like RA and Dec) but the third is in units of Angstroms for
+wavelength. Because of this different nature (which also affects the
+processing), it may be necessary for the kernel to have a different extent
+in that direction@footnote{So far, only cases where the third dimension
+extent must be less than the first two have been useful for us. If you need
+it to be larger, then please contact us and we will implement it.}. For
+example, let's have a look at the two examples above but in 3D:
+
+@table @option
+@item --kernel=moffat-3d,3,2.8,5,0.5
+An ellipsoid Moffat kernel with FWHM of 3 pixels, @mymath{\beta=2.8} which
+is truncated at 5 times the FWHM. The ellipsoid is circular in the first
+two dimensions, but in the third dimension its extent is half the first
+two.
+
+@item --kernel=gaussian-3d,2,3,1
+A spherical Gaussian kernel with FWHM of 2 pixels and truncated at 3 times
+the FWHM.
 @end table
 
 @item -x INT,INT
 @itemx --naxis=INT,INT
-The number of pixels along each dimension axis of the output in FITS
+The number of pixels along each dimension axis of f the output in FITS
 order. This is before over-sampling. For example if you call MakeProfiles
 with @option{--naxis=100,150 --oversample=5} (assuming no shift due for
 later convolution), then the final image size along the first axis will be
@@ -19448,7 +19479,10 @@ gal_fits_key_read_from_ptr(fptr, keysll, 0, 0);
 
 /* Use the values as you like... */
 
-/* Free all the allocated spaces. */
+/* Free all the allocated spaces. Note that `name' wasn't allocated
+   in this example, so we should explicitly set it to NULL before
+   calling `gal_data_array_free'.                                   */
+for(i=0;i<N-2;++i) keysll[i].name=NULL;
 gal_data_array_free(keysll, N, 1);
 @end example
 
diff --git a/tests/Makefile.am b/tests/Makefile.am
index 3df0714..a341a91 100644
--- a/tests/Makefile.am
+++ b/tests/Makefile.am
@@ -125,19 +125,21 @@ if COND_MKNOISE
   MAYBE_MKNOISE_TESTS = mknoise/addnoise.sh mknoise/addnoise-3d.sh
 
   mknoise/addnoise.sh: warp/warp_scale.sh.log
-  mknoise/addnoise-3d.sh: mkprof/3dcat.sh.log
+  mknoise/addnoise-3d.sh: mkprof/3d-cat.sh.log
 endif
 if COND_MKPROF
   MAYBE_MKPROF_TESTS = mkprof/mosaic1.sh mkprof/mosaic2.sh         \
   mkprof/mosaic3.sh mkprof/mosaic4.sh mkprof/radeccat.sh           \
-  mkprof/ellipticalmasks.sh mkprof/clearcanvas.sh mkprof/3dcat.sh
+  mkprof/ellipticalmasks.sh mkprof/clearcanvas.sh mkprof/3d-cat.sh \
+  mkprof/3d-kernel.sh
 
-  mkprof/3dcat.sh: prepconf.sh.log
+  mkprof/3d-cat.sh: prepconf.sh.log
   mkprof/mosaic1.sh: prepconf.sh.log
   mkprof/mosaic2.sh: prepconf.sh.log
   mkprof/mosaic3.sh: prepconf.sh.log
   mkprof/mosaic4.sh: prepconf.sh.log
   mkprof/radeccat.sh: prepconf.sh.log
+  mkprof/3d-kernel.sh: prepconf.sh.log
   mkprof/ellipticalmasks.sh: mknoise/addnoise.sh.log
   mkprof/clearcanvas.sh: mknoise/addnoise.sh.log
 endif
diff --git a/tests/mknoise/addnoise-3d.sh b/tests/mknoise/addnoise-3d.sh
index b04e193..ec010f4 100755
--- a/tests/mknoise/addnoise-3d.sh
+++ b/tests/mknoise/addnoise-3d.sh
@@ -28,7 +28,7 @@
 # those that do: for example SubtractSky and NoiseChisel will be better
 # tested on a larger image.
 prog=mknoise
-img=3dcat.fits
+img=3d-cat.fits
 execname=../bin/$prog/ast$prog
 
 
diff --git a/tests/mkprof/3dcat.sh b/tests/mkprof/3d-cat.sh
similarity index 97%
copy from tests/mkprof/3dcat.sh
copy to tests/mkprof/3d-cat.sh
index 4d39abd..459ce02 100755
--- a/tests/mkprof/3dcat.sh
+++ b/tests/mkprof/3d-cat.sh
@@ -21,7 +21,7 @@
 # file exists (basicchecks.sh is in the source tree).
 prog=mkprof
 execname=../bin/$prog/ast$prog
-cat=$topsrc/tests/$prog/3dcat.txt
+cat=$topsrc/tests/$prog/3d-cat.txt
 
 
 
diff --git a/tests/mkprof/3dcat.txt b/tests/mkprof/3d-cat.txt
similarity index 100%
rename from tests/mkprof/3dcat.txt
rename to tests/mkprof/3d-cat.txt
diff --git a/tests/mkprof/3dcat.sh b/tests/mkprof/3d-kernel.sh
similarity index 83%
rename from tests/mkprof/3dcat.sh
rename to tests/mkprof/3d-kernel.sh
index 4d39abd..d89f25e 100755
--- a/tests/mkprof/3dcat.sh
+++ b/tests/mkprof/3d-kernel.sh
@@ -1,4 +1,4 @@
-# Create a 3D image with MakeProfiles
+# Create a 3D kernel with MakeProfiles.
 #
 # See the Tests subsection of the manual for a complete explanation
 # (in the Installing gnuastro section).
@@ -21,7 +21,6 @@
 # file exists (basicchecks.sh is in the source tree).
 prog=mkprof
 execname=../bin/$prog/ast$prog
-cat=$topsrc/tests/$prog/3dcat.txt
 
 
 
@@ -37,7 +36,6 @@ cat=$topsrc/tests/$prog/3dcat.txt
 #
 #   - Catalog doesn't exist (problem in tarball release).
 if [ ! -f $execname ]; then echo "$execname not created."; exit 77; fi
-if [ ! -f $cat      ]; then echo "$cat does not exist.";   exit 77; fi
 
 
 
@@ -45,4 +43,5 @@ if [ ! -f $cat      ]; then echo "$cat does not exist.";   
exit 77; fi
 
 # Actual test script
 # ==================
-$execname $cat --config=.gnuastro/astmkprof-3d.conf --oversample=1
+$execname $cat --kernel=gaussian-3d,2,3,0.5 --oversample=1    \
+          --output=3d-kernel.fits