NcmModel tutorial

Curve fitting tutorial

Author

NumCosmo developers

Overview

This tutorial defines an abstract curve model in C and fits concrete implementations.

NumCosmo separates an abstract model, which declares only the fundamental quantities of a phenomenon, from the implementations that choose a functional form and a parametrization for them. The separation gives four things:

  • Analysis code can be written against the abstract interface, without reference to any particular implementation.
  • The fundamental quantities of the model stay distinct from the choices made when implementing them.
  • An implementation can be changed or replaced without affecting the rest of the codebase.
  • Separate groups can develop implementations independently.

Models can also be written in Python; that is covered in a separate tutorial.

Our abstract model

The phenomenon here is a curve \(f(x)\). Its only fundamental quantities are the curve itself and its domain \([x_l,\;x_u]\): no functional form and no parametrization, since those belong to the implementations.

See the GObject tutorial for the object framework.

The abstract NcCurve type subclasses NcmModel. nc_curve_f evaluates \(f(x)\); xl and xu define the domain.

Abstract model header

The header defines the virtual method, class structure, properties, and registration functions:

Virtual method type

This line determines the prototype of our virtual function. While this step is not strictly necessary – since one could use this prototype directly in the class structure – defining a type improves readability.

typedef gdouble (*NcCurveF) (NcCurve *curve, const gdouble x);

Object class

In the GObject framework, virtual functions are defined as function pointers in the class structure. In this simple example, there is a single virtual function, which we call f using the type defined above. The class structure is defined as follows:

struct _NcCurveClass
{
  /*< private >*/
  NcmModelClass parent_class;
  NcCurveF f;
};

Implementation enumerator

The NcmModel framework is designed to support complex models that may involve several virtual methods. Some of these methods may not be required by every implementation. For this reason in introduced the implementation enumerator, this way, the object that use the model may ask if that model implements a specific virtual method. In our example we have a single virtual method, so the enumerator is given by:

typedef enum _NcCurveImpl
{
  NC_CURVE_IMPL_f = 0,
} NcCurveImpl;

The enumerator above states that the first bit in the implementation enumerator controls if the child implements this method.

Defining the class id prototype

In order to be able to analyze different models, we need an object to aggregate all necessary models. In NumCosmo this task is performed by the NcmMSet (model set) object. We use the NCM_MSET_MODEL_DECLARE_ID macro to declare the id function nc_curve_id in the NcmMSet framework.

NCM_MSET_MODEL_DECLARE_ID (nc_curve);

Method prototypes

Finally we declare the two methods necessary to close the implementation (apart from the default ones).

void nc_curve_set_f_impl (NcCurveClass *curve_class, NcCurveF f);
gdouble nc_curve_f (NcCurve *curve, const gdouble x);

The first one just declare the class method that will be used to set the virtual method f. It is necessary because it has two jobs, it sets the function pointer to the member f of the class structure and updates the implementation flags to assert that child implements this method.

Abstract model source

We include the full C file below, however, first we are going to highlight the main steps:

Private properties

The first declaration in our code is the definition of the private struct

struct _NcCurvePrivate
{
  gdouble xl;
  gdouble xu;
};

In principle we could have added these members directly to the NcCurve struct, however, this would allow users to modify it directly. By hiding it in the compilation unit we avoid this, the only way to modify these values will be through the exported methods. This is important since these methods will perform consistency checks and call all necessary hooks that update the model.

As a rule of thumb, if the access to the struct member is not part of computationally heavy part of the code (here these members are likely to be set just once in the beginning and stay constant throughout all computations), then it just me hidden in the private struct for the reasons described above. Nevertheless, if the member will be accessed constantly during the computations we would attribute it to the object struct and include a inline method to access it, this technique will be discussed in detail in a future tutorial.

Property enumerator

The property enumerator is the part of the GObject framework used to define properties. Here the PROP_0 member is required by the GObject implementation and do not correspond to any actual property. The following members PROP_XL and PROP_XU correspond respectively to \(x_l\) and \(x_u\). Finally, PROP_SIZE* is a convenience member (that we use by convention in NumCosmo) used to count the number of properties, in this case PROP_SIZE == 3 (an enumerator always assign 0 to its first member and +1 to the subsequent, this can be changes by explicitly assigning values).

enum
{
  PROP_0,
  PROP_XL,
  PROP_XU,
  PROP_SIZE,
};

Object defining macro

The first step in creating the object is to call the GObject macro G_DEFINE_ABSTRACT_TYPE. This macro creates an abstract object (that cannot be instantiated) that is a child of NcmModel.

G_DEFINE_ABSTRACT_TYPE (NcCurve, nc_curve, NCM_TYPE_MODEL);

This macro assumes that two functions will be implemented: - nc_curve_init: this function should make all necessary initialization to the object that happens before any property is set. This usually consists of assign null or zero value for all struct members. - nc_curve_class_init: the class method must define all object’s properties and virtual methods.

Object initialization

In this example, the initialization function nc_curve_init has two tasks, getting the pointer to its private struct and assign the initial values for the properties. Our convention is that if both xl and xu are equal to zero then the object is in the initial state.

static void
nc_curve_init (NcCurve *curve)
{
  curve->priv = G_TYPE_INSTANCE_GET_PRIVATE (curve, NC_TYPE_CURVE, NcCurvePrivate);
  curve->priv->xl = 0.0;
  curve->priv->xu = 0.0;
}

Set/Get hooks

In the GObject framework all properties are modified through the set/get hooks. Here we have a two properties, for both we delegate the set/get job to methods implemented below. We could implement it directly inside the enumerators below, however, implementing it in separated functions improve readability.

static void
_nc_curve_set_property (GObject *object, guint prop_id, const GValue *value, GParamSpec *pspec)
{
  NcCurve *curve = NC_CURVE (object);
  g_return_if_fail (NC_IS_CURVE (object));

  switch (prop_id)
  {
    case PROP_XL:
      nc_curve_set_xl (curve, g_value_get_double (value));
      break;
    case PROP_XU:
      nc_curve_set_xu (curve, g_value_get_double (value));
      break;
    default:
      G_OBJECT_WARN_INVALID_PROPERTY_ID (object, prop_id, pspec);
      break;
  }
}

static void
_nc_curve_get_property (GObject *object, guint prop_id, GValue *value, GParamSpec *pspec)
{
  NcCurve *curve = NC_CURVE (object);
  g_return_if_fail (NC_IS_CURVE (object));

  switch (prop_id)
  {
    case PROP_XL:
      g_value_set_double (value, nc_curve_get_xl (curve));
      break;
    case PROP_XU:
      g_value_set_double (value, nc_curve_get_xu (curve));
      break;
    default:
      G_OBJECT_WARN_INVALID_PROPERTY_ID (object, prop_id, pspec);
      break;
  }
}

Defining the id function

The NcmMSet framework provides the macro NCM_MSET_MODEL_REGISTER_ID that should be used in the object code to create the nc_curve_id function. This function will return the unique ID defined in the NcmMSet type list.

NCM_MSET_MODEL_REGISTER_ID (nc_curve, NC_TYPE_CURVE);

The class initialization hook

Attributing the Set/Get hooks

In the class initialization hook one must define the object related hooks, as the set/get functions. Since we are using the NcmModel framework, we need to attribute the set/get hooks to the NcmModelClass struct. This is necessary since the framework takes care in organizing all model parameters and only the additional properties are handled by the set/get hooks defined here.

  GObjectClass* object_class = G_OBJECT_CLASS (klass);
  NcmModelClass *model_class = NCM_MODEL_CLASS (klass);

  ...
  
  model_class->set_property = &_nc_curve_set_property;
  model_class->get_property = &_nc_curve_get_property;

The next step is to initialize the model parameters.

ncm_model_class_set_name_nick (model_class, "Curve-f", "NcCurve");
ncm_model_class_add_params (model_class, 0, 0, PROP_SIZE);

These functions set the model name and nick and tells the framework that it requires 0 scalar parameters and 0 vector parameters and it has additionally PROP_SIZE properties. It is usually the case that the abstract model has no parameters, this happens because the parameters are usually connected to the specific implementations and not to the abstract definitions.

Registering the ID in NcmMSet

  ncm_mset_model_register_id (model_class,
                              "NcCurve",
                              "Curve model.",
                              NULL,
                              FALSE,
                              NCM_MSET_MODEL_MAIN);

This functions call takes care in registering the NcCurve type in the NcmMSet type list. This is necessary to set/get objects of this type from a NcmMSet.

Registering the properties

As discussed above, the properties \(x_l\) and \(x_u\) will not be part of the model parameters and will be taken as fixed numbers, in this case, we include them as usual GObject properties.

  g_object_class_install_property (object_class,
                                   PROP_XL,
                                   g_param_spec_double ("xl",
                                                        NULL,
                                                        "x lower bound",
                                                        -G_MAXDOUBLE, G_MAXDOUBLE, 0.0,
                                                        G_PARAM_READWRITE | G_PARAM_CONSTRUCT | G_PARAM_STATIC_NAME | G_PARAM_STATIC_BLURB));
  g_object_class_install_property (object_class,
                                   PROP_XU,
                                   g_param_spec_double ("xu",
                                                        NULL,
                                                        "x upper bound",
                                                        -G_MAXDOUBLE, G_MAXDOUBLE, 1.0,
                                                        G_PARAM_READWRITE | G_PARAM_CONSTRUCT | G_PARAM_STATIC_NAME | G_PARAM_STATIC_BLURB));

Checking parameters implementation

The function call below asserts that all required scalar and vector parameters were implemented correctly. In this case, none were required.

  ncm_model_class_check_params_info (model_class);

Setting a default implementation for the virtual function

This line sets the default implementation for the virtual function f.

  curve_class = &_nc_curve_f;

Methods

Our implementation of the virtual function f here is just a placeholder that raises an error if called.

static gdouble 
_nc_curve_f (NcCurve *curve, const gdouble x)
{
  g_error ("nc_curve_f: model `%s' does not implement this function.", 
           G_OBJECT_TYPE_NAME (curve)); 
  return 0.0;
}

We use the NCM_MODEL_SET_IMPL_FUNC to generate the code for the implementation function nc_curve_set_f_impl. This class method, when called by the child implementations will attribute the function to the class struct and enable the implementation flag for this method.

NCM_MODEL_SET_IMPL_FUNC (NC_CURVE, NcCurve, nc_curve, NcCurveF, f)

Example sources

Source

/***************************************************************************
 *            nc_curve.c
 *
 *  Tue July 11 16:30:26 2017
 *  Copyright  2017  Sandro Dias Pinto Vitenti
 *  <vitenti@uel.br>
 ****************************************************************************/
/*
 * nc_curve.c
 * Copyright (C) 2017 Sandro Dias Pinto Vitenti <vitenti@uel.br>
 *
 * numcosmo is free software: you can redistribute it and/or modify it
 * under the terms of the GNU General Public License as published by the
 * Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * numcosmo is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 * See the GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License along
 * with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

/**
 * NcCurve:
 *
 * Abstract curve model with domain bounds and a virtual f(x) method.
 */

#ifdef HAVE_CONFIG_H
#  include "config.h"
#endif /* HAVE_CONFIG_H */
#include "build_cfg.h"

#include "nc_curve.h"

/* Domain bounds. */
typedef struct _NcCurvePrivate
{
  gdouble xl;
  gdouble xu;
} NcCurvePrivate;

/* Object properties. */
enum
{
  PROP_0,
  PROP_XL,
  PROP_XU,
  PROP_SIZE,
};

/*
 * Define the abstract GType.
 * everything necessary to define an object in the GLib type system. The
 * ABSTRACT version of the macro creates a GType that cannot be instantiated
 * this means that to use this object we *must* define a child.
 *
 * The last argument provides the GType of the parent object.
 *
 */
G_DEFINE_ABSTRACT_TYPE (NcCurve, nc_curve, NCM_TYPE_MODEL);

static void
nc_curve_init (NcCurve *curve)
{
  /*
   * The first step is the creation of the private structure, all allocation and
   * de-allocation is automatically performed by the GObject framework.
   */
  curve->priv = G_TYPE_INSTANCE_GET_PRIVATE (curve, NC_TYPE_CURVE, NcCurvePrivate);

  /*
   * Here we initialize all structure member to null/zero.
   *
   */
  curve->priv->xl = 0.0;
  curve->priv->xu = 0.0;
}

/*
 * Here we call all property accessors in order to set property values.
 *
 */
static void
_nc_curve_set_property (GObject *object, guint prop_id, const GValue *value, GParamSpec *pspec)
{
  NcCurve *curve = NC_CURVE (object);

  g_return_if_fail (NC_IS_CURVE (object));

  switch (prop_id)
  {
    case PROP_XL:
      nc_curve_set_xl (curve, g_value_get_double (value));
      break;
    case PROP_XU:
      nc_curve_set_xu (curve, g_value_get_double (value));
      break;
    default:
      G_OBJECT_WARN_INVALID_PROPERTY_ID (object, prop_id, pspec);
      break;
  }
}

/*
 * Here we call all property accessors in order to get property values.
 *
 */
static void
_nc_curve_get_property (GObject *object, guint prop_id, GValue *value, GParamSpec *pspec)
{
  NcCurve *curve = NC_CURVE (object);

  g_return_if_fail (NC_IS_CURVE (object));

  switch (prop_id)
  {
    case PROP_XL:
      g_value_set_double (value, nc_curve_get_xl (curve));
      break;
    case PROP_XU:
      g_value_set_double (value, nc_curve_get_xu (curve));
      break;
    default:
      G_OBJECT_WARN_INVALID_PROPERTY_ID (object, prop_id, pspec);
      break;
  }
}

/* Release GObject-owned resources. */
static void
_nc_curve_dispose (GObject *object)
{
  G_OBJECT_CLASS (nc_curve_parent_class)->dispose (object);
}

/* Release non-GObject resources. */
static void
_nc_curve_finalize (GObject *object)
{
  G_OBJECT_CLASS (nc_curve_parent_class)->finalize (object);
}

/*
 * Registry the ID in NumCosmo model system.
 *
 */
NCM_MSET_MODEL_REGISTER_ID (nc_curve, NC_TYPE_CURVE);

/*
 * Prototype of the default implementation.
 *
 */
static gdouble _nc_curve_f (NcCurve *curve, const gdouble x);

/*
 * Define properties, parameters, and the default virtual method.
 */
static void
nc_curve_class_init (NcCurveClass *klass)
{
  GObjectClass *object_class = G_OBJECT_CLASS (klass);
  NcmModelClass *model_class = NCM_MODEL_CLASS (klass);

  /*
   * Tells the type system that we have a private struct
   *
   */
  g_type_class_add_private (klass, sizeof (NcCurvePrivate));

  /*
   * The NcmModel class takes cares of the parameters, thus,
   * the set/get functions above must be set in NcmModelClass
   * structure.
   *
   */
  model_class->set_property = &_nc_curve_set_property;
  model_class->get_property = &_nc_curve_get_property;

  /*
   * The other functions are the usual GObject hooks and
   * must be assigned to GObjectClass.
   *
   */
  object_class->dispose  = &_nc_curve_dispose;
  object_class->finalize = &_nc_curve_finalize;

  /*
   * First we set the model nick `Curve-f' and its name
   *
   */
  ncm_model_class_set_name_nick (model_class, "Curve-f", "NcCurve");

  /*
   * Now we inform that we have 0 scalar parameters and 0 vector parameters.
   * Note that here we should include parameters that *all* implementations
   * would share, here we have none. Finallt, with the last argument we assert
   * that we have PROP_SIZE == 2 additional properties.
   *
   */
  ncm_model_class_add_params (model_class, 0, 0, PROP_SIZE);

  /*
   * Next step is to register the object in NumCosmo's model set object class.
   * We don't include a long description (first NULL below). The last arguments
   * FALSE and NCM_MSET_MODEL_MAIN determines that this object in non-stackable
   * and it is a MAIN model, both concepts will be discussed in an advanced example.
   *
   */
  ncm_mset_model_register_id (model_class,
                              "NcCurve",
                              "Curve model.",
                              NULL,
                              FALSE,
                              NCM_MSET_MODEL_MAIN);


  /*
   * The property PROP_XL is allowed in range [-G_MAXDOUBLE, G_MAXDOUBLE]
   * and its default value is 0.0. This property can be set only during
   * the object construction.
   *
   */
  g_object_class_install_property (object_class,
                                   PROP_XL,
                                   g_param_spec_double ("xl",
                                                        NULL,
                                                        "x lower bound",
                                                        -G_MAXDOUBLE, G_MAXDOUBLE, 0.0,
                                                        G_PARAM_READWRITE | G_PARAM_CONSTRUCT | G_PARAM_STATIC_NAME | G_PARAM_STATIC_BLURB));

  /*
   * The property PROP_XU is allowed in range [-G_MAXDOUBLE, G_MAXDOUBLE]
   * and its default value is 1.0. This property can be set only during
   * the object construction.
   *
   */
  g_object_class_install_property (object_class,
                                   PROP_XU,
                                   g_param_spec_double ("xu",
                                                        NULL,
                                                        "x upper bound",
                                                        -G_MAXDOUBLE, G_MAXDOUBLE, 1.0,
                                                        G_PARAM_READWRITE | G_PARAM_CONSTRUCT | G_PARAM_STATIC_NAME | G_PARAM_STATIC_BLURB));

  /*
   * Check for errors in parameters initialization.
   *
   */
  ncm_model_class_check_params_info (model_class);

  klass->f = &_nc_curve_f;
}

/*
 * Default implementation. It raises an error if the method is called
 * and the child didn't implement it.
 *
 */
static gdouble
_nc_curve_f (NcCurve *curve, const gdouble x)
{
  g_error ("nc_curve_f: model `%s' does not implement this function.",
           G_OBJECT_TYPE_NAME (curve));

  return 0.0;
}

/*
 * This creates a function to be used by the children to
 * implement this method.
 *
 */
/**
 * nc_curve_set_f_impl: (skip)
 * @curve_class: a #NcCurveClass
 * @f: function $f(x)$
 *
 * Sets the implementation of the curve $f(x)$.
 *
 */
NCM_MODEL_SET_IMPL_FUNC (NC_CURVE, NcCurve, nc_curve, NcCurveF, f)

/*
 * This defines a generic constructor for the subclasses
 *
 */
/**
 * nc_curve_new_from_name:
 * @curve_name: #NcCurve child type name
 *
 * Creates a new #NcCurve of the type described by @curve_name.
 *
 * Returns: (transfer full): a new #NcCurve
 */
NcCurve *
nc_curve_new_from_name (const gchar *curve_name)
{
  /*
   * We use the serialization object to transform a string into an instance.
   * This also allows us to chose the parameters through the string.
   *
   */
  GObject *obj = ncm_serialize_global_from_string (curve_name);

  /*
   * This gets the GType of this new instance.
   *
   */
  GType curve_type = G_OBJECT_TYPE (obj);

  /*
   * Check if the string represent an actual child of #NcCurve.
   *
   */
  if (!g_type_is_a (curve_type, NC_TYPE_CLUSTER_MASS))
    g_error ("nc_curve_new_from_name: NcCurve `%s' do not descend from `%s'.",
             curve_name, g_type_name (NC_TYPE_CURVE));

  /*
   * Returns the correct cast.
   *
   */
  return NC_CURVE (obj);
}

/*
 * The reference increasing function
 *
 */

/**
 * nc_curve_ref:
 * @curve: a #NcCurve
 *
 * Increase reference count by one.
 *
 * Returns: (transfer full): @curve.
 */
NcCurve *
nc_curve_ref (NcCurve *curve)
{
  return g_object_ref (curve);
}

/*
 * The reference decreasing function
 *
 */

/**
 * nc_curve_free:
 * @curve: a #NcCurve
 *
 * Decrease reference count by one.
 *
 */
void
nc_curve_free (NcCurve *curve)
{
  g_object_unref (curve);
}

/*
 * This function decreases the reference count
 * by one only if *curve != NULL, and in that case
 * it sets *curve to NULL after decreasing the
 * reference count. It is useful to use these
 * functions in the dispose hooks.
 *
 *
 */

/**
 * nc_curve_clear:
 * @curve: a #NcCurve
 *
 * Decrease reference count by one if *@curve != NULL
 * and sets @curve to NULL.
 *
 */
void
nc_curve_clear (NcCurve **curve)
{
  g_clear_object (curve);
}

/*
 * Below we implement the accessor functions.
 *
 */

/**
 * nc_curve_set_xl:
 * @curve: a #NcCurve
 * @xl: new $x$ lower bound
 *
 * Sets $x$ lower bound to @xl.
 *
 */
void
nc_curve_set_xl (NcCurve *curve, const gdouble xl)
{
  /*
   * Checking if the object is still in the unintitalized state,
   * if that's the case just assign the value to the private
   * struct.
   *
   */
  if ((curve->priv->xl == 0.0) && (curve->priv->xl == curve->priv->xu))
  {
    curve->priv->xl = xl;
  }
  else
  {
    /*
     * Otherwise assert that xl will be less than xu.
     *
     */
    g_assert_cmpfloat (xl, <, curve->priv->xu);
    curve->priv->xl = xl;
  }
}

/**
 * nc_curve_set_xu:
 * @curve: a #NcCurve
 * @xu: new $x$ upper bound
 *
 * Sets $x$ upper bound to @xu.
 *
 */
void
nc_curve_set_xu (NcCurve *curve, const gdouble xu)
{
  /*
   * Checking if the object is still in the unintitalized state,
   * if that's the case just assign the value to the private
   * struct.
   *
   */
  if ((curve->priv->xl == 0.0) && (curve->priv->xl == curve->priv->xu))
  {
    curve->priv->xu = xu;
  }
  else
  {
    /*
     * Otherwise assert that xl will be less than xu.
     *
     */
    g_assert_cmpfloat (curve->priv->xl, <, xu);
    curve->priv->xu = xu;
  }
}

/**
 * nc_curve_get_xl:
 * @curve: a #NcCurve
 *
 * Sets $x$ lower bound.
 *
 * Returns: $x_l$.
 */
gdouble
nc_curve_get_xl (NcCurve *curve)
{
  return curve->priv->xl;
}

/**
 * nc_curve_get_xu:
 * @curve: a #NcCurve
 *
 * Sets $x$ lower bound.
 *
 * Returns: $x_l$.
 */
gdouble
nc_curve_get_xu (NcCurve *curve)
{
  return curve->priv->xu;
}

/*
/**
 * nc_curve_f: (virtual f)
 * @curve: a #NcCurve
 * @x: $x$
 *
 * Computes $f(x)$.
 *
 * Returns: the value of $f(x)$.
 */
gdouble
nc_curve_f (NcCurve *curve, const gdouble x)
{
  return NC_CURVE_GET_CLASS (curve)->f (curve, x);
}