My goal is to use MathLink to access the functions in sofa (Standards of Fundamental Astronomy)

So far I have had mixed success. I have very little C programing experience so I am hoping to get access to these algorithms via MathLink.

I have a very simple C code file that loads the sofa libraries:

#include "mathlink.h"
#include <stdio.h>
#include "sofa.h"
#include "sofam.h"

int main(int argc, char* argv[])
return MLMain(argc, argv);

Then, from the sofa documentation I find function definitions like:

double iauAnp(double a)

Which I then turn into inputs, in the sofa.tm file

:Function:       iauAnp
:Pattern:        iauAnp[i_Real]
:Arguments:      { i }
:ArgumentTypes:  { Real64 }
:ReturnType:     Real64

:Evaluate: iauAnp::usage = "iauAnp[x] Normalize angle into the range 0 <= a < 2pi."

I'm then able to compile this and load it into Mathematica via the Install[] function. So I know I can turn the crank once for a simple case.

The question: How do I write the lines for the sofa.tm file for the following function?

void iauC2i00a(double date1, double date2, double rc2i[3][3])

My understanding of the workings of this function is that date1 and date2 are inputs and the 2D array, rc2i is returned by reference.

I would like the Mathematica function to look like this:

rc2i = iauC2i00a[date1, date2]

1 Answer 1


Caution: All code below is untested written from the top of my head. It should only serve as help for the explanation and is most likely not compilable.

In your last case you need to do what is usually done with all library-calls: you write a wrapper function. This is basically a function which takes the input, prepares everything (if there is something), calls your library function, collects the result and sends everything back to Mathematica.

Your logic should work as follows. You define a template which calls your wrapper function. What you want is to send two double values in and get one double array out.

:Function:       MLiauC2i00a
:Pattern:        MLiauC2i00a[date1_?NumericQ, date2_?NumericQ]
:Arguments:      {N[date1], N[date2]}
:ArgumentTypes:  { Real64, Real64 }
:ReturnType:     Manual

Note that I have put Manual as return type which says: I will do the work of putting the result back to Mathematica by myself using MathLink functions.

void MLiauC2i00a(double date1, double date2) {
  double rc2i[3][3];
  iauC2i00a(date1, date2, rc2i);

  MLPutFunction(stdlink, "List", 3);
    MLPutReal64List(stdlink, rc2i[0], 3);
    MLPutReal64List(stdlink, rc2i[1], 3);
    MLPutReal64List(stdlink, rc2i[2], 3);

In the first part I define rc2i and call the library-function. Afterwards the result is stored in rc2i. To send it back to Mathematica, I build up the expression exactly as you would do it in Mathematica itself:

(* {rc2i[0],rc2i[1],rc2i[2]} *)

where the rc2i are then again lists of length 3 of real numbers. Now you can probably guess that the last argument to MLPutFunction is the number of list-elements we want to supply.

Needless to say, that you could call the MLPutReal64List functions in a loop if you have more than 3.

  • $\begingroup$ Thank you this works! (Except in the :Pattern: the second data should be a data2) I tried to edit it but an edit must be at least 6 letters. $\endgroup$
    – c186282
    Commented Jul 15, 2013 at 13:58
  • $\begingroup$ @c186282 I changed it. Nice that it helped you. $\endgroup$
    – halirutan
    Commented Jul 15, 2013 at 16:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.