source: issm/trunk-jpl/src/modules/Xy2ll/Xy2ll.cpp@ 12032

/*\file Xy2ll.c
 *\brief: x/y to lat/long coordinate mex module.
 */
#include "./Xy2ll.h"

void mexFunction( int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]){

        int i,verbose=1;

        /*input: */
        double  *x=NULL,*y=NULL;
        int     nx,ny,ncoord;
        int     sgn;

        Options* options=NULL;
        double   cm=0.,sp=0.;

        /* output: */
        double  *lat=NULL,*lon=NULL;
        int     iret=0;

        /*Boot module: */
        MODULEBOOT();

        /*checks on arguments on the matlab side: */
        if (nlhs > NLHS) {
                Xy2llUsage();
                _error_("Xy2ll usage error");
        }
        if (nrhs < NRHS) {
                Xy2llUsage();
                _error_("Xy2ll usage error");
        }

        /*Input datasets: */
        if (verbose) printf("Fetching inputs:\n");
        FetchData(&x,&nx,X_IN);
        if (verbose)
                if   (nx == 1) printf("  x=%g\n",x[0]);
                else           printf("  x=[%d values]\n",nx);
//      for (i=0; i<nx; i++) printf("  x[%d]=%g\n",i,x[i]);
        FetchData(&y,&ny,Y_IN);
        if (verbose)
                if   (ny == 1) printf("  y=%g\n",y[0]);
                else           printf("  y=[%d values]\n",ny);
//      for (i=0; i<ny; i++) printf("  y[%d]=%g\n",i,y[i]);
        FetchData(&sgn,SGN_IN);
        if (verbose) printf("  sgn=%d\n",sgn);

        if (verbose) printf("Parsing options:\n");
        options=new Options(NRHS,nrhs,prhs);
        if (options->Size()) for(i=0;i<options->Size();i++) ((Option*)options->GetObjectByOffset(i))->DeepEcho();
        /*  defaults are in Xy2lldef, so don't duplicate them here, and only use user values if both have been specified  */
        if (options->GetOption("central_meridian") || options->GetOption("standard_parallel")) {
                options->Get(&cm,"central_meridian");
                if (verbose) printf("  cm=%g\n",cm);
                options->Get(&sp,"standard_parallel");
                if (verbose) printf("  sp=%g\n",sp);
        }

        /*some checks*/
        if (verbose) printf("Checking inputs:\n");

        if   (nx != ny) _error_("Must have same number of x[%d] and y[%d] coordinates.",nx,ny);
        else            ncoord=nx;
        if (sgn != +1 && sgn != -1) _error_("Hemisphere sgn=%d must be +1 (north) or -1 (south).",sgn);
        if (fabs(cm)      > 180.) _error_("Central meridian cm=%g must be between -180 (west) and +180 (east) degrees.",cm);
        if (sp < 0. || sp >  90.) _error_("Standard parallel sp=%g must be between 0 and 90 degrees (in specified hemisphere).",sp);

        lat=(double *)xmalloc(ncoord*sizeof(double));
        lon=(double *)xmalloc(ncoord*sizeof(double));

        /* Run core computations: */
        if (verbose) printf("Calling core:\n");
        if (options->GetOption("central_meridian") && options->GetOption("standard_parallel"))
                iret=Xy2llx(lat,lon,
                                        x,y,ncoord,
                                        sgn,cm,sp);
        else
                iret=Xy2llx(lat,lon,
                                        x,y,ncoord,
                                        sgn);
        if (verbose) printf("  iret=%d\n",iret);

        /*Write data: */
        WriteData(LAT_OUT,lat,ncoord);
        WriteData(LON_OUT,lon,ncoord);

        /*Clean-up*/
        delete options;

        /*end module: */
        MODULEEND();
}

void Xy2llUsage(void)
{
        _printf_(true,"Xy2ll - x/y to lat/long coordinate transformation module:\n");
        _printf_(true,"\n");
        _printf_(true,"   This module transforms x/y to lat/long coordinates.\n");
        _printf_(true,"\n");
        _printf_(true,"   Usage:\n");
        _printf_(true,"      [lat,lon]=Xy2ll(x,y,sgn,'param name',param,...);\n");
        _printf_(true,"\n");
        _printf_(true,"      x           x coordinates (double vector)\n");
        _printf_(true,"      y           y coordinates (double vector)\n");
        _printf_(true,"      sgn         sign for hemisphere (double, +1 (north) or -1 (south))\n");
        _printf_(true,"\n");
        _printf_(true,"      central_meridian     central meridian (double, optional, but must specify with sp)\n");
        _printf_(true,"      standard_parallel    standard parallel (double, optional, but must specify with cm)\n");
        _printf_(true,"\n");
        _printf_(true,"      lat         latitude coordinates (double vector)\n");
        _printf_(true,"      lon         longitude coordinates (double vector)\n");
        _printf_(true,"\n");
        _printf_(true,"   Examples:\n");
        _printf_(true,"      [lat,lon]=Xy2ll(x,y, 1);\n");
        _printf_(true,"      [lat,lon]=Xy2ll(x,y, 1,'central_meridian',45,'standard_parallel',70);\n");
        _printf_(true,"      [lat,lon]=Xy2ll(x,y,-1,'central_meridian', 0,'standard_parallel',71);\n");
        _printf_(true,"\n");
}