static char mc_version[] = "MATLAB Compiler 1.2 Jan 17 1998 infun";
/*
 *  MATLAB Compiler: 1.2
 *  Date: Jan 17 1998
 *  Arguments: vision_engine_b 
 */
#ifndef ARRAY_ACCESS_INLINING
#error You must use the -inline option when compiling MATLAB compiler generated code with MEX or MBUILD
#endif
#ifndef MATLAB_COMPILER_GENERATED_CODE
#define MATLAB_COMPILER_GENERATED_CODE
#endif

#include <math.h>
#include "mex.h"
#include "mcc.h"

static void vision_engine_b_color_b(mxArray *, mxArray *);
/* static array S0_ (3 x 3) int, line 74 */
static int S0r_[] =
{
	2, 1, 8, 3, 0, 7, 4, 5, 
	6, 
};
static mxArray S0_ = mccCINIT( mccINT, 3, 3, S0r_, 0 );
/* static array S1_ (1 x 3) int, line 76 */
static int S1r_[] =
{
	-1, 0, 1, 
};
static mxArray S1_ = mccCINIT( mccINT, 1, 3, S1r_, 0 );
/* static array S2_ (8 x 2) int, line 78 */
static int S2r_[] =
{
	0, -1, -1, -1, 
	0, 1, 1, 1, -1, -1, 0, 1, 
	1, 1, 0, -1, 
};
static mxArray S2_ = mccCINIT( mccINT, 8, 2, S2r_, 0 );
/* static array S3_ (1 x 2) int, line 89 */
static int S3r_[] =
{
	1, 1, 
};
static mxArray S3_ = mccCINIT( mccINT, 1, 2, S3r_, 0 );
/* static array S4_ (1 x 5) text, line 12: 'close' */
static unsigned short S4__r_[] =
{
         99,  108,  111,  115,  101,
};
static mxArray S4_ = mccCINIT( mccTEXT,  1, 5, S4__r_, 0);
/* static array S5_ (1 x 4) text, line 13: 'fill' */
static unsigned short S5__r_[] =
{
        102,  105,  108,  108,
};
static mxArray S5_ = mccCINIT( mccTEXT,  1, 4, S5__r_, 0);
/* static array S6_ (1 x 4) text, line 14: 'spur' */
static unsigned short S6__r_[] =
{
        115,  112,  117,  114,
};
static mxArray S6_ = mccCINIT( mccTEXT,  1, 4, S6__r_, 0);
/* static array S7_ (1 x 5) text, line 15: 'clean' */
static unsigned short S7__r_[] =
{
         99,  108,  101,   97,  110,
};
static mxArray S7_ = mccCINIT( mccTEXT,  1, 5, S7__r_, 0);
/* static array S8_ (1 x 50) text, line 38: 'number of topological com...' */
static unsigned short S8__r_[] =
{
        110,  117,  109,   98,  101,  114,   32,  111,
        102,   32,  116,  111,  112,  111,  108,  111,
        103,  105,   99,   97,  108,   32,   99,  111,
        109,  112,  111,  110,  101,  110,  116,  115,
         32,  105,  110,   32,  111,   98,  115,  101,
        114,  118,  101,  100,   32,  105,  109,   97,
        103,  101,
};
static mxArray S8_ = mccCINIT( mccTEXT,  1, 50, S8__r_, 0);
/* static array S9_ (1 x 4) text, line 49: 'spur' */
static unsigned short S9__r_[] =
{
        115,  112,  117,  114,
};
static mxArray S9_ = mccCINIT( mccTEXT,  1, 4, S9__r_, 0);
/* static array S10_ (1 x 5) text, line 54: 'image' */
static unsigned short S10__r_[] =
{
        105,  109,   97,  103,  101,
};
static mxArray S10_ = mccCINIT( mccTEXT,  1, 5, S10__r_, 0);
/* static array S11_ (1 x 14) text, line 55: 'OBSERVED IMAGE' */
static unsigned short S11__r_[] =
{
         79,   66,   83,   69,   82,   86,   69,   68,
         32,   73,   77,   65,   71,   69,
};
static mxArray S11_ = mccCINIT( mccTEXT,  1, 14, S11__r_, 0);
/* static array S12_ (1 x 5) text, line 59: 'image' */
static unsigned short S12__r_[] =
{
        105,  109,   97,  103,  101,
};
static mxArray S12_ = mccCINIT( mccTEXT,  1, 5, S12__r_, 0);
/* static array S13_ (1 x 13) text, line 60: 'CURRENT FOVEA' */
static unsigned short S13__r_[] =
{
         67,   85,   82,   82,   69,   78,   84,   32,
         70,   79,   86,   69,   65,
};
static mxArray S13_ = mccCINIT( mccTEXT,  1, 13, S13__r_, 0);
/* static array S14_ (1 x 22) text, line 64: 'Please wait a while...' */
static unsigned short S14__r_[] =
{
         80,  108,  101,   97,  115,  101,   32,  119,
         97,  105,  116,   32,   97,   32,  119,  104,
        105,  108,  101,   46,   46,   46,
};
static mxArray S14_ = mccCINIT( mccTEXT,  1, 22, S14__r_, 0);
/* static array S15_ (1 x 23) text, line 117: 'chaincoding boundary...' */
static unsigned short S15__r_[] =
{
         99,  104,   97,  105,  110,   99,  111,  100,
        105,  110,  103,   32,   98,  111,  117,  110,
        100,   97,  114,  121,   46,   46,   46,
};
static mxArray S15_ = mccCINIT( mccTEXT,  1, 23, S15__r_, 0);
/* static array S16_ (1 x 22) text, line 166: 'computing curvature...' */
static unsigned short S16__r_[] =
{
         99,  111,  109,  112,  117,  116,  105,  110,
        103,   32,   99,  117,  114,  118,   97,  116,
        117,  114,  101,   46,   46,   46,
};
static mxArray S16_ = mccCINIT( mccTEXT,  1, 22, S16__r_, 0);
/* static array S17_ (1 x 33) text, line 206: 'dynamic programming, back...' */
static unsigned short S17__r_[] =
{
        100,  121,  110,   97,  109,  105,   99,   32,
        112,  114,  111,  103,  114,   97,  109,  109,
        105,  110,  103,   44,   32,   98,   97,   99,
        107,  119,   97,  114,  100,  115,   46,   46,
         46,
};
static mxArray S17_ = mccCINIT( mccTEXT,  1, 33, S17__r_, 0);
/* static array S18_ (1 x 148) text, line 228: 'Function zeros being call...' */
static unsigned short S18__r_[] =
{
         70,  117,  110,   99,  116,  105,  111,  110,
         32,  122,  101,  114,  111,  115,   32,   98,
        101,  105,  110,  103,   32,   99,   97,  108,
        108,  101,  100,   32,  105,  110,   99,  111,
        114,  114,  101,   99,  116,  108,  121,   46,
         32,   32,   69,  105,  116,  104,  101,  114,
         32,  116,  111,  111,   32,  109,   97,  110,
        121,   32,  105,  110,  112,  117,  116,   32,
        111,  114,   32,  111,  117,  116,  112,  117,
        116,   32,   97,  114,  103,  117,  109,  101,
        110,  116,  115,   32,   97,  114,  101,   32,
         98,  101,  105,  110,  103,   32,  112,   97,
        115,  115,  101,  100,   32,  111,  114,   32,
        116,  104,  101,   32,   97,  114,  103,  117,
        109,  101,  110,  116,  115,   32,   97,  114,
        101,   32,  105,  110,  118,   97,  108,  105,
        100,   32,  102,  111,  114,   32,  116,  104,
        105,  115,   32,  102,  117,  110,   99,  116,
        105,  111,  110,   46,
};
static mxArray S18_ = mccCINIT( mccTEXT,  1, 148, S18__r_, 0);
/* static array S19_ (1 x 28) text, line 267: 'dynamic programming, forw...' */
static unsigned short S19__r_[] =
{
        100,  121,  110,   97,  109,  105,   99,   32,
        112,  114,  111,  103,  114,   97,  109,  109,
        105,  110,  103,   44,   32,  102,  111,  114,
        119,   97,  114,  100,
};
static mxArray S19_ = mccCINIT( mccTEXT,  1, 28, S19__r_, 0);
/* static array S20_ (1 x 15) text, line 295: 'Recognized arcs' */
static unsigned short S20__r_[] =
{
         82,  101,   99,  111,  103,  110,  105,  122,
        101,  100,   32,   97,  114,   99,  115,
};
static mxArray S20_ = mccCINIT( mccTEXT,  1, 15, S20__r_, 0);
/* static array S21_ (1 x 5) text, line 296: 'image' */
static unsigned short S21__r_[] =
{
        105,  109,   97,  103,  101,
};
static mxArray S21_ = mccCINIT( mccTEXT,  1, 5, S21__r_, 0);
/* static array S22_ (1 x 5) text, line 302: 'image' */
static unsigned short S22__r_[] =
{
        105,  109,   97,  103,  101,
};
static mxArray S22_ = mccCINIT( mccTEXT,  1, 5, S22__r_, 0);
/* static array S23_ (1 x 29) text, line 303: 'Estimated deformation of ...' */
static unsigned short S23__r_[] =
{
         69,  115,  116,  105,  109,   97,  116,  101,
        100,   32,  100,  101,  102,  111,  114,  109,
         97,  116,  105,  111,  110,   32,  111,  102,
         32,   97,  114,   99,  115,
};
static mxArray S23_ = mccCINIT( mccTEXT,  1, 29, S23__r_, 0);
/* static array S24_ (3 x 3) int, line 74 */
static int S24r_[] =
{
	2, 1, 
	8, 3, 0, 7, 4, 5, 6, 
};
static mxArray S24_ = mccCINIT( mccINT, 3, 3, S24r_, 0 );
/* static array S25_ (1 x 3) int, line 76 */
static int S25r_[] =
{
	-1, 
	0, 1, 
};
static mxArray S25_ = mccCINIT( mccINT, 1, 3, S25r_, 0 );
/* static array S26_ (8 x 2) int, line 78 */
static int S26r_[] =
{
	0, -1, -1, -1, 0, 1, 
	1, 1, -1, -1, 0, 1, 1, 1, 
	0, -1, 
};
static mxArray S26_ = mccCINIT( mccINT, 8, 2, S26r_, 0 );
/* static array S27_ (1 x 2) int, line 89 */
static int S27r_[] =
{
	1, 1, 
};
static mxArray S27_ = mccCINIT( mccINT, 1, 2, S27r_, 0 );
/* static array S28_ (1 x 5) text, line 12: 'close' */
static unsigned short S28__r_[] =
{
         99,  108,  111,  115,  101,
};
static mxArray S28_ = mccCINIT( mccTEXT,  1, 5, S28__r_, 0);
/* static array S29_ (1 x 4) text, line 13: 'fill' */
static unsigned short S29__r_[] =
{
        102,  105,  108,  108,
};
static mxArray S29_ = mccCINIT( mccTEXT,  1, 4, S29__r_, 0);
/* static array S30_ (1 x 4) text, line 14: 'spur' */
static unsigned short S30__r_[] =
{
        115,  112,  117,  114,
};
static mxArray S30_ = mccCINIT( mccTEXT,  1, 4, S30__r_, 0);
/* static array S31_ (1 x 5) text, line 15: 'clean' */
static unsigned short S31__r_[] =
{
         99,  108,  101,   97,  110,
};
static mxArray S31_ = mccCINIT( mccTEXT,  1, 5, S31__r_, 0);
/* static array S32_ (1 x 50) text, line 38: 'number of topological com...' */
static unsigned short S32__r_[] =
{
        110,  117,  109,   98,  101,  114,   32,  111,
        102,   32,  116,  111,  112,  111,  108,  111,
        103,  105,   99,   97,  108,   32,   99,  111,
        109,  112,  111,  110,  101,  110,  116,  115,
         32,  105,  110,   32,  111,   98,  115,  101,
        114,  118,  101,  100,   32,  105,  109,   97,
        103,  101,
};
static mxArray S32_ = mccCINIT( mccTEXT,  1, 50, S32__r_, 0);
/* static array S33_ (1 x 4) text, line 49: 'spur' */
static unsigned short S33__r_[] =
{
        115,  112,  117,  114,
};
static mxArray S33_ = mccCINIT( mccTEXT,  1, 4, S33__r_, 0);
/* static array S34_ (1 x 5) text, line 54: 'image' */
static unsigned short S34__r_[] =
{
        105,  109,   97,  103,  101,
};
static mxArray S34_ = mccCINIT( mccTEXT,  1, 5, S34__r_, 0);
/* static array S35_ (1 x 14) text, line 55: 'OBSERVED IMAGE' */
static unsigned short S35__r_[] =
{
         79,   66,   83,   69,   82,   86,   69,   68,
         32,   73,   77,   65,   71,   69,
};
static mxArray S35_ = mccCINIT( mccTEXT,  1, 14, S35__r_, 0);
/* static array S36_ (1 x 5) text, line 59: 'image' */
static unsigned short S36__r_[] =
{
        105,  109,   97,  103,  101,
};
static mxArray S36_ = mccCINIT( mccTEXT,  1, 5, S36__r_, 0);
/* static array S37_ (1 x 13) text, line 60: 'CURRENT FOVEA' */
static unsigned short S37__r_[] =
{
         67,   85,   82,   82,   69,   78,   84,   32,
         70,   79,   86,   69,   65,
};
static mxArray S37_ = mccCINIT( mccTEXT,  1, 13, S37__r_, 0);
/* static array S38_ (1 x 22) text, line 64: 'Please wait a while...' */
static unsigned short S38__r_[] =
{
         80,  108,  101,   97,  115,  101,   32,  119,
         97,  105,  116,   32,   97,   32,  119,  104,
        105,  108,  101,   46,   46,   46,
};
static mxArray S38_ = mccCINIT( mccTEXT,  1, 22, S38__r_, 0);
/* static array S39_ (1 x 23) text, line 117: 'chaincoding boundary...' */
static unsigned short S39__r_[] =
{
         99,  104,   97,  105,  110,   99,  111,  100,
        105,  110,  103,   32,   98,  111,  117,  110,
        100,   97,  114,  121,   46,   46,   46,
};
static mxArray S39_ = mccCINIT( mccTEXT,  1, 23, S39__r_, 0);
/* static array S40_ (1 x 22) text, line 166: 'computing curvature...' */
static unsigned short S40__r_[] =
{
         99,  111,  109,  112,  117,  116,  105,  110,
        103,   32,   99,  117,  114,  118,   97,  116,
        117,  114,  101,   46,   46,   46,
};
static mxArray S40_ = mccCINIT( mccTEXT,  1, 22, S40__r_, 0);
/* static array S41_ (1 x 33) text, line 206: 'dynamic programming, back...' */
static unsigned short S41__r_[] =
{
        100,  121,  110,   97,  109,  105,   99,   32,
        112,  114,  111,  103,  114,   97,  109,  109,
        105,  110,  103,   44,   32,   98,   97,   99,
        107,  119,   97,  114,  100,  115,   46,   46,
         46,
};
static mxArray S41_ = mccCINIT( mccTEXT,  1, 33, S41__r_, 0);
/* static array S42_ (1 x 148) text, line 228: 'Function zeros being call...' */
static unsigned short S42__r_[] =
{
         70,  117,  110,   99,  116,  105,  111,  110,
         32,  122,  101,  114,  111,  115,   32,   98,
        101,  105,  110,  103,   32,   99,   97,  108,
        108,  101,  100,   32,  105,  110,   99,  111,
        114,  114,  101,   99,  116,  108,  121,   46,
         32,   32,   69,  105,  116,  104,  101,  114,
         32,  116,  111,  111,   32,  109,   97,  110,
        121,   32,  105,  110,  112,  117,  116,   32,
        111,  114,   32,  111,  117,  116,  112,  117,
        116,   32,   97,  114,  103,  117,  109,  101,
        110,  116,  115,   32,   97,  114,  101,   32,
         98,  101,  105,  110,  103,   32,  112,   97,
        115,  115,  101,  100,   32,  111,  114,   32,
        116,  104,  101,   32,   97,  114,  103,  117,
        109,  101,  110,  116,  115,   32,   97,  114,
        101,   32,  105,  110,  118,   97,  108,  105,
        100,   32,  102,  111,  114,   32,  116,  104,
        105,  115,   32,  102,  117,  110,   99,  116,
        105,  111,  110,   46,
};
static mxArray S42_ = mccCINIT( mccTEXT,  1, 148, S42__r_, 0);
/* static array S43_ (1 x 28) text, line 267: 'dynamic programming, forw...' */
static unsigned short S43__r_[] =
{
        100,  121,  110,   97,  109,  105,   99,   32,
        112,  114,  111,  103,  114,   97,  109,  109,
        105,  110,  103,   44,   32,  102,  111,  114,
        119,   97,  114,  100,
};
static mxArray S43_ = mccCINIT( mccTEXT,  1, 28, S43__r_, 0);
/* static array S44_ (1 x 15) text, line 295: 'Recognized arcs' */
static unsigned short S44__r_[] =
{
         82,  101,   99,  111,  103,  110,  105,  122,
        101,  100,   32,   97,  114,   99,  115,
};
static mxArray S44_ = mccCINIT( mccTEXT,  1, 15, S44__r_, 0);
/* static array S45_ (1 x 5) text, line 296: 'image' */
static unsigned short S45__r_[] =
{
        105,  109,   97,  103,  101,
};
static mxArray S45_ = mccCINIT( mccTEXT,  1, 5, S45__r_, 0);
/* static array S46_ (1 x 5) text, line 302: 'image' */
static unsigned short S46__r_[] =
{
        105,  109,   97,  103,  101,
};
static mxArray S46_ = mccCINIT( mccTEXT,  1, 5, S46__r_, 0);
/* static array S47_ (1 x 29) text, line 303: 'Estimated deformation of ...' */
static unsigned short S47__r_[] =
{
         69,  115,  116,  105,  109,   97,  116,  101,
        100,   32,  100,  101,  102,  111,  114,  109,
         97,  116,  105,  111,  110,   32,  111,  102,
         32,   97,  114,   99,  115,
};
static mxArray S47_ = mccCINIT( mccTEXT,  1, 29, S47__r_, 0);
/* static array S48_ (1 x 2) text, line 392: 'on' */
static unsigned short S48__r_[] =
{
        111,  110,
};
static mxArray S48_ = mccCINIT( mccTEXT,  1, 2, S48__r_, 0);
/* static array S49_ (1 x 3) text, line 394: 'off' */
static unsigned short S49__r_[] =
{
        111,  102,  102,
};
static mxArray S49_ = mccCINIT( mccTEXT,  1, 3, S49__r_, 0);
/***************** Compiler Assumptions ****************
 *
 *       B           	complex vector/matrix
 *       BC          	complex vector/matrix
 *       BM0_        	boolean vector/matrix temporary
 *       BM1_        	boolean vector/matrix temporary
 *       BM2_        	boolean vector/matrix temporary
 *       BM3_        	boolean vector/matrix temporary
 *       Bframe      	complex vector/matrix
 *       C0_         	complex scalar temporary
 *       C1_         	complex scalar temporary
 *       C2_         	complex scalar temporary
 *       C3_         	complex scalar temporary
 *       CM0_        	complex vector/matrix temporary
 *       CM1_        	complex vector/matrix temporary
 *       CM2_        	complex vector/matrix temporary
 *       CM3_        	complex vector/matrix temporary
 *       CM4_        	complex vector/matrix temporary
 *       CM5_        	complex vector/matrix temporary
 *       CM6_        	complex vector/matrix temporary
 *       I           	boolean vector/matrix
 *       I0_         	integer scalar temporary
 *       I1_         	integer scalar temporary
 *       I2_         	integer scalar temporary
 *       I3_         	integer scalar temporary
 *       ID          	real vector/matrix
 *       ID          	complex vector/matrix  => ID_1
 *       IM0_        	integer vector/matrix temporary
 *       IM1_        	integer vector/matrix temporary
 *       IM2_        	integer vector/matrix temporary
 *       IM3_        	integer vector/matrix temporary
 *       I           	complex vector/matrix  => I_1
 *       Itemps      	complex vector/matrix
 *       M           	integer vector/matrix
 *       M           	complex vector/matrix  => M_1
 *       Neigh       	complex vector/matrix
 *       R0_         	real scalar temporary
 *       R1_         	real scalar temporary
 *       RM0_        	real vector/matrix temporary
 *       RM1_        	real vector/matrix temporary
 *       RM2_        	real vector/matrix temporary
 *       RM3_        	real vector/matrix temporary
 *       RM4_        	real vector/matrix temporary
 *       RM5_        	real vector/matrix temporary
 *       S0_         	<constant>
 *       S1_         	<constant>
 *       S2_         	<constant>
 *       S3_         	<constant>
 *       U           	complex vector/matrix
 *       V           	complex vector/matrix
 *       W           	integer vector/matrix
 *       Y           	integer vector/matrix
 *       Z           	complex vector/matrix
 *       a           	real scalar
 *       a           	real vector/matrix  => a_1
 *       abs         	<function>
 *       alpha       	integer scalar
 *       alphas      	integer vector/matrix
 *       analysis    	integer vector/matrix
 *       any         	<function>
 *       asin        	<function>
 *       attention   	integer vector/matrix
 *       axis        	<function>
 *       b           	real vector/matrix
 *       bwlabel     	<function>
 *       bwmorph     	<function>
 *       bwperim     	<function>
 *       c1          	complex vector/matrix
 *       chaincode   	complex vector/matrix
 *       close       	<function>
 *       code        	complex vector/matrix
 *       colindex    	integer vector/matrix
 *       compnumber  	integer scalar
 *       cost        	complex vector/matrix
 *       cost1       	complex vector/matrix
 *       cuml        	real vector/matrix
 *       cumsum      	<function>
 *       d           	integer scalar
 *       diffs       	real vector/matrix
 *       dl          	integer scalar
 *       double      	<function>
 *       dphis       	complex vector/matrix
 *       elastc      	real scalar
 *       error       	<function>
 *       fft         	<function>
 *       floor       	<function>
 *       fovea       	boolean vector/matrix
 *       i           	integer scalar
 *       i           	integer vector/matrix  => i_1
 *       ifft        	<function>
 *       image       	<function>
 *       incr        	integer vector/matrix
 *       jumpc       	real scalar
 *       kron        	<function>
 *       l           	integer scalar
 *       l1          	integer scalar
 *       l2          	integer scalar
 *       labels      	complex vector/matrix
 *       length      	<function>
 *       lengths     	real vector/matrix
 *       m           	real scalar
 *       m           	integer vector/matrix  => m_1
 *       max         	<function>
 *       maxatt      	integer vector/matrix
 *       min         	<function>
 *       mod         	<function>
 *       n           	integer scalar
 *       n           	integer vector/matrix  => n_1
 *       nalpha      	integer scalar
 *       ncomp       	integer scalar
 *       next        	integer vector/matrix
 *       norm        	<function>
 *       off         	integer vector/matrix
 *       ones        	<function>
 *       p           	real scalar
 *       pause       	<function>
 *       pi          	<function>
 *       plot        	<function>
 *       r           	complex scalar
 *       real        	<function>
 *       reshape     	<function>
 *       rhos        	complex vector/matrix
 *       rowindex    	integer vector/matrix
 *       size        	<function>
 *       spline      	<function>
 *       sqrt        	<function>
 *       start       	integer vector/matrix
 *       sum         	<function>
 *       t           	integer scalar
 *       t1          	integer scalar
 *       templs      	complex vector/matrix
 *       title       	<function>
 *       tmax        	integer scalar
 *       ts          	integer vector/matrix
 *       u           	complex vector/matrix
 *       u           	integer vector/matrix  => u_1
 *       uhat        	complex vector/matrix
 *       v           	boolean vector/matrix
 *       v           	real vector/matrix  => v_1
 *       v           	integer vector/matrix  => v_2
 *       vision_engine_b	<function being defined>
 *       vision_engine_b/color_b	<function>
 *       w           	complex vector/matrix
 *       x           	real scalar
 *       x1shat      	complex vector/matrix
 *       x2shat      	complex vector/matrix
 *       x           	complex vector/matrix  => x_1
 *       xs          	real vector/matrix
 *       y           	real scalar
 *       y1          	real vector/matrix
 *       y2          	real vector/matrix
 *       y           	complex scalar  => y_1
 *       y           	complex vector/matrix  => y_2
 *       ys          	real vector/matrix
 *       zeros       	<function>
 *******************************************************/

void
mexFunction(
    int nlhs_,
    mxArray *plhs_[],
    int nrhs_,
    const mxArray *prhs_[]
)
{
   int ci_;
   int i_;
   mxArray *Mplhs_[3];
   mxArray *Mprhs_[17];
   

   if (nrhs_ > 3 )
   {
      mexErrMsgTxt( "Too many input arguments." );
   }

   if (nlhs_ > 6 )
   {
      mexErrMsgTxt( "Too many output arguments." );
   }

   for (ci_=i_=0; i_<nrhs_; ++i_)
   {
      if ( mccPI(prhs_[i_]) )
      {
         ci_ = 1;
         break;
      }
   }
   if (ci_)
   {
      mxArray B;
      mxArray BC;
      mxArray rhos;
      mxArray analysis;
      mxArray xs;
      mxArray ys;
      mxArray ID;
      mxArray Itemps;
      mxArray templs;
      mxArray I;
      mxArray I_1;
      mxArray labels;
      int ncomp = 0;
      mxArray attention;
      int i = 0;
      int compnumber = 0;
      mxArray maxatt;
      mxArray i_1;
      mxArray fovea;
      mxArray v;
      mxArray Y;
      mxArray colindex;
      mxArray Z;
      mxArray rowindex;
      mxArray start;
      mxArray W;
      mxArray off;
      mxArray incr;
      int l1 = 0;
      int l2 = 0;
      mxArray Bframe;
      mxArray Neigh;
      mxArray chaincode;
      mxArray code;
      mxArray next;
      int t = 0;
      int l = 0;
      int l_flags_ = 0;
      int d = 0;
      double m = 0.0;
      double a = 0.0;
      double p = 0.0;
      mxArray v_1;
      mxArray w;
      mxArray u;
      mxArray uhat;
      mxArray x1shat;
      mxArray x2shat;
      mxArray lengths;
      mxArray cuml;
      mxArray dphis;
      mxArray a_1;
      mxArray b;
      double x = 0.0;
      double y = 0.0;
      double r_r=0.0, r_i=0.0;
      int dl = 0;
      int n = 0;
      mxArray M_y1;
      mxArray y2;
      double jumpc = 0.0;
      double elastc = 0.0;
      int nalpha = 0;
      int nalpha_set_ = 0;
      int tmax = 0;
      int tmax_set_ = 0;
      mxArray alphas;
      mxArray ts;
      mxArray M;
      mxArray ID_1;
      mxArray u_1;
      mxArray v_2;
      int alpha = 0;
      mxArray diffs;
      mxArray M_1;
      mxArray x_1;
      mxArray m_1;
      double y_1_r=0.0, y_1_i=0.0;
      mxArray n_1;
      mxArray cost;
      int t1 = 0;
      mxArray cost1;
      mxArray c1;
      mxArray U;
      mxArray V;
      mxArray y_2;
      mxArray CM0_;
      mxArray CM1_;
      int I0_ = 0;
      mxArray BM0_;
      mxArray IM0_;
      mxArray IM1_;
      int I1_ = 0;
      mxArray CM2_;
      mxArray IM2_;
      mxArray IM3_;
      mxArray RM0_;
      mxArray RM1_;
      mxArray CM3_;
      double R0_ = 0.0;
      double C0__r=0.0, C0__i=0.0;
      double C1__r=0.0, C1__i=0.0;
      double C2__r=0.0, C2__i=0.0;
      mxArray CM4_;
      double C3__r=0.0, C3__i=0.0;
      mxArray RM2_;
      mxArray RM3_;
      mxArray RM4_;
      mxArray RM5_;
      double R1_ = 0.0;
      int I2_ = 0;
      mxArray BM1_;
      mxArray BM2_;
      mxArray BM3_;
      mxArray CM5_;
      mxArray CM6_;
      int I3_ = 0;
      
      mccRealInit(ID);
      mccImport(&ID, ((nrhs_>0) ? prhs_[0] : 0), 0, 0);
      mccComplexInit(Itemps);
      mccImport(&Itemps, ((nrhs_>1) ? prhs_[1] : 0), 0, 0);
      mccComplexInit(templs);
      mccImport(&templs, ((nrhs_>2) ? prhs_[2] : 0), 0, 0);
      mccComplexInit(B);
      mccComplexInit(BC);
      mccComplexInit(rhos);
      mccIntInit(analysis);
      mccRealInit(xs);
      mccRealInit(ys);
      mccBoolInit(I);
      mccComplexInit(I_1);
      mccComplexInit(labels);
      mccIntInit(attention);
      mccIntInit(maxatt);
      mccIntInit(i_1);
      mccBoolInit(fovea);
      mccBoolInit(v);
      mccIntInit(Y);
      mccIntInit(colindex);
      mccComplexInit(Z);
      mccIntInit(rowindex);
      mccIntInit(start);
      mccIntInit(W);
      mccIntInit(off);
      mccIntInit(incr);
      mccComplexInit(Bframe);
      mccComplexInit(Neigh);
      mccComplexInit(chaincode);
      mccComplexInit(code);
      mccIntInit(next);
      mccRealInit(v_1);
      mccComplexInit(w);
      mccComplexInit(u);
      mccComplexInit(uhat);
      mccComplexInit(x1shat);
      mccComplexInit(x2shat);
      mccRealInit(lengths);
      mccRealInit(cuml);
      mccComplexInit(dphis);
      mccRealInit(a_1);
      mccRealInit(b);
      mccRealInit(M_y1);
      mccRealInit(y2);
      mccIntInit(alphas);
      mccIntInit(ts);
      mccIntInit(M);
      mccComplexInit(ID_1);
      mccIntInit(u_1);
      mccIntInit(v_2);
      mccRealInit(diffs);
      mccComplexInit(M_1);
      mccComplexInit(x_1);
      mccIntInit(m_1);
      mccIntInit(n_1);
      mccComplexInit(cost);
      mccComplexInit(cost1);
      mccComplexInit(c1);
      mccComplexInit(U);
      mccComplexInit(V);
      mccComplexInit(y_2);
      mccComplexInit(CM0_);
      mccComplexInit(CM1_);
      mccBoolInit(BM0_);
      mccIntInit(IM0_);
      mccIntInit(IM1_);
      mccComplexInit(CM2_);
      mccIntInit(IM2_);
      mccIntInit(IM3_);
      mccRealInit(RM0_);
      mccRealInit(RM1_);
      mccComplexInit(CM3_);
      mccComplexInit(CM4_);
      mccRealInit(RM2_);
      mccRealInit(RM3_);
      mccRealInit(RM4_);
      mccRealInit(RM5_);
      mccBoolInit(BM1_);
      mccBoolInit(BM2_);
      mccBoolInit(BM3_);
      mccComplexInit(CM5_);
      mccComplexInit(CM6_);
      
      /* % executes inference engine in Strategies for Seeing, Experiment 2 */
      /* % where generators are arcs */
      /* % ID=oberved image */
      /* % Itemps =matrix of templates extended periodically to length tmax */
      /* % templs =vector of template lengths */
      
      /* % 1.ENTER. */
      /* % make image binary by thresholding */
      /* I=ID>200; */
      if(mccNOTSET(&ID))
      {
         mexErrMsgTxt( "variable ID undefined, line 10" );
      }
      {
         int i_, j_;
         int m_=1, n_=1, cx_ = 0;
         unsigned short *p_I;
         int I_I=1;
         double *p_ID;
         int I_ID=1;
         m_ = mcmCalcResultSize(m_, &n_, mccM(&ID), mccN(&ID));
         mccAllocateMatrix(&I, m_, n_);
         I_I = (mccM(&I) != 1 || mccN(&I) != 1);
         p_I = mccSPR(&I);
         I_ID = (mccM(&ID) != 1 || mccN(&ID) != 1);
         p_ID = mccPR(&ID);
         if (m_ != 0)
         {
            for (j_=0; j_<n_; ++j_)
            {
               for (i_=0; i_<m_; ++i_, p_I+=I_I, p_ID+=I_ID)
               {
                  *p_I = ( (*p_ID > 200) && !mccREL_NAN(*p_ID) );
               }
            }
         }
      }
      /* % regularize image */
      /* I=bwmorph(I,'close'); */
      Mprhs_[0] = &I;
      Mprhs_[1] = &S4_;
      Mplhs_[0] = &I_1;
      mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "bwmorph", 12);
      /* I=bwmorph(I,'fill'); */
      Mprhs_[0] = &I_1;
      Mprhs_[1] = &S5_;
      Mplhs_[0] = &I_1;
      mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "bwmorph", 13);
      /* I=bwmorph(I,'spur'); */
      Mprhs_[0] = &I_1;
      Mprhs_[1] = &S6_;
      Mplhs_[0] = &I_1;
      mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "bwmorph", 14);
      /* I=bwmorph(I,'clean'); */
      Mprhs_[0] = &I_1;
      Mprhs_[1] = &S7_;
      Mplhs_[0] = &I_1;
      mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "bwmorph", 15);
      /* I=double(I); */
      Mprhs_[0] = &I_1;
      Mplhs_[0] = &I_1;
      mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "double", 16);
      
      /* % opening frames */
      /* % page1 */
      /* % pause(1) */
      /* % close */
      /* % page2 */
      /* % pause(1) */
      /* % close */
      /* % page3b */
      /* % pause(1) */
      /* % close */
      
      
      /* % 2. SACCADIC SEARCH FOR FOVEAS */
      /* % find the ncomp topological components and their sizes */
      /* labels=bwlabel(I); */
      Mprhs_[0] = &I_1;
      Mplhs_[0] = &labels;
      mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "bwlabel", 32);
      /* ncomp=max(max(labels)); */
      Mprhs_[0] = &labels;
      Mplhs_[0] = &CM0_;
      mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 33);
      Mprhs_[0] = &CM0_;
      Mplhs_[0] = &CM1_;
      mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 33);
      ncomp = (mccGetRealVectorElement(&CM1_, mccRint(1)));
      /* attention=zeros(1,ncomp); */
      mccZerosMN(&attention, 1, ncomp);
      /* for i=1:ncomp */
      for (I0_ = 1; I0_ <= ncomp; I0_ = I0_ + 1)
      {
         i = I0_;
         /* attention(i)=sum(sum(labels==i)); */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            unsigned short *p_BM0_;
            int I_BM0_=1;
            double *p_labels;
            int I_labels=1;
            double *q_labels;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&labels), mccN(&labels));
            mccAllocateMatrix(&BM0_, m_, n_);
            I_BM0_ = (mccM(&BM0_) != 1 || mccN(&BM0_) != 1);
            p_BM0_ = mccSPR(&BM0_);
            I_labels = (mccM(&labels) != 1 || mccN(&labels) != 1);
            p_labels = mccPR(&labels);
            q_labels = mccPI(&labels);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_BM0_+=I_BM0_, p_labels+=I_labels, q_labels+=I_labels)
                  {
                     *p_BM0_ = ( ( (*p_labels == i) && !mccREL_NAN(*p_labels) ) && (*q_labels == 0.));
                  }
               }
            }
         }
         mccSum(&IM0_, &BM0_);
         mccSum(&IM1_, &IM0_);
         I1_ = ((int)mccGetRealVectorElement(&IM1_, mccRint(1)));
         mccSetIntVectorElement(&attention, mccRint(i), (double)I1_);
         /* end  */
      }
      /* 'number of topological components in observed image' */
      mccPrint (&S8_, 0);
      /* ncomp */
      mccPrint (mccTempMatrix(ncomp, 0., mccINT, 0 ), "ncomp");
      /* attention */
      mccPrint (&attention, "attention");
      
      /* compnumber=1; */
      compnumber = 1;
      /* [maxatt,i]=max(attention); */
      Mprhs_[0] = &attention;
      Mplhs_[0] = &maxatt;
      Mplhs_[1] = &i_1;
      mccCallMATLAB(2, Mplhs_, 1, Mprhs_, "max", 43);
      /* while maxatt>100 */
      while (1)
      {
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            unsigned short *p_BM0_;
            int I_BM0_=1;
            int *p_maxatt;
            int I_maxatt=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&maxatt), mccN(&maxatt));
            mccAllocateMatrix(&BM0_, m_, n_);
            I_BM0_ = (mccM(&BM0_) != 1 || mccN(&BM0_) != 1);
            p_BM0_ = mccSPR(&BM0_);
            I_maxatt = (mccM(&maxatt) != 1 || mccN(&maxatt) != 1);
            p_maxatt = mccIPR(&maxatt);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_BM0_+=I_BM0_, p_maxatt+=I_maxatt)
                  {
                     *p_BM0_ = (((int)*p_maxatt) > 100);
                  }
               }
            }
         }
         if (!(mccIfCondition(&BM0_)))
            break;
         /* fovea=labels==i; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            unsigned short *p_fovea;
            int I_fovea=1;
            double *p_labels;
            int I_labels=1;
            double *q_labels;
            int *p_i_1;
            int I_i_1=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&labels), mccN(&labels));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&i_1), mccN(&i_1));
            mccAllocateMatrix(&fovea, m_, n_);
            I_fovea = (mccM(&fovea) != 1 || mccN(&fovea) != 1);
            p_fovea = mccSPR(&fovea);
            I_labels = (mccM(&labels) != 1 || mccN(&labels) != 1);
            p_labels = mccPR(&labels);
            q_labels = mccPI(&labels);
            I_i_1 = (mccM(&i_1) != 1 || mccN(&i_1) != 1);
            p_i_1 = mccIPR(&i_1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_fovea+=I_fovea, p_labels+=I_labels, q_labels+=I_labels, p_i_1+=I_i_1)
                  {
                     *p_fovea = ( ( (*p_labels == ((int)*p_i_1)) && !mccREL_NAN(*p_labels) ) && (*q_labels == 0.));
                  }
               }
            }
         }
         /* %find boundary B for fovea and chaincode BC */
         /* B=bwperim(fovea); */
         Mprhs_[0] = &fovea;
         Mplhs_[0] = &B;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "bwperim", 47);
         /* %regularize boundary */
         /* B=bwmorph(B,'spur'); */
         Mprhs_[0] = &B;
         Mprhs_[1] = &S9_;
         Mplhs_[0] = &B;
         mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "bwmorph", 49);
         /* B=double(B); */
         Mprhs_[0] = &B;
         Mplhs_[0] = &B;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "double", 50);
         
         /* %display observed image and current fovea */
         /* image(30.*I) */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_CM2_;
            int I_CM2_=1;
            double *q_CM2_;
            double *p_I_1;
            int I_I_1=1;
            double *q_I_1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&I_1), mccN(&I_1));
            mccAllocateMatrix(&CM2_, m_, n_);
            I_CM2_ = (mccM(&CM2_) != 1 || mccN(&CM2_) != 1);
            p_CM2_ = mccPR(&CM2_);
            q_CM2_ = mccPI(&CM2_);
            I_I_1 = (mccM(&I_1) != 1 || mccN(&I_1) != 1);
            p_I_1 = mccPR(&I_1);
            q_I_1 = mccPI(&I_1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_CM2_+=I_CM2_, q_CM2_+=I_CM2_, p_I_1+=I_I_1, q_I_1+=I_I_1)
                  {
                     {
                        double t_ = (30 * *p_I_1);
                        *q_CM2_ = (30 * *q_I_1);
                        *p_CM2_ = t_;
                     }
                  }
               }
            }
         }
         mccSTRING(&CM2_) = 0;
         Mprhs_[0] = &CM2_;
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 1, Mprhs_, "image", 53);
         mccPrint(Mplhs_[0], 0);
         /* axis('image') */
         Mprhs_[0] = &S10_;
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 1, Mprhs_, "axis", 54);
         mccPrint(Mplhs_[0], 0);
         /* title('OBSERVED IMAGE') */
         Mprhs_[0] = &S11_;
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 1, Mprhs_, "title", 55);
         /* pause */
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 0, Mprhs_, "pause", 56);
         /* close */
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 0, Mprhs_, "close", 57);
         /* image(30.*B) */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_CM2_;
            int I_CM2_=1;
            double *q_CM2_;
            double *p_B;
            int I_B=1;
            double *q_B;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&B), mccN(&B));
            mccAllocateMatrix(&CM2_, m_, n_);
            I_CM2_ = (mccM(&CM2_) != 1 || mccN(&CM2_) != 1);
            p_CM2_ = mccPR(&CM2_);
            q_CM2_ = mccPI(&CM2_);
            I_B = (mccM(&B) != 1 || mccN(&B) != 1);
            p_B = mccPR(&B);
            q_B = mccPI(&B);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_CM2_+=I_CM2_, q_CM2_+=I_CM2_, p_B+=I_B, q_B+=I_B)
                  {
                     {
                        double t_ = (30 * *p_B);
                        *q_CM2_ = (30 * *q_B);
                        *p_CM2_ = t_;
                     }
                  }
               }
            }
         }
         mccSTRING(&CM2_) = 0;
         Mprhs_[0] = &CM2_;
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 1, Mprhs_, "image", 58);
         mccPrint(Mplhs_[0], 0);
         /* axis('image') */
         Mprhs_[0] = &S12_;
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 1, Mprhs_, "axis", 59);
         mccPrint(Mplhs_[0], 0);
         /* title('CURRENT FOVEA') */
         Mprhs_[0] = &S13_;
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 1, Mprhs_, "title", 60);
         /* pause */
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 0, Mprhs_, "pause", 61);
         /* close */
         Mplhs_[0] = 0;
         mccCallMATLAB(0, Mplhs_, 0, Mprhs_, "close", 62);
         
         /* 'Please wait a while...' */
         mccPrint (&S14_, 0);
         
         
         /* %3. PREPARE  */
         /* %find starting pixel */
         /* v=any(B); */
         mccAny(&v, &B);
         /* [Y,colindex]=max(v); */
         Mprhs_[0] = &v;
         Mplhs_[0] = &Y;
         Mplhs_[1] = &colindex;
         mccCallMATLAB(2, Mplhs_, 1, Mprhs_, "max", 70);
         /* [Z,rowindex]=max(B(:,colindex)); */
         mccFindIndex(&IM1_, &colindex);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_CM2_;
            int I_CM2_=1;
            double *q_CM2_;
            double *p_B;
            int I_B=1;
            double *q_B;
            int *p_IM1_;
            int I_IM1_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&B), (mccM(&IM1_) * mccN(&IM1_)));
            mccAllocateMatrix(&CM2_, m_, n_);
            mccCheckMatrixSize(&B, mccM(&B), mccGetMaxIndex(&IM1_ ,mccN(&B)));
            I_CM2_ = (mccM(&CM2_) != 1 || mccN(&CM2_) != 1);
            p_CM2_ = mccPR(&CM2_);
            q_CM2_ = mccPI(&CM2_);
            I_IM1_ = (mccM(&IM1_) != 1 || mccN(&IM1_) != 1);
            p_IM1_ = mccIPR(&IM1_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_IM1_ += I_IM1_)
               {
                  p_B = mccPR(&B) + mccM(&B) * ((int)(*p_IM1_ - .5)) + 0;
                  q_B = mccPI(&B) + mccM(&B) * ((int)(*p_IM1_ - .5)) + 0;
                  for (i_=0; i_<m_; ++i_, p_CM2_+=I_CM2_, q_CM2_+=I_CM2_, p_B+=I_B, q_B+=I_B)
                  {
                     *p_CM2_ = *p_B;
                     *q_CM2_ = *q_B;
                  }
               }
            }
         }
         Mprhs_[0] = &CM2_;
         Mplhs_[0] = &Z;
         Mplhs_[1] = &rowindex;
         mccCallMATLAB(2, Mplhs_, 1, Mprhs_, "max", 71);
         /* start=[rowindex,colindex]; */
         mccCatenateColumns(&start, &rowindex, &colindex);
         /* %introduce window W of 8-neighborhood */
         /* W=[2 3 4;1 0 5; 8 7 6]; */
         mccCopy(&W, &S0_);
         /* %and vector of offsets */
         /* off=[-1 0 1]; */
         mccCopy(&off, &S1_);
         /* %also 2-column matrix with x and y increments for given code */
         /* incr=[0 -1;-1 -1;-1 0;-1 1;0 1;1 1;1 0; 1 -1]; */
         mccCopy(&incr, &S2_);
         
         /* %4. INITIALIZE */
         /* %initializes two column coordinate matrix */
         /* xs=[]; */
         mccCreateEmpty(&xs);
         /* %put frame of 0's around B */
         /* [l1,l2]=size(B); */
         mccGetMatrixSize(&l1,&l2, &B);
         /* Bframe=zeros(l1+2,l2+2); */
         mccZerosMN(&Bframe, (l1+2), (l2+2));
         /* Bframe(2:l1+1,2:l2+1)=B; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_Bframe;
            int I_Bframe=1, J_Bframe;
            double *q_Bframe;
            double *p_B;
            int I_B=1;
            double *q_B;
            m_ = mcmCalcResultSize(m_, &n_, ((int)((l1 + 1) - 2) + 1), ((int)((l2 + 1) - 2) + 1));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&B), mccN(&B));
            mccGrowMatrix(&Bframe, (l1 + 1), (l2 + 1));
            if (mccM(&Bframe) == 1 && mccN(&Bframe) == 1) { I_Bframe = J_Bframe = 0;}
            else { I_Bframe = 1; J_Bframe=mccM(&Bframe)-m_; }
            p_Bframe = mccPR(&Bframe) + ((int)(2 - .5)) + mccM(&Bframe) * ((int)(2 - .5));
            q_Bframe = mccPI(&Bframe) + ((int)(2 - .5)) + mccM(&Bframe) * ((int)(2 - .5));
            I_B = (mccM(&B) != 1 || mccN(&B) != 1);
            p_B = mccPR(&B);
            q_B = mccPI(&B);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_Bframe += J_Bframe, q_Bframe += J_Bframe)
               {
                  for (i_=0; i_<m_; ++i_, p_Bframe+=I_Bframe, q_Bframe+=I_Bframe, p_B+=I_B, q_B+=I_B)
                  {
                     *p_Bframe = *p_B;
                     *q_Bframe = *q_B;
                  }
               }
            }
         }
         
         /* %remember frame! */
         /* start=start+[1 1]; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_start;
            int I_start=1;
            int *p_1start;
            int I_1start=1;
            int *p_S3_;
            int I_S3_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&start), mccN(&start));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&S3_), mccN(&S3_));
            mccAllocateMatrix(&start, m_, n_);
            I_start = (mccM(&start) != 1 || mccN(&start) != 1);
            p_start = mccIPR(&start);
            I_1start = (mccM(&start) != 1 || mccN(&start) != 1);
            p_1start = mccIPR(&start);
            I_S3_ = (mccM(&S3_) != 1 || mccN(&S3_) != 1);
            p_S3_ = mccIPR(&S3_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_start+=I_start, p_1start+=I_1start, p_S3_+=I_S3_)
                  {
                     *p_start = (((int)*p_1start) + ((int)*p_S3_));
                  }
               }
            }
         }
         /* Bframe(start(1),start(2))=0; */
         mccSetMatrixElement(&Bframe, mccRint(((int)mccGetRealVectorElement(&start, mccRint(1)))), mccRint(((int)mccGetRealVectorElement(&start, mccRint(2)))), 0, 0.);
         
         /* %initialize for chain coding */
         /* Neigh=Bframe(start(1)+off,start(2)+off); */
         I0_ = ((int)mccGetRealVectorElement(&start, mccRint(1)));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_IM1_;
            int I_IM1_=1;
            int *p_off;
            int I_off=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
            mccAllocateMatrix(&IM1_, m_, n_);
            I_IM1_ = (mccM(&IM1_) != 1 || mccN(&IM1_) != 1);
            p_IM1_ = mccIPR(&IM1_);
            I_off = (mccM(&off) != 1 || mccN(&off) != 1);
            p_off = mccIPR(&off);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_IM1_+=I_IM1_, p_off+=I_off)
                  {
                     *p_IM1_ = (I0_ + ((int)*p_off));
                  }
               }
            }
         }
         mccFindIndex(&IM0_, &IM1_);
         I1_ = ((int)mccGetRealVectorElement(&start, mccRint(2)));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_IM2_;
            int I_IM2_=1;
            int *p_off;
            int I_off=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
            mccAllocateMatrix(&IM2_, m_, n_);
            I_IM2_ = (mccM(&IM2_) != 1 || mccN(&IM2_) != 1);
            p_IM2_ = mccIPR(&IM2_);
            I_off = (mccM(&off) != 1 || mccN(&off) != 1);
            p_off = mccIPR(&off);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_IM2_+=I_IM2_, p_off+=I_off)
                  {
                     *p_IM2_ = (I1_ + ((int)*p_off));
                  }
               }
            }
         }
         mccFindIndex(&IM3_, &IM2_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_Neigh;
            int I_Neigh=1;
            double *q_Neigh;
            double *p_Bframe;
            double *q_Bframe;
            int *p_IM0_;
            int I_IM0_=1;
            int *p_IM3_;
            int I_IM3_=1;
            m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM0_) * mccN(&IM0_)), (mccM(&IM3_) * mccN(&IM3_)));
            mccAllocateMatrix(&Neigh, m_, n_);
            mccCheckMatrixSize(&Bframe, mccGetMaxIndex(&IM0_ ,mccM(&Bframe)), mccGetMaxIndex(&IM3_ ,mccN(&Bframe)));
            I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
            p_Neigh = mccPR(&Neigh);
            q_Neigh = mccPI(&Neigh);
            I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
            I_IM3_ = (mccM(&IM3_) != 1 || mccN(&IM3_) != 1);
            p_IM3_ = mccIPR(&IM3_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_IM3_ += I_IM3_)
               {
                  p_Bframe = mccPR(&Bframe) + mccM(&Bframe) * ((int)(*p_IM3_ - .5));
                  q_Bframe = mccPI(&Bframe) + mccM(&Bframe) * ((int)(*p_IM3_ - .5));
                  p_IM0_ = mccIPR(&IM0_);
                  for (i_=0; i_<m_; ++i_, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_IM0_+=I_IM0_)
                  {
                     *p_Neigh = p_Bframe[((int)(*p_IM0_ - .5))];
                     *q_Neigh = q_Bframe[((int)(*p_IM0_ - .5))];
                  }
               }
            }
         }
         /* chaincode=Neigh.*W; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_chaincode;
            int I_chaincode=1;
            double *q_chaincode;
            double *p_Neigh;
            int I_Neigh=1;
            double *q_Neigh;
            int *p_W;
            int I_W=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&Neigh), mccN(&Neigh));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&W), mccN(&W));
            mccAllocateMatrix(&chaincode, m_, n_);
            I_chaincode = (mccM(&chaincode) != 1 || mccN(&chaincode) != 1);
            p_chaincode = mccPR(&chaincode);
            q_chaincode = mccPI(&chaincode);
            I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
            p_Neigh = mccPR(&Neigh);
            q_Neigh = mccPI(&Neigh);
            I_W = (mccM(&W) != 1 || mccN(&W) != 1);
            p_W = mccIPR(&W);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_chaincode+=I_chaincode, q_chaincode+=I_chaincode, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_W+=I_W)
                  {
                     *p_chaincode = (*p_Neigh * ((int)*p_W));
                     *q_chaincode = (*q_Neigh * ((int)*p_W));
                  }
               }
            }
         }
         /* code=max(max(chaincode)); */
         Mprhs_[0] = &chaincode;
         Mplhs_[0] = &CM2_;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 95);
         Mprhs_[0] = &CM2_;
         Mplhs_[0] = &code;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 95);
         /* BC=code; */
         mccCopy(&BC, &code);
         /* next=start+incr(code,:); */
         mccFindIndex(&IM3_, &code);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_next;
            int I_next=1;
            int *p_start;
            int I_start=1;
            int *p_incr;
            int I_incr=1, J_incr;
            int *p_IM3_;
            int I_IM3_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&start), mccN(&start));
            m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM3_) * mccN(&IM3_)), mccN(&incr));
            mccAllocateMatrix(&next, m_, n_);
            mccCheckMatrixSize(&incr, mccGetMaxIndex(&IM3_ ,mccM(&incr)), mccN(&incr));
            I_next = (mccM(&next) != 1 || mccN(&next) != 1);
            p_next = mccIPR(&next);
            I_start = (mccM(&start) != 1 || mccN(&start) != 1);
            p_start = mccIPR(&start);
            J_incr = ((mccM(&incr) != 1 || mccN(&incr) != 1) ? mccM(&incr) : 0);
            p_incr = mccIPR(&incr) + mccM(&incr) * 0;
            I_IM3_ = (mccM(&IM3_) != 1 || mccN(&IM3_) != 1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_incr += J_incr)
               {
                  p_IM3_ = mccIPR(&IM3_);
                  for (i_=0; i_<m_; ++i_, p_next+=I_next, p_start+=I_start, p_IM3_+=I_IM3_)
                  {
                     *p_next = (((int)*p_start) + ((int)p_incr[((int)(*p_IM3_ - .5))]));
                  }
               }
            }
         }
         /* Bframe(next(1),next(2))=0; */
         mccSetMatrixElement(&Bframe, mccRint(((int)mccGetRealVectorElement(&next, mccRint(1)))), mccRint(((int)mccGetRealVectorElement(&next, mccRint(2)))), 0, 0.);
         /* xs=[xs;next]; */
         mccCatenateRows(&xs, &xs, &next);
         
         /* Neigh=Bframe(next(1)+off,next(2)+off); */
         I1_ = ((int)mccGetRealVectorElement(&next, mccRint(1)));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_IM3_;
            int I_IM3_=1;
            int *p_off;
            int I_off=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
            mccAllocateMatrix(&IM3_, m_, n_);
            I_IM3_ = (mccM(&IM3_) != 1 || mccN(&IM3_) != 1);
            p_IM3_ = mccIPR(&IM3_);
            I_off = (mccM(&off) != 1 || mccN(&off) != 1);
            p_off = mccIPR(&off);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_IM3_+=I_IM3_, p_off+=I_off)
                  {
                     *p_IM3_ = (I1_ + ((int)*p_off));
                  }
               }
            }
         }
         mccFindIndex(&IM2_, &IM3_);
         I0_ = ((int)mccGetRealVectorElement(&next, mccRint(2)));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_IM0_;
            int I_IM0_=1;
            int *p_off;
            int I_off=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
            mccAllocateMatrix(&IM0_, m_, n_);
            I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
            p_IM0_ = mccIPR(&IM0_);
            I_off = (mccM(&off) != 1 || mccN(&off) != 1);
            p_off = mccIPR(&off);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_IM0_+=I_IM0_, p_off+=I_off)
                  {
                     *p_IM0_ = (I0_ + ((int)*p_off));
                  }
               }
            }
         }
         mccFindIndex(&IM1_, &IM0_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_Neigh;
            int I_Neigh=1;
            double *q_Neigh;
            double *p_Bframe;
            double *q_Bframe;
            int *p_IM2_;
            int I_IM2_=1;
            int *p_IM1_;
            int I_IM1_=1;
            m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM2_) * mccN(&IM2_)), (mccM(&IM1_) * mccN(&IM1_)));
            mccAllocateMatrix(&Neigh, m_, n_);
            mccCheckMatrixSize(&Bframe, mccGetMaxIndex(&IM2_ ,mccM(&Bframe)), mccGetMaxIndex(&IM1_ ,mccN(&Bframe)));
            I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
            p_Neigh = mccPR(&Neigh);
            q_Neigh = mccPI(&Neigh);
            I_IM2_ = (mccM(&IM2_) != 1 || mccN(&IM2_) != 1);
            I_IM1_ = (mccM(&IM1_) != 1 || mccN(&IM1_) != 1);
            p_IM1_ = mccIPR(&IM1_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_IM1_ += I_IM1_)
               {
                  p_Bframe = mccPR(&Bframe) + mccM(&Bframe) * ((int)(*p_IM1_ - .5));
                  q_Bframe = mccPI(&Bframe) + mccM(&Bframe) * ((int)(*p_IM1_ - .5));
                  p_IM2_ = mccIPR(&IM2_);
                  for (i_=0; i_<m_; ++i_, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_IM2_+=I_IM2_)
                  {
                     *p_Neigh = p_Bframe[((int)(*p_IM2_ - .5))];
                     *q_Neigh = q_Bframe[((int)(*p_IM2_ - .5))];
                  }
               }
            }
         }
         /* chaincode=Neigh.*W; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_chaincode;
            int I_chaincode=1;
            double *q_chaincode;
            double *p_Neigh;
            int I_Neigh=1;
            double *q_Neigh;
            int *p_W;
            int I_W=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&Neigh), mccN(&Neigh));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&W), mccN(&W));
            mccAllocateMatrix(&chaincode, m_, n_);
            I_chaincode = (mccM(&chaincode) != 1 || mccN(&chaincode) != 1);
            p_chaincode = mccPR(&chaincode);
            q_chaincode = mccPI(&chaincode);
            I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
            p_Neigh = mccPR(&Neigh);
            q_Neigh = mccPI(&Neigh);
            I_W = (mccM(&W) != 1 || mccN(&W) != 1);
            p_W = mccIPR(&W);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_chaincode+=I_chaincode, q_chaincode+=I_chaincode, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_W+=I_W)
                  {
                     *p_chaincode = (*p_Neigh * ((int)*p_W));
                     *q_chaincode = (*q_Neigh * ((int)*p_W));
                  }
               }
            }
         }
         /* code=max(max(chaincode)); */
         Mprhs_[0] = &chaincode;
         Mplhs_[0] = &CM2_;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 103);
         Mprhs_[0] = &CM2_;
         Mplhs_[0] = &code;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 103);
         /* BC=[BC;code]; */
         mccCatenateRows(&BC, &BC, &code);
         /* Bframe(next(1),next(2))=0; */
         mccSetMatrixElement(&Bframe, mccRint(((int)mccGetRealVectorElement(&next, mccRint(1)))), mccRint(((int)mccGetRealVectorElement(&next, mccRint(2)))), 0, 0.);
         /* next=next+incr(code,:); */
         mccFindIndex(&IM1_, &code);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_next;
            int I_next=1;
            int *p_1next;
            int I_1next=1;
            int *p_incr;
            int I_incr=1, J_incr;
            int *p_IM1_;
            int I_IM1_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&next), mccN(&next));
            m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM1_) * mccN(&IM1_)), mccN(&incr));
            mccAllocateMatrix(&next, m_, n_);
            mccCheckMatrixSize(&incr, mccGetMaxIndex(&IM1_ ,mccM(&incr)), mccN(&incr));
            I_next = (mccM(&next) != 1 || mccN(&next) != 1);
            p_next = mccIPR(&next);
            I_1next = (mccM(&next) != 1 || mccN(&next) != 1);
            p_1next = mccIPR(&next);
            J_incr = ((mccM(&incr) != 1 || mccN(&incr) != 1) ? mccM(&incr) : 0);
            p_incr = mccIPR(&incr) + mccM(&incr) * 0;
            I_IM1_ = (mccM(&IM1_) != 1 || mccN(&IM1_) != 1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_incr += J_incr)
               {
                  p_IM1_ = mccIPR(&IM1_);
                  for (i_=0; i_<m_; ++i_, p_next+=I_next, p_1next+=I_1next, p_IM1_+=I_IM1_)
                  {
                     *p_next = (((int)*p_1next) + ((int)p_incr[((int)(*p_IM1_ - .5))]));
                  }
               }
            }
         }
         /* xs=[xs;next]; */
         mccCatenateRows(&xs, &xs, &next);
         
         /* Neigh=Bframe(next(1)+off,next(2)+off); */
         I0_ = ((int)mccGetRealVectorElement(&next, mccRint(1)));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_IM1_;
            int I_IM1_=1;
            int *p_off;
            int I_off=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
            mccAllocateMatrix(&IM1_, m_, n_);
            I_IM1_ = (mccM(&IM1_) != 1 || mccN(&IM1_) != 1);
            p_IM1_ = mccIPR(&IM1_);
            I_off = (mccM(&off) != 1 || mccN(&off) != 1);
            p_off = mccIPR(&off);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_IM1_+=I_IM1_, p_off+=I_off)
                  {
                     *p_IM1_ = (I0_ + ((int)*p_off));
                  }
               }
            }
         }
         mccFindIndex(&IM0_, &IM1_);
         I1_ = ((int)mccGetRealVectorElement(&next, mccRint(2)));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_IM2_;
            int I_IM2_=1;
            int *p_off;
            int I_off=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
            mccAllocateMatrix(&IM2_, m_, n_);
            I_IM2_ = (mccM(&IM2_) != 1 || mccN(&IM2_) != 1);
            p_IM2_ = mccIPR(&IM2_);
            I_off = (mccM(&off) != 1 || mccN(&off) != 1);
            p_off = mccIPR(&off);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_IM2_+=I_IM2_, p_off+=I_off)
                  {
                     *p_IM2_ = (I1_ + ((int)*p_off));
                  }
               }
            }
         }
         mccFindIndex(&IM3_, &IM2_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_Neigh;
            int I_Neigh=1;
            double *q_Neigh;
            double *p_Bframe;
            double *q_Bframe;
            int *p_IM0_;
            int I_IM0_=1;
            int *p_IM3_;
            int I_IM3_=1;
            m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM0_) * mccN(&IM0_)), (mccM(&IM3_) * mccN(&IM3_)));
            mccAllocateMatrix(&Neigh, m_, n_);
            mccCheckMatrixSize(&Bframe, mccGetMaxIndex(&IM0_ ,mccM(&Bframe)), mccGetMaxIndex(&IM3_ ,mccN(&Bframe)));
            I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
            p_Neigh = mccPR(&Neigh);
            q_Neigh = mccPI(&Neigh);
            I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
            I_IM3_ = (mccM(&IM3_) != 1 || mccN(&IM3_) != 1);
            p_IM3_ = mccIPR(&IM3_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_IM3_ += I_IM3_)
               {
                  p_Bframe = mccPR(&Bframe) + mccM(&Bframe) * ((int)(*p_IM3_ - .5));
                  q_Bframe = mccPI(&Bframe) + mccM(&Bframe) * ((int)(*p_IM3_ - .5));
                  p_IM0_ = mccIPR(&IM0_);
                  for (i_=0; i_<m_; ++i_, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_IM0_+=I_IM0_)
                  {
                     *p_Neigh = p_Bframe[((int)(*p_IM0_ - .5))];
                     *q_Neigh = q_Bframe[((int)(*p_IM0_ - .5))];
                  }
               }
            }
         }
         /* chaincode=Neigh.*W; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_chaincode;
            int I_chaincode=1;
            double *q_chaincode;
            double *p_Neigh;
            int I_Neigh=1;
            double *q_Neigh;
            int *p_W;
            int I_W=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&Neigh), mccN(&Neigh));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&W), mccN(&W));
            mccAllocateMatrix(&chaincode, m_, n_);
            I_chaincode = (mccM(&chaincode) != 1 || mccN(&chaincode) != 1);
            p_chaincode = mccPR(&chaincode);
            q_chaincode = mccPI(&chaincode);
            I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
            p_Neigh = mccPR(&Neigh);
            q_Neigh = mccPI(&Neigh);
            I_W = (mccM(&W) != 1 || mccN(&W) != 1);
            p_W = mccIPR(&W);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_chaincode+=I_chaincode, q_chaincode+=I_chaincode, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_W+=I_W)
                  {
                     *p_chaincode = (*p_Neigh * ((int)*p_W));
                     *q_chaincode = (*q_Neigh * ((int)*p_W));
                  }
               }
            }
         }
         /* code=max(max(chaincode)); */
         Mprhs_[0] = &chaincode;
         Mplhs_[0] = &CM2_;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 111);
         Mprhs_[0] = &CM2_;
         Mplhs_[0] = &code;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 111);
         /* BC=[BC;code]; */
         mccCatenateRows(&BC, &BC, &code);
         mccSTRING(&BC) = 0;
         /* Bframe(next(1),next(2))=0; */
         mccSetMatrixElement(&Bframe, mccRint(((int)mccGetRealVectorElement(&next, mccRint(1)))), mccRint(((int)mccGetRealVectorElement(&next, mccRint(2)))), 0, 0.);
         /* next=next+incr(code,:); */
         mccFindIndex(&IM3_, &code);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_next;
            int I_next=1;
            int *p_1next;
            int I_1next=1;
            int *p_incr;
            int I_incr=1, J_incr;
            int *p_IM3_;
            int I_IM3_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&next), mccN(&next));
            m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM3_) * mccN(&IM3_)), mccN(&incr));
            mccAllocateMatrix(&next, m_, n_);
            mccCheckMatrixSize(&incr, mccGetMaxIndex(&IM3_ ,mccM(&incr)), mccN(&incr));
            I_next = (mccM(&next) != 1 || mccN(&next) != 1);
            p_next = mccIPR(&next);
            I_1next = (mccM(&next) != 1 || mccN(&next) != 1);
            p_1next = mccIPR(&next);
            J_incr = ((mccM(&incr) != 1 || mccN(&incr) != 1) ? mccM(&incr) : 0);
            p_incr = mccIPR(&incr) + mccM(&incr) * 0;
            I_IM3_ = (mccM(&IM3_) != 1 || mccN(&IM3_) != 1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_incr += J_incr)
               {
                  p_IM3_ = mccIPR(&IM3_);
                  for (i_=0; i_<m_; ++i_, p_next+=I_next, p_1next+=I_1next, p_IM3_+=I_IM3_)
                  {
                     *p_next = (((int)*p_1next) + ((int)p_incr[((int)(*p_IM3_ - .5))]));
                  }
               }
            }
         }
         /* xs=[xs;next]; */
         mccCatenateRows(&xs, &xs, &next);
         
         /* 'chaincoding boundary...' */
         mccPrint (&S15_, 0);
         
         /* %5. START CHAINCODING LOOP */
         /* t=3; */
         t = 3;
         /* %start iterating along boundary */
         /* %use sup norm distance for termination criterion */
         /* while max(abs(next-start))>1 */
         while (1)
         {
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               int *p_IM3_;
               int I_IM3_=1;
               int *p_next;
               int I_next=1;
               int *p_start;
               int I_start=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&next), mccN(&next));
               m_ = mcmCalcResultSize(m_, &n_, mccM(&start), mccN(&start));
               mccAllocateMatrix(&IM3_, m_, n_);
               I_IM3_ = (mccM(&IM3_) != 1 || mccN(&IM3_) != 1);
               p_IM3_ = mccIPR(&IM3_);
               I_next = (mccM(&next) != 1 || mccN(&next) != 1);
               p_next = mccIPR(&next);
               I_start = (mccM(&start) != 1 || mccN(&start) != 1);
               p_start = mccIPR(&start);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_IM3_+=I_IM3_, p_next+=I_next, p_start+=I_start)
                     {
                        *p_IM3_ = (((int)*p_next) - ((int)*p_start));
                     }
                  }
               }
            }
            mccAbs(&IM2_, &IM3_);
            mccMax(&IM0_, &IM2_);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               unsigned short *p_BM0_;
               int I_BM0_=1;
               int *p_IM0_;
               int I_IM0_=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&IM0_), mccN(&IM0_));
               mccAllocateMatrix(&BM0_, m_, n_);
               I_BM0_ = (mccM(&BM0_) != 1 || mccN(&BM0_) != 1);
               p_BM0_ = mccSPR(&BM0_);
               I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
               p_IM0_ = mccIPR(&IM0_);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_BM0_+=I_BM0_, p_IM0_+=I_IM0_)
                     {
                        *p_BM0_ = (((int)*p_IM0_) > 1);
                     }
                  }
               }
            }
            if (!(mccIfCondition(&BM0_)))
               break;
            /* Neigh=Bframe(next(1)+off,next(2)+off); */
            I1_ = ((int)mccGetRealVectorElement(&next, mccRint(1)));
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               int *p_IM0_;
               int I_IM0_=1;
               int *p_off;
               int I_off=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
               mccAllocateMatrix(&IM0_, m_, n_);
               I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
               p_IM0_ = mccIPR(&IM0_);
               I_off = (mccM(&off) != 1 || mccN(&off) != 1);
               p_off = mccIPR(&off);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_IM0_+=I_IM0_, p_off+=I_off)
                     {
                        *p_IM0_ = (I1_ + ((int)*p_off));
                     }
                  }
               }
            }
            mccFindIndex(&IM2_, &IM0_);
            I0_ = ((int)mccGetRealVectorElement(&next, mccRint(2)));
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               int *p_IM3_;
               int I_IM3_=1;
               int *p_off;
               int I_off=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
               mccAllocateMatrix(&IM3_, m_, n_);
               I_IM3_ = (mccM(&IM3_) != 1 || mccN(&IM3_) != 1);
               p_IM3_ = mccIPR(&IM3_);
               I_off = (mccM(&off) != 1 || mccN(&off) != 1);
               p_off = mccIPR(&off);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_IM3_+=I_IM3_, p_off+=I_off)
                     {
                        *p_IM3_ = (I0_ + ((int)*p_off));
                     }
                  }
               }
            }
            mccFindIndex(&IM1_, &IM3_);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_Neigh;
               int I_Neigh=1;
               double *q_Neigh;
               double *p_Bframe;
               double *q_Bframe;
               int *p_IM2_;
               int I_IM2_=1;
               int *p_IM1_;
               int I_IM1_=1;
               m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM2_) * mccN(&IM2_)), (mccM(&IM1_) * mccN(&IM1_)));
               mccAllocateMatrix(&Neigh, m_, n_);
               mccCheckMatrixSize(&Bframe, mccGetMaxIndex(&IM2_ ,mccM(&Bframe)), mccGetMaxIndex(&IM1_ ,mccN(&Bframe)));
               I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
               p_Neigh = mccPR(&Neigh);
               q_Neigh = mccPI(&Neigh);
               I_IM2_ = (mccM(&IM2_) != 1 || mccN(&IM2_) != 1);
               I_IM1_ = (mccM(&IM1_) != 1 || mccN(&IM1_) != 1);
               p_IM1_ = mccIPR(&IM1_);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_, p_IM1_ += I_IM1_)
                  {
                     p_Bframe = mccPR(&Bframe) + mccM(&Bframe) * ((int)(*p_IM1_ - .5));
                     q_Bframe = mccPI(&Bframe) + mccM(&Bframe) * ((int)(*p_IM1_ - .5));
                     p_IM2_ = mccIPR(&IM2_);
                     for (i_=0; i_<m_; ++i_, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_IM2_+=I_IM2_)
                     {
                        *p_Neigh = p_Bframe[((int)(*p_IM2_ - .5))];
                        *q_Neigh = q_Bframe[((int)(*p_IM2_ - .5))];
                     }
                  }
               }
            }
            /* chaincode=Neigh.*W; */
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_chaincode;
               int I_chaincode=1;
               double *q_chaincode;
               double *p_Neigh;
               int I_Neigh=1;
               double *q_Neigh;
               int *p_W;
               int I_W=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&Neigh), mccN(&Neigh));
               m_ = mcmCalcResultSize(m_, &n_, mccM(&W), mccN(&W));
               mccAllocateMatrix(&chaincode, m_, n_);
               I_chaincode = (mccM(&chaincode) != 1 || mccN(&chaincode) != 1);
               p_chaincode = mccPR(&chaincode);
               q_chaincode = mccPI(&chaincode);
               I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
               p_Neigh = mccPR(&Neigh);
               q_Neigh = mccPI(&Neigh);
               I_W = (mccM(&W) != 1 || mccN(&W) != 1);
               p_W = mccIPR(&W);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_chaincode+=I_chaincode, q_chaincode+=I_chaincode, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_W+=I_W)
                     {
                        *p_chaincode = (*p_Neigh * ((int)*p_W));
                        *q_chaincode = (*q_Neigh * ((int)*p_W));
                     }
                  }
               }
            }
            /* code=max(max(chaincode)); */
            Mprhs_[0] = &chaincode;
            Mplhs_[0] = &CM2_;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 126);
            Mprhs_[0] = &CM2_;
            Mplhs_[0] = &code;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 126);
            /* BC=[BC;code]; */
            mccCatenateRows(&BC, &BC, &code);
            /* Bframe(next(1),next(2))=0; */
            mccSetMatrixElement(&Bframe, mccRint(((int)mccGetRealVectorElement(&next, mccRint(1)))), mccRint(((int)mccGetRealVectorElement(&next, mccRint(2)))), 0, 0.);
            /* next=next+incr(code,:); */
            mccFindIndex(&IM1_, &code);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               int *p_next;
               int I_next=1;
               int *p_1next;
               int I_1next=1;
               int *p_incr;
               int I_incr=1, J_incr;
               int *p_IM1_;
               int I_IM1_=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&next), mccN(&next));
               m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM1_) * mccN(&IM1_)), mccN(&incr));
               mccAllocateMatrix(&next, m_, n_);
               mccCheckMatrixSize(&incr, mccGetMaxIndex(&IM1_ ,mccM(&incr)), mccN(&incr));
               I_next = (mccM(&next) != 1 || mccN(&next) != 1);
               p_next = mccIPR(&next);
               I_1next = (mccM(&next) != 1 || mccN(&next) != 1);
               p_1next = mccIPR(&next);
               J_incr = ((mccM(&incr) != 1 || mccN(&incr) != 1) ? mccM(&incr) : 0);
               p_incr = mccIPR(&incr) + mccM(&incr) * 0;
               I_IM1_ = (mccM(&IM1_) != 1 || mccN(&IM1_) != 1);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_, p_incr += J_incr)
                  {
                     p_IM1_ = mccIPR(&IM1_);
                     for (i_=0; i_<m_; ++i_, p_next+=I_next, p_1next+=I_1next, p_IM1_+=I_IM1_)
                     {
                        *p_next = (((int)*p_1next) + ((int)p_incr[((int)(*p_IM1_ - .5))]));
                     }
                  }
               }
            }
            /* xs=[xs;next]; */
            mccCatenateRows(&xs, &xs, &next);
            /* t=t+1; */
            t = (t + 1);
            /* end       */
         }
         /* Bframe(start(1),start(2))=1; */
         mccSetMatrixElement(&Bframe, mccRint(((int)mccGetRealVectorElement(&start, mccRint(1)))), mccRint(((int)mccGetRealVectorElement(&start, mccRint(2)))), 1, 0.);
         /* Neigh=Bframe(next(1)+off,next(2)+off); */
         I0_ = ((int)mccGetRealVectorElement(&next, mccRint(1)));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_IM1_;
            int I_IM1_=1;
            int *p_off;
            int I_off=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
            mccAllocateMatrix(&IM1_, m_, n_);
            I_IM1_ = (mccM(&IM1_) != 1 || mccN(&IM1_) != 1);
            p_IM1_ = mccIPR(&IM1_);
            I_off = (mccM(&off) != 1 || mccN(&off) != 1);
            p_off = mccIPR(&off);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_IM1_+=I_IM1_, p_off+=I_off)
                  {
                     *p_IM1_ = (I0_ + ((int)*p_off));
                  }
               }
            }
         }
         mccFindIndex(&IM3_, &IM1_);
         I1_ = ((int)mccGetRealVectorElement(&next, mccRint(2)));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_IM2_;
            int I_IM2_=1;
            int *p_off;
            int I_off=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&off), mccN(&off));
            mccAllocateMatrix(&IM2_, m_, n_);
            I_IM2_ = (mccM(&IM2_) != 1 || mccN(&IM2_) != 1);
            p_IM2_ = mccIPR(&IM2_);
            I_off = (mccM(&off) != 1 || mccN(&off) != 1);
            p_off = mccIPR(&off);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_IM2_+=I_IM2_, p_off+=I_off)
                  {
                     *p_IM2_ = (I1_ + ((int)*p_off));
                  }
               }
            }
         }
         mccFindIndex(&IM0_, &IM2_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_Neigh;
            int I_Neigh=1;
            double *q_Neigh;
            double *p_Bframe;
            double *q_Bframe;
            int *p_IM3_;
            int I_IM3_=1;
            int *p_IM0_;
            int I_IM0_=1;
            m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM3_) * mccN(&IM3_)), (mccM(&IM0_) * mccN(&IM0_)));
            mccAllocateMatrix(&Neigh, m_, n_);
            mccCheckMatrixSize(&Bframe, mccGetMaxIndex(&IM3_ ,mccM(&Bframe)), mccGetMaxIndex(&IM0_ ,mccN(&Bframe)));
            I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
            p_Neigh = mccPR(&Neigh);
            q_Neigh = mccPI(&Neigh);
            I_IM3_ = (mccM(&IM3_) != 1 || mccN(&IM3_) != 1);
            I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
            p_IM0_ = mccIPR(&IM0_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_IM0_ += I_IM0_)
               {
                  p_Bframe = mccPR(&Bframe) + mccM(&Bframe) * ((int)(*p_IM0_ - .5));
                  q_Bframe = mccPI(&Bframe) + mccM(&Bframe) * ((int)(*p_IM0_ - .5));
                  p_IM3_ = mccIPR(&IM3_);
                  for (i_=0; i_<m_; ++i_, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_IM3_+=I_IM3_)
                  {
                     *p_Neigh = p_Bframe[((int)(*p_IM3_ - .5))];
                     *q_Neigh = q_Bframe[((int)(*p_IM3_ - .5))];
                  }
               }
            }
         }
         /* chaincode=Neigh.*W; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_chaincode;
            int I_chaincode=1;
            double *q_chaincode;
            double *p_Neigh;
            int I_Neigh=1;
            double *q_Neigh;
            int *p_W;
            int I_W=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&Neigh), mccN(&Neigh));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&W), mccN(&W));
            mccAllocateMatrix(&chaincode, m_, n_);
            I_chaincode = (mccM(&chaincode) != 1 || mccN(&chaincode) != 1);
            p_chaincode = mccPR(&chaincode);
            q_chaincode = mccPI(&chaincode);
            I_Neigh = (mccM(&Neigh) != 1 || mccN(&Neigh) != 1);
            p_Neigh = mccPR(&Neigh);
            q_Neigh = mccPI(&Neigh);
            I_W = (mccM(&W) != 1 || mccN(&W) != 1);
            p_W = mccIPR(&W);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_chaincode+=I_chaincode, q_chaincode+=I_chaincode, p_Neigh+=I_Neigh, q_Neigh+=I_Neigh, p_W+=I_W)
                  {
                     *p_chaincode = (*p_Neigh * ((int)*p_W));
                     *q_chaincode = (*q_Neigh * ((int)*p_W));
                  }
               }
            }
         }
         /* code=max(max(chaincode)); */
         Mprhs_[0] = &chaincode;
         Mplhs_[0] = &CM2_;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 136);
         Mprhs_[0] = &CM2_;
         Mplhs_[0] = &code;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "max", 136);
         /* BC=[BC;code]; */
         mccCatenateRows(&BC, &BC, &code);
         
         
         
         /* % 6. COMPUTE CURVATURE      */
         /* %compute curvature at equidistant points of ID */
         /* %separated by distance dl */
         /* %computes curvature function rhos */
         /* % for coordinate matrix ys.  */
         /* l=length(xs(:,1)); */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM0_;
            int I_RM0_=1;
            double *p_xs;
            int I_xs=1, J_xs;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&xs), 1);
            mccAllocateMatrix(&RM0_, m_, n_);
            mccCheckMatrixSize(&xs, mccM(&xs), 1);
            I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
            p_RM0_ = mccPR(&RM0_);
            if (mccM(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + 0 + mccM(&xs) * (1-1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_RM0_+=I_RM0_, p_xs+=I_xs)
                  {
                     *p_RM0_ = *p_xs;
                  }
               }
            }
         }
         l = mccGetLength(&RM0_);
         /* % decimate boundary data by a factor d */
         /* d=5; */
         d = 5;
         /* m=floor(l/d); */
         m = floor((l / (double) d));
         /* % extend perimeter by one period */
         /* xs=[xs;xs(1,:)]; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM0_;
            int I_RM0_=1;
            double *p_xs;
            int I_xs=1, J_xs;
            m_ = mcmCalcResultSize(m_, &n_, 1, mccN(&xs));
            mccAllocateMatrix(&RM0_, m_, n_);
            mccCheckMatrixSize(&xs, 1, mccN(&xs));
            I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
            p_RM0_ = mccPR(&RM0_);
            if (mccN(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + (1-1) + mccM(&xs) * 0;
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_RM0_+=I_RM0_, p_xs+=I_xs)
                  {
                     *p_RM0_ = *p_xs;
                  }
               }
            }
         }
         mccSTRING(&RM0_) = mccSTRING(&xs);
         mccCatenateRows(&xs, &xs, &RM0_);
         /* xs=[xs;xs]; */
         mccCatenateRows(&xs, &xs, &xs);
         /* xs=xs(d.*[1:m],:); */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM0_;
            int I_RM0_=1;
            m_ = mcmCalcResultSize(m_, &n_, 1, ((int)(m - 1) + 1));
            mccAllocateMatrix(&RM0_, m_, n_);
            I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
            p_RM0_ = mccPR(&RM0_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_RM0_+=I_RM0_)
                  {
                     *p_RM0_ = (d * (double) ((int)(1 + i_+j_*m_)));
                  }
               }
            }
         }
         mccFindIndex(&IM0_, &RM0_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM1_;
            int I_RM1_=1;
            double *p_xs;
            int I_xs=1, J_xs;
            int *p_IM0_;
            int I_IM0_=1;
            m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM0_) * mccN(&IM0_)), mccN(&xs));
            mccAllocateMatrix(&RM1_, m_, n_);
            mccCheckMatrixSize(&xs, mccGetMaxIndex(&IM0_ ,mccM(&xs)), mccN(&xs));
            I_RM1_ = (mccM(&RM1_) != 1 || mccN(&RM1_) != 1);
            p_RM1_ = mccPR(&RM1_);
            J_xs = ((mccM(&xs) != 1 || mccN(&xs) != 1) ? mccM(&xs) : 0);
            p_xs = mccPR(&xs) + mccM(&xs) * 0;
            I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  p_IM0_ = mccIPR(&IM0_);
                  for (i_=0; i_<m_; ++i_, p_RM1_+=I_RM1_, p_IM0_+=I_IM0_)
                  {
                     *p_RM1_ = p_xs[((int)(*p_IM0_ - .5))];
                  }
               }
            }
         }
         mccSTRING(&RM1_) = mccSTRING(&xs);
         mccCopy(&xs, &RM1_);
         mccSTRING(&xs) = mccSTRING(&RM1_);
         
         /* % 6. REGULARIZE BOUNDARY */
         /* % regularize discrete boundary by cyclic convolution */
         /* l=length(xs(:,1)); */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM1_;
            int I_RM1_=1;
            double *p_xs;
            int I_xs=1, J_xs;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&xs), 1);
            mccAllocateMatrix(&RM1_, m_, n_);
            mccCheckMatrixSize(&xs, mccM(&xs), 1);
            I_RM1_ = (mccM(&RM1_) != 1 || mccN(&RM1_) != 1);
            p_RM1_ = mccPR(&RM1_);
            if (mccM(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + 0 + mccM(&xs) * (1-1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_RM1_+=I_RM1_, p_xs+=I_xs)
                  {
                     *p_RM1_ = *p_xs;
                  }
               }
            }
         }
         l = mccGetLength(&RM1_);
         /* a=.4;p=floor(l/4);v=a.^[0:p-1];w=a.^[p:-1:1]; */
         a = .4;
         p = floor((l / (double) 4));
         mccColon2(&IM0_, (double)0, (p - 1));
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_v_1;
            int I_v_1=1;
            int *p_IM0_;
            int I_IM0_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&IM0_), mccN(&IM0_));
            mccAllocateMatrix(&v_1, m_, n_);
            I_v_1 = (mccM(&v_1) != 1 || mccN(&v_1) != 1);
            p_v_1 = mccPR(&v_1);
            I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
            p_IM0_ = mccIPR(&IM0_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_v_1+=I_v_1, p_IM0_+=I_IM0_)
                  {
                     *p_v_1 = mcmRealPowerInt(a, ((int)*p_IM0_));
                  }
               }
            }
         }
         mccColon(&RM1_, p, (double) -1, (double)1);
         mccArrayPower(&w, mccTempMatrix(a, 0., mccREAL, 0 ), &RM1_);
         /* u=[v,zeros(1,l-2*p),w];u=u./sum(u); */
         mccZerosMN(&IM0_, 1, ((int)(l - (2 * (double) p))));
         mccCatenateColumns(&RM1_, &v_1, &IM0_);
         mccCatenateColumns(&u, &RM1_, &w);
         Mprhs_[0] = &u;
         Mplhs_[0] = &CM2_;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "sum", 159);
         mccArrayRightDivide(&u, &u, &CM2_);
         /* uhat=fft(u);x1shat=fft(xs(:,1)');x2shat=fft(xs(:,2)'); */
         Mprhs_[0] = &u;
         Mplhs_[0] = &uhat;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "fft", 160);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM1_;
            int I_RM1_=1;
            double *p_xs;
            int I_xs=1, J_xs;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&xs), 1);
            mccAllocateMatrix(&RM1_, m_, n_);
            mccCheckMatrixSize(&xs, mccM(&xs), 1);
            I_RM1_ = (mccM(&RM1_) != 1 || mccN(&RM1_) != 1);
            p_RM1_ = mccPR(&RM1_);
            if (mccM(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + 0 + mccM(&xs) * (1-1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_RM1_+=I_RM1_, p_xs+=I_xs)
                  {
                     *p_RM1_ = *p_xs;
                  }
               }
            }
         }
         mccConjTrans(&RM0_, &RM1_);
         Mprhs_[0] = &RM0_;
         Mplhs_[0] = &x1shat;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "fft", 160);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM0_;
            int I_RM0_=1;
            double *p_xs;
            int I_xs=1, J_xs;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&xs), 1);
            mccAllocateMatrix(&RM0_, m_, n_);
            mccCheckMatrixSize(&xs, mccM(&xs), 2);
            I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
            p_RM0_ = mccPR(&RM0_);
            if (mccM(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + 0 + mccM(&xs) * (2-1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_RM0_+=I_RM0_, p_xs+=I_xs)
                  {
                     *p_RM0_ = *p_xs;
                  }
               }
            }
         }
         mccConjTrans(&RM1_, &RM0_);
         Mprhs_[0] = &RM1_;
         Mplhs_[0] = &x2shat;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "fft", 160);
         /* % execute convolution by Fourier transform */
         /* xs(:,1)=real(ifft(uhat.*x1shat))';xs(:,2)=real(ifft(uhat.*x2shat))'; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_CM2_;
            int I_CM2_=1;
            double *q_CM2_;
            double *p_uhat;
            int I_uhat=1;
            double *q_uhat;
            double *p_x1shat;
            int I_x1shat=1;
            double *q_x1shat;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&uhat), mccN(&uhat));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&x1shat), mccN(&x1shat));
            mccAllocateMatrix(&CM2_, m_, n_);
            I_CM2_ = (mccM(&CM2_) != 1 || mccN(&CM2_) != 1);
            p_CM2_ = mccPR(&CM2_);
            q_CM2_ = mccPI(&CM2_);
            I_uhat = (mccM(&uhat) != 1 || mccN(&uhat) != 1);
            p_uhat = mccPR(&uhat);
            q_uhat = mccPI(&uhat);
            I_x1shat = (mccM(&x1shat) != 1 || mccN(&x1shat) != 1);
            p_x1shat = mccPR(&x1shat);
            q_x1shat = mccPI(&x1shat);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_CM2_+=I_CM2_, q_CM2_+=I_CM2_, p_uhat+=I_uhat, q_uhat+=I_uhat, p_x1shat+=I_x1shat, q_x1shat+=I_x1shat)
                  {
                     {
                        double t_ = (*p_uhat * *p_x1shat - *q_uhat * *q_x1shat);
                        *q_CM2_ = (*p_uhat * *q_x1shat + *q_uhat * *p_x1shat);
                        *p_CM2_ = t_;
                     }
                  }
               }
            }
         }
         Mprhs_[0] = &CM2_;
         Mplhs_[0] = &CM3_;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "ifft", 162);
         mccReal(&RM1_, &CM3_);
         mccConjTrans(&RM0_, &RM1_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_xs;
            int I_xs=1, J_xs;
            double *p_RM0_;
            int I_RM0_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&RM0_), mccN(&RM0_));
            mccGrowMatrix(&xs, m_, 1);
            if (mccM(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + 0 + mccM(&xs) * (1-1);
            I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
            p_RM0_ = mccPR(&RM0_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_xs+=I_xs, p_RM0_+=I_RM0_)
                  {
                     *p_xs = *p_RM0_;
                  }
               }
            }
         }
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_CM3_;
            int I_CM3_=1;
            double *q_CM3_;
            double *p_uhat;
            int I_uhat=1;
            double *q_uhat;
            double *p_x2shat;
            int I_x2shat=1;
            double *q_x2shat;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&uhat), mccN(&uhat));
            m_ = mcmCalcResultSize(m_, &n_, mccM(&x2shat), mccN(&x2shat));
            mccAllocateMatrix(&CM3_, m_, n_);
            I_CM3_ = (mccM(&CM3_) != 1 || mccN(&CM3_) != 1);
            p_CM3_ = mccPR(&CM3_);
            q_CM3_ = mccPI(&CM3_);
            I_uhat = (mccM(&uhat) != 1 || mccN(&uhat) != 1);
            p_uhat = mccPR(&uhat);
            q_uhat = mccPI(&uhat);
            I_x2shat = (mccM(&x2shat) != 1 || mccN(&x2shat) != 1);
            p_x2shat = mccPR(&x2shat);
            q_x2shat = mccPI(&x2shat);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_CM3_+=I_CM3_, q_CM3_+=I_CM3_, p_uhat+=I_uhat, q_uhat+=I_uhat, p_x2shat+=I_x2shat, q_x2shat+=I_x2shat)
                  {
                     {
                        double t_ = (*p_uhat * *p_x2shat - *q_uhat * *q_x2shat);
                        *q_CM3_ = (*p_uhat * *q_x2shat + *q_uhat * *p_x2shat);
                        *p_CM3_ = t_;
                     }
                  }
               }
            }
         }
         Mprhs_[0] = &CM3_;
         Mplhs_[0] = &CM2_;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "ifft", 162);
         mccReal(&RM0_, &CM2_);
         mccConjTrans(&RM1_, &RM0_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_xs;
            int I_xs=1, J_xs;
            double *p_RM1_;
            int I_RM1_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&RM1_), mccN(&RM1_));
            mccGrowMatrix(&xs, m_, 2);
            if (mccM(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + 0 + mccM(&xs) * (2-1);
            I_RM1_ = (mccM(&RM1_) != 1 || mccN(&RM1_) != 1);
            p_RM1_ = mccPR(&RM1_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_xs+=I_xs, p_RM1_+=I_RM1_)
                  {
                     *p_xs = *p_RM1_;
                  }
               }
            }
         }
         /* % extend boundary data one period */
         /* xs=[xs;xs]; */
         mccCatenateRows(&xs, &xs, &xs);
         
         /* 'computing curvature...' */
         mccPrint (&S16_, 0);
         
         /* % 7. GET ARCLENGTH, CURVATURES */
         /* % compute arclengths */
         /* lengths=zeros(1,m); */
         mccZerosMN(&lengths, 1, ((int)m));
         /* for t=1:m */
         for (I1_ = 1; I1_ <= m; I1_ = I1_ + 1)
         {
            t = I1_;
            /* lengths(t)=norm(xs(t+1,:)-xs(t,:)); */
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_RM1_;
               int I_RM1_=1;
               double *p_xs;
               int I_xs=1, J_xs;
               double *p_1xs;
               int I_1xs=1, J_1xs;
               m_ = mcmCalcResultSize(m_, &n_, 1, mccN(&xs));
               m_ = mcmCalcResultSize(m_, &n_, 1, mccN(&xs));
               mccAllocateMatrix(&RM1_, m_, n_);
               mccCheckMatrixSize(&xs, (t+1), mccN(&xs));
               mccCheckMatrixSize(&xs, t, mccN(&xs));
               I_RM1_ = (mccM(&RM1_) != 1 || mccN(&RM1_) != 1);
               p_RM1_ = mccPR(&RM1_);
               if (mccN(&xs) == 1) { I_xs = J_xs = 0;}
               else { I_xs = 1; J_xs=mccM(&xs)-m_; }
               p_xs = mccPR(&xs) + ((t+1)-1) + mccM(&xs) * 0;
               if (mccN(&xs) == 1) { I_1xs = J_1xs = 0;}
               else { I_1xs = 1; J_1xs=mccM(&xs)-m_; }
               p_1xs = mccPR(&xs) + (t-1) + mccM(&xs) * 0;
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_, p_xs += J_xs, p_1xs += J_1xs)
                  {
                     for (i_=0; i_<m_; ++i_, p_RM1_+=I_RM1_, p_xs+=I_xs, p_1xs+=I_1xs)
                     {
                        *p_RM1_ = (*p_xs - *p_1xs);
                     }
                  }
               }
            }
            Mprhs_[0] = &RM1_;
            Mplhs_[0] = 0;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "norm", 172);
            R0_ = mccImportReal(0, 0, Mplhs_[ 0 ], 1, " (vision_engine_b, line 172): R0_");
            mccSetRealVectorElement(&lengths, mccRint(t), R0_);
            /* end */
         }
         
         /* cuml=cumsum(lengths); */
         Mprhs_[0] = &lengths;
         Mplhs_[0] = &cuml;
         mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "cumsum", 175);
         /* % now get discretized curvatures rhos */
         /* dphis=zeros(1,m);rhos=zeros(1,m); */
         mccZerosMN(&dphis, 1, ((int)m));
         mccZerosMN(&rhos, 1, ((int)m));
         
         /* for t=1:m */
         for (I1_ = 1; I1_ <= m; I1_ = I1_ + 1)
         {
            t = I1_;
            
            /* a=xs(mod(t+3,m)+1,:)-xs(t,:);a=a./norm(a); */
            Mprhs_[0] = mccTempMatrix((t + 3), 0., mccINT, 0 );
            Mprhs_[1] = mccTempMatrix(m, 0., mccREAL, 0 );
            Mplhs_[0] = &RM1_;
            mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "mod", 181);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_RM0_;
               int I_RM0_=1;
               double *p_RM1_;
               int I_RM1_=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&RM1_), mccN(&RM1_));
               mccAllocateMatrix(&RM0_, m_, n_);
               I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
               p_RM0_ = mccPR(&RM0_);
               I_RM1_ = (mccM(&RM1_) != 1 || mccN(&RM1_) != 1);
               p_RM1_ = mccPR(&RM1_);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_RM0_+=I_RM0_, p_RM1_+=I_RM1_)
                     {
                        *p_RM0_ = (*p_RM1_ + 1);
                     }
                  }
               }
            }
            mccFindIndex(&IM0_, &RM0_);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_a_1;
               int I_a_1=1;
               double *p_xs;
               int I_xs=1, J_xs;
               int *p_IM0_;
               int I_IM0_=1;
               double *p_1xs;
               int I_1xs=1, J_1xs;
               m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM0_) * mccN(&IM0_)), mccN(&xs));
               m_ = mcmCalcResultSize(m_, &n_, 1, mccN(&xs));
               mccAllocateMatrix(&a_1, m_, n_);
               mccCheckMatrixSize(&xs, mccGetMaxIndex(&IM0_ ,mccM(&xs)), mccN(&xs));
               mccCheckMatrixSize(&xs, t, mccN(&xs));
               I_a_1 = (mccM(&a_1) != 1 || mccN(&a_1) != 1);
               p_a_1 = mccPR(&a_1);
               J_xs = ((mccM(&xs) != 1 || mccN(&xs) != 1) ? mccM(&xs) : 0);
               p_xs = mccPR(&xs) + mccM(&xs) * 0;
               I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
               if (mccN(&xs) == 1) { I_1xs = J_1xs = 0;}
               else { I_1xs = 1; J_1xs=mccM(&xs)-m_; }
               p_1xs = mccPR(&xs) + (t-1) + mccM(&xs) * 0;
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_, p_xs += J_xs, p_1xs += J_1xs)
                  {
                     p_IM0_ = mccIPR(&IM0_);
                     for (i_=0; i_<m_; ++i_, p_a_1+=I_a_1, p_IM0_+=I_IM0_, p_1xs+=I_1xs)
                     {
                        *p_a_1 = (p_xs[((int)(*p_IM0_ - .5))] - *p_1xs);
                     }
                  }
               }
            }
            Mprhs_[0] = &a_1;
            Mplhs_[0] = 0;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "norm", 181);
            R0_ = mccImportReal(0, 0, Mplhs_[ 0 ], 1, " (vision_engine_b, line 181): R0_");
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_a_1;
               int I_a_1=1;
               double *p_1a_1;
               int I_1a_1=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&a_1), mccN(&a_1));
               mccAllocateMatrix(&a_1, m_, n_);
               I_a_1 = (mccM(&a_1) != 1 || mccN(&a_1) != 1);
               p_a_1 = mccPR(&a_1);
               I_1a_1 = (mccM(&a_1) != 1 || mccN(&a_1) != 1);
               p_1a_1 = mccPR(&a_1);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_a_1+=I_a_1, p_1a_1+=I_1a_1)
                     {
                        *p_a_1 = (*p_1a_1 / (double) R0_);
                     }
                  }
               }
            }
            /* b=xs(t,:)-xs(mod(t-5,m)+1,:);b=b./norm(b); */
            Mprhs_[0] = mccTempMatrix((t - 5), 0., mccINT, 0 );
            Mprhs_[1] = mccTempMatrix(m, 0., mccREAL, 0 );
            Mplhs_[0] = &RM0_;
            mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "mod", 182);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_RM1_;
               int I_RM1_=1;
               double *p_RM0_;
               int I_RM0_=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&RM0_), mccN(&RM0_));
               mccAllocateMatrix(&RM1_, m_, n_);
               I_RM1_ = (mccM(&RM1_) != 1 || mccN(&RM1_) != 1);
               p_RM1_ = mccPR(&RM1_);
               I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
               p_RM0_ = mccPR(&RM0_);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_RM1_+=I_RM1_, p_RM0_+=I_RM0_)
                     {
                        *p_RM1_ = (*p_RM0_ + 1);
                     }
                  }
               }
            }
            mccFindIndex(&IM0_, &RM1_);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_b;
               int I_b=1;
               double *p_xs;
               int I_xs=1, J_xs;
               double *p_1xs;
               int I_1xs=1, J_1xs;
               int *p_IM0_;
               int I_IM0_=1;
               m_ = mcmCalcResultSize(m_, &n_, 1, mccN(&xs));
               m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM0_) * mccN(&IM0_)), mccN(&xs));
               mccAllocateMatrix(&b, m_, n_);
               mccCheckMatrixSize(&xs, t, mccN(&xs));
               mccCheckMatrixSize(&xs, mccGetMaxIndex(&IM0_ ,mccM(&xs)), mccN(&xs));
               I_b = (mccM(&b) != 1 || mccN(&b) != 1);
               p_b = mccPR(&b);
               if (mccN(&xs) == 1) { I_xs = J_xs = 0;}
               else { I_xs = 1; J_xs=mccM(&xs)-m_; }
               p_xs = mccPR(&xs) + (t-1) + mccM(&xs) * 0;
               J_1xs = ((mccM(&xs) != 1 || mccN(&xs) != 1) ? mccM(&xs) : 0);
               p_1xs = mccPR(&xs) + mccM(&xs) * 0;
               I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_, p_xs += J_xs, p_1xs += J_1xs)
                  {
                     p_IM0_ = mccIPR(&IM0_);
                     for (i_=0; i_<m_; ++i_, p_b+=I_b, p_xs+=I_xs, p_IM0_+=I_IM0_)
                     {
                        *p_b = (*p_xs - p_1xs[((int)(*p_IM0_ - .5))]);
                     }
                  }
               }
            }
            Mprhs_[0] = &b;
            Mplhs_[0] = 0;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "norm", 182);
            R0_ = mccImportReal(0, 0, Mplhs_[ 0 ], 1, " (vision_engine_b, line 182): R0_");
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_b;
               int I_b=1;
               double *p_1b;
               int I_1b=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&b), mccN(&b));
               mccAllocateMatrix(&b, m_, n_);
               I_b = (mccM(&b) != 1 || mccN(&b) != 1);
               p_b = mccPR(&b);
               I_1b = (mccM(&b) != 1 || mccN(&b) != 1);
               p_1b = mccPR(&b);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_b+=I_b, p_1b+=I_1b)
                     {
                        *p_b = (*p_1b / (double) R0_);
                     }
                  }
               }
            }
            /* %change to image (matrix) coordinates */
            /* a=[a(2), -a(1)];b=[b(2),-b(1)]; */
            mccCatenateColumns(&a_1, mccTempVectorElement(&a_1, 2), mccTempMatrix((-(mccGetRealVectorElement(&a_1, mccRint(1)))), 0., mccREAL, 0 ));
            mccCatenateColumns(&b, mccTempVectorElement(&b, 2), mccTempMatrix((-(mccGetRealVectorElement(&b, mccRint(1)))), 0., mccREAL, 0 ));
            
            
            /* x=a(1)*b(1)+a(2)*b(2); */
            x = (((mccGetRealVectorElement(&a_1, mccRint(1))) * (double) (mccGetRealVectorElement(&b, mccRint(1)))) + ((mccGetRealVectorElement(&a_1, mccRint(2))) * (double) (mccGetRealVectorElement(&b, mccRint(2)))));
            /* y=a(1)*b(2)-a(2)*b(1); */
            y = (((mccGetRealVectorElement(&a_1, mccRint(1))) * (double) (mccGetRealVectorElement(&b, mccRint(2)))) - ((mccGetRealVectorElement(&a_1, mccRint(2))) * (double) (mccGetRealVectorElement(&b, mccRint(1)))));
            /* r=sqrt(x^2+y^2); */
            Mprhs_[0] = mccTempMatrix((mcmRealPowerInt(x, 2) + mcmRealPowerInt(y, 2)), 0., mccREAL, 0 );
            Mplhs_[0] = &CM2_;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "sqrt", 189);
            r_r = (mccGetRealVectorElement(&CM2_, mccRint(1)));
            r_i = mccGetImagVectorElement(&CM2_, mccRint(1));
            /* %dphis(t)=asin(y/r); */
            /* dphis(t)=(x>=0)*asin(y/r)+(x<0)*(pi-asin(y/r)); */
            {
               double t_ = mcmDivideRealpart(y, 0., r_r, r_i);
               C0__i = mcmDivideImagPart(y, 0., r_r, r_i);
               C0__r = t_;
            }
            Mprhs_[0] = mccTempMatrix(C0__r, C0__i, mccCX, 0 );
            Mplhs_[0] = &CM3_;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "asin", 191);
            C1__r = (mccGetRealVectorElement(&CM3_, mccRint(1)));
            C1__i = mccGetImagVectorElement(&CM3_, mccRint(1));
            R0_ = mcmPi();
            {
               double t_ = mcmDivideRealpart(y, 0., r_r, r_i);
               C2__i = mcmDivideImagPart(y, 0., r_r, r_i);
               C2__r = t_;
            }
            Mprhs_[0] = mccTempMatrix(C2__r, C2__i, mccCX, 0 );
            Mplhs_[0] = &CM4_;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "asin", 191);
            C3__r = (mccGetRealVectorElement(&CM4_, mccRint(1)));
            C3__i = mccGetImagVectorElement(&CM4_, mccRint(1));
            mccSetVectorElement(&dphis, mccRint(t), ((( (x >= 0) && !mccREL_NAN(x) ) * C1__r) + (( (x < 0) && !mccREL_NAN(x) ) * (R0_ - C3__r))), ((( (x >= 0) && !mccREL_NAN(x) ) * C1__i) + (( (x < 0) && !mccREL_NAN(x) ) * (0. - C3__i))));
            /* rhos(t)=2*dphis(t)/(norm(xs(mod(t+3,m)+1,:)-xs(t,:))+norm(xs(t,:)-xs(mod(t-5,m)+1,:))); */
            Mprhs_[0] = mccTempMatrix((t + 3), 0., mccINT, 0 );
            Mprhs_[1] = mccTempMatrix(m, 0., mccREAL, 0 );
            Mplhs_[0] = &RM1_;
            mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "mod", 192);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_RM0_;
               int I_RM0_=1;
               double *p_RM1_;
               int I_RM1_=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&RM1_), mccN(&RM1_));
               mccAllocateMatrix(&RM0_, m_, n_);
               I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
               p_RM0_ = mccPR(&RM0_);
               I_RM1_ = (mccM(&RM1_) != 1 || mccN(&RM1_) != 1);
               p_RM1_ = mccPR(&RM1_);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_RM0_+=I_RM0_, p_RM1_+=I_RM1_)
                     {
                        *p_RM0_ = (*p_RM1_ + 1);
                     }
                  }
               }
            }
            mccFindIndex(&IM0_, &RM0_);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_RM2_;
               int I_RM2_=1;
               double *p_xs;
               int I_xs=1, J_xs;
               int *p_IM0_;
               int I_IM0_=1;
               double *p_1xs;
               int I_1xs=1, J_1xs;
               m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM0_) * mccN(&IM0_)), mccN(&xs));
               m_ = mcmCalcResultSize(m_, &n_, 1, mccN(&xs));
               mccAllocateMatrix(&RM2_, m_, n_);
               mccCheckMatrixSize(&xs, mccGetMaxIndex(&IM0_ ,mccM(&xs)), mccN(&xs));
               mccCheckMatrixSize(&xs, t, mccN(&xs));
               I_RM2_ = (mccM(&RM2_) != 1 || mccN(&RM2_) != 1);
               p_RM2_ = mccPR(&RM2_);
               J_xs = ((mccM(&xs) != 1 || mccN(&xs) != 1) ? mccM(&xs) : 0);
               p_xs = mccPR(&xs) + mccM(&xs) * 0;
               I_IM0_ = (mccM(&IM0_) != 1 || mccN(&IM0_) != 1);
               if (mccN(&xs) == 1) { I_1xs = J_1xs = 0;}
               else { I_1xs = 1; J_1xs=mccM(&xs)-m_; }
               p_1xs = mccPR(&xs) + (t-1) + mccM(&xs) * 0;
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_, p_xs += J_xs, p_1xs += J_1xs)
                  {
                     p_IM0_ = mccIPR(&IM0_);
                     for (i_=0; i_<m_; ++i_, p_RM2_+=I_RM2_, p_IM0_+=I_IM0_, p_1xs+=I_1xs)
                     {
                        *p_RM2_ = (p_xs[((int)(*p_IM0_ - .5))] - *p_1xs);
                     }
                  }
               }
            }
            Mprhs_[0] = &RM2_;
            Mplhs_[0] = 0;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "norm", 192);
            R0_ = mccImportReal(0, 0, Mplhs_[ 0 ], 1, " (vision_engine_b, line 192): R0_");
            Mprhs_[0] = mccTempMatrix((t - 5), 0., mccINT, 0 );
            Mprhs_[1] = mccTempMatrix(m, 0., mccREAL, 0 );
            Mplhs_[0] = &RM3_;
            mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "mod", 192);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_RM4_;
               int I_RM4_=1;
               double *p_RM3_;
               int I_RM3_=1;
               m_ = mcmCalcResultSize(m_, &n_, mccM(&RM3_), mccN(&RM3_));
               mccAllocateMatrix(&RM4_, m_, n_);
               I_RM4_ = (mccM(&RM4_) != 1 || mccN(&RM4_) != 1);
               p_RM4_ = mccPR(&RM4_);
               I_RM3_ = (mccM(&RM3_) != 1 || mccN(&RM3_) != 1);
               p_RM3_ = mccPR(&RM3_);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_)
                  {
                     for (i_=0; i_<m_; ++i_, p_RM4_+=I_RM4_, p_RM3_+=I_RM3_)
                     {
                        *p_RM4_ = (*p_RM3_ + 1);
                     }
                  }
               }
            }
            mccFindIndex(&IM2_, &RM4_);
            {
               int i_, j_;
               int m_=1, n_=1, cx_ = 0;
               double *p_RM5_;
               int I_RM5_=1;
               double *p_xs;
               int I_xs=1, J_xs;
               double *p_1xs;
               int I_1xs=1, J_1xs;
               int *p_IM2_;
               int I_IM2_=1;
               m_ = mcmCalcResultSize(m_, &n_, 1, mccN(&xs));
               m_ = mcmCalcResultSize(m_, &n_, (mccM(&IM2_) * mccN(&IM2_)), mccN(&xs));
               mccAllocateMatrix(&RM5_, m_, n_);
               mccCheckMatrixSize(&xs, t, mccN(&xs));
               mccCheckMatrixSize(&xs, mccGetMaxIndex(&IM2_ ,mccM(&xs)), mccN(&xs));
               I_RM5_ = (mccM(&RM5_) != 1 || mccN(&RM5_) != 1);
               p_RM5_ = mccPR(&RM5_);
               if (mccN(&xs) == 1) { I_xs = J_xs = 0;}
               else { I_xs = 1; J_xs=mccM(&xs)-m_; }
               p_xs = mccPR(&xs) + (t-1) + mccM(&xs) * 0;
               J_1xs = ((mccM(&xs) != 1 || mccN(&xs) != 1) ? mccM(&xs) : 0);
               p_1xs = mccPR(&xs) + mccM(&xs) * 0;
               I_IM2_ = (mccM(&IM2_) != 1 || mccN(&IM2_) != 1);
               if (m_ != 0)
               {
                  for (j_=0; j_<n_; ++j_, p_xs += J_xs, p_1xs += J_1xs)
                  {
                     p_IM2_ = mccIPR(&IM2_);
                     for (i_=0; i_<m_; ++i_, p_RM5_+=I_RM5_, p_xs+=I_xs, p_IM2_+=I_IM2_)
                     {
                        *p_RM5_ = (*p_xs - p_1xs[((int)(*p_IM2_ - .5))]);
                     }
                  }
               }
            }
            Mprhs_[0] = &RM5_;
            Mplhs_[0] = 0;
            mccCallMATLAB(1, Mplhs_, 1, Mprhs_, "norm", 192);
            R1_ = mccImportReal(0, 0, Mplhs_[ 0 ], 1, " (vision_engine_b, line 192): R1_");
            mccSetVectorElement(&rhos, mccRint(t), ((2 * (mccGetRealVectorElement(&dphis, mccRint(t)))) / (double) (R0_ + R1_)), ((2 * mccGetImagVectorElement(&dphis, mccRint(t))) / (double) (R0_ + R1_)));
            
            /* end */
         }
         /* %arc element length=dl */
         
         /* dl=5; */
         dl = 5;
         /* n=floor(cuml(m)/dl); */
         n = floor(((mccGetRealVectorElement(&cuml, mccRint(m))) / (double) dl));
         /* %interpolate curvature function to equi-distant points */
         /* rhos=spline([0,cuml(1:m-1)],rhos,dl.*[0:n-1]);  */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM5_;
            int I_RM5_=1;
            double *p_cuml;
            int I_cuml=1;
            m_ = mccCalcSubscriptDimensions(m_, &n_, 1, ((int)((m - 1) - 1) + 1), &cuml);
            mccAllocateMatrix(&RM5_, m_, n_);
            mccCheckVectorSize(&cuml, (m - 1));
            I_RM5_ = (mccM(&RM5_) != 1 || mccN(&RM5_) != 1);
            p_RM5_ = mccPR(&RM5_);
            I_cuml = (mccM(&cuml) != 1 || mccN(&cuml) != 1);
            p_cuml = mccPR(&cuml) + ((int)(1 - .5));
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_RM5_+=I_RM5_, p_cuml+=I_cuml)
                  {
                     *p_RM5_ = *p_cuml;
                  }
               }
            }
         }
         mccCatenateColumns(&RM4_, mccTempMatrix(0, 0., mccINT, 0 ), &RM5_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM3_;
            int I_RM3_=1;
            m_ = mcmCalcResultSize(m_, &n_, 1, ((int)((n - 1) - 0) + 1));
            mccAllocateMatrix(&RM3_, m_, n_);
            I_RM3_ = (mccM(&RM3_) != 1 || mccN(&RM3_) != 1);
            p_RM3_ = mccPR(&RM3_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_RM3_+=I_RM3_)
                  {
                     *p_RM3_ = (dl * (double) ((int)(0 + i_+j_*m_)));
                  }
               }
            }
         }
         Mprhs_[0] = &RM4_;
         Mprhs_[1] = &rhos;
         Mprhs_[2] = &RM3_;
         Mplhs_[0] = &rhos;
         mccCallMATLAB(1, Mplhs_, 3, Mprhs_, "spline", 200);
         /* ys=zeros(n,2); */
         mccZerosMN(&ys, n, 2);
         /* y1=spline([0,cuml(1:m-1)],xs(1:m,1)',dl.*[0:n-1]);  */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM3_;
            int I_RM3_=1;
            double *p_cuml;
            int I_cuml=1;
            m_ = mccCalcSubscriptDimensions(m_, &n_, 1, ((int)((m - 1) - 1) + 1), &cuml);
            mccAllocateMatrix(&RM3_, m_, n_);
            mccCheckVectorSize(&cuml, (m - 1));
            I_RM3_ = (mccM(&RM3_) != 1 || mccN(&RM3_) != 1);
            p_RM3_ = mccPR(&RM3_);
            I_cuml = (mccM(&cuml) != 1 || mccN(&cuml) != 1);
            p_cuml = mccPR(&cuml) + ((int)(1 - .5));
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_RM3_+=I_RM3_, p_cuml+=I_cuml)
                  {
                     *p_RM3_ = *p_cuml;
                  }
               }
            }
         }
         mccCatenateColumns(&RM4_, mccTempMatrix(0, 0., mccINT, 0 ), &RM3_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM5_;
            int I_RM5_=1;
            double *p_xs;
            int I_xs=1, J_xs;
            m_ = mcmCalcResultSize(m_, &n_, ((int)(m - 1) + 1), 1);
            mccAllocateMatrix(&RM5_, m_, n_);
            mccCheckMatrixSize(&xs, (((int)(1 + .5)) + ((int)(m - 1))), 1);
            I_RM5_ = (mccM(&RM5_) != 1 || mccN(&RM5_) != 1);
            p_RM5_ = mccPR(&RM5_);
            if (mccM(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + ((int)(1 - .5)) + mccM(&xs) * (1-1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_RM5_+=I_RM5_, p_xs+=I_xs)
                  {
                     *p_RM5_ = *p_xs;
                  }
               }
            }
         }
         mccConjTrans(&RM2_, &RM5_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM0_;
            int I_RM0_=1;
            m_ = mcmCalcResultSize(m_, &n_, 1, ((int)((n - 1) - 0) + 1));
            mccAllocateMatrix(&RM0_, m_, n_);
            I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
            p_RM0_ = mccPR(&RM0_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_RM0_+=I_RM0_)
                  {
                     *p_RM0_ = (dl * (double) ((int)(0 + i_+j_*m_)));
                  }
               }
            }
         }
         Mprhs_[0] = &RM4_;
         Mprhs_[1] = &RM2_;
         Mprhs_[2] = &RM0_;
         Mplhs_[0] = &M_y1;
         mccCallMATLAB(1, Mplhs_, 3, Mprhs_, "spline", 202);
         /* y2=spline([0,cuml(1:m-1)],xs(1:m,2)',dl.*[0:n-1]);  */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM0_;
            int I_RM0_=1;
            double *p_cuml;
            int I_cuml=1;
            m_ = mccCalcSubscriptDimensions(m_, &n_, 1, ((int)((m - 1) - 1) + 1), &cuml);
            mccAllocateMatrix(&RM0_, m_, n_);
            mccCheckVectorSize(&cuml, (m - 1));
            I_RM0_ = (mccM(&RM0_) != 1 || mccN(&RM0_) != 1);
            p_RM0_ = mccPR(&RM0_);
            I_cuml = (mccM(&cuml) != 1 || mccN(&cuml) != 1);
            p_cuml = mccPR(&cuml) + ((int)(1 - .5));
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_RM0_+=I_RM0_, p_cuml+=I_cuml)
                  {
                     *p_RM0_ = *p_cuml;
                  }
               }
            }
         }
         mccCatenateColumns(&RM2_, mccTempMatrix(0, 0., mccINT, 0 ), &RM0_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM5_;
            int I_RM5_=1;
            double *p_xs;
            int I_xs=1, J_xs;
            m_ = mcmCalcResultSize(m_, &n_, ((int)(m - 1) + 1), 1);
            mccAllocateMatrix(&RM5_, m_, n_);
            mccCheckMatrixSize(&xs, (((int)(1 + .5)) + ((int)(m - 1))), 2);
            I_RM5_ = (mccM(&RM5_) != 1 || mccN(&RM5_) != 1);
            p_RM5_ = mccPR(&RM5_);
            if (mccM(&xs) == 1) { I_xs = J_xs = 0;}
            else { I_xs = 1; J_xs=mccM(&xs)-m_; }
            p_xs = mccPR(&xs) + ((int)(1 - .5)) + mccM(&xs) * (2-1);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_xs += J_xs)
               {
                  for (i_=0; i_<m_; ++i_, p_RM5_+=I_RM5_, p_xs+=I_xs)
                  {
                     *p_RM5_ = *p_xs;
                  }
               }
            }
         }
         mccConjTrans(&RM4_, &RM5_);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_RM3_;
            int I_RM3_=1;
            m_ = mcmCalcResultSize(m_, &n_, 1, ((int)((n - 1) - 0) + 1));
            mccAllocateMatrix(&RM3_, m_, n_);
            I_RM3_ = (mccM(&RM3_) != 1 || mccN(&RM3_) != 1);
            p_RM3_ = mccPR(&RM3_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_RM3_+=I_RM3_)
                  {
                     *p_RM3_ = (dl * (double) ((int)(0 + i_+j_*m_)));
                  }
               }
            }
         }
         Mprhs_[0] = &RM2_;
         Mprhs_[1] = &RM4_;
         Mprhs_[2] = &RM3_;
         Mplhs_[0] = &y2;
         mccCallMATLAB(1, Mplhs_, 3, Mprhs_, "spline", 203);
         /* ys(:,1)=y1';ys(:,2)=y2'; */
         mccConjTrans(&RM3_, &M_y1);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_ys;
            int I_ys=1, J_ys;
            double *p_RM3_;
            int I_RM3_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&RM3_), mccN(&RM3_));
            mccGrowMatrix(&ys, m_, 1);
            if (mccM(&ys) == 1) { I_ys = J_ys = 0;}
            else { I_ys = 1; J_ys=mccM(&ys)-m_; }
            p_ys = mccPR(&ys) + 0 + mccM(&ys) * (1-1);
            I_RM3_ = (mccM(&RM3_) != 1 || mccN(&RM3_) != 1);
            p_RM3_ = mccPR(&RM3_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_ys += J_ys)
               {
                  for (i_=0; i_<m_; ++i_, p_ys+=I_ys, p_RM3_+=I_RM3_)
                  {
                     *p_ys = *p_RM3_;
                  }
               }
            }
         }
         mccConjTrans(&RM3_, &y2);
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            double *p_ys;
            int I_ys=1, J_ys;
            double *p_RM3_;
            int I_RM3_=1;
            m_ = mcmCalcResultSize(m_, &n_, mccM(&RM3_), mccN(&RM3_));
            mccGrowMatrix(&ys, m_, 2);
            if (mccM(&ys) == 1) { I_ys = J_ys = 0;}
            else { I_ys = 1; J_ys=mccM(&ys)-m_; }
            p_ys = mccPR(&ys) + 0 + mccM(&ys) * (2-1);
            I_RM3_ = (mccM(&RM3_) != 1 || mccN(&RM3_) != 1);
            p_RM3_ = mccPR(&RM3_);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_, p_ys += J_ys)
               {
                  for (i_=0; i_<m_; ++i_, p_ys+=I_ys, p_RM3_+=I_RM3_)
                  {
                     *p_ys = *p_RM3_;
                  }
               }
            }
         }
         
         /* 'dynamic programming, backwards...' */
         mccPrint (&S17_, 0);
         
         /* %8. INITIALIZE FOR DYNAMIC PROGRAMMING      */
         /* %prepare for  pattern analysis by dynamic programming */
         /* % for state variables given as  */
         /* % 2-vectors (alpha,t); alpha for arctype and t for arclength */
         /* % alpha=1,2...nalpha;t=1,2,..tmax. %jumpc=cost for jumping from arctype to arctype */
         
         /* jumpc=.020; */
         jumpc = .020;
         /* % elastc=cost of elastic change with fixed arctype */
         /* elastc=.0005; */
         elastc = .0005;
         /* [nalpha,tmax]=size(Itemps); */
         if(mccNOTSET(&Itemps))
         {
            mexErrMsgTxt( "variable Itemps undefined, line 217" );
         }
         mccGetMatrixSize(&nalpha,&tmax, &Itemps);
         /* alphas=1:nalpha;ts=1:tmax; */
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_alphas;
            int I_alphas=1;
            m_ = mcmCalcResultSize(m_, &n_, 1, ((int)(nalpha - 1) + 1));
            mccAllocateMatrix(&alphas, m_, n_);
            I_alphas = (mccM(&alphas) != 1 || mccN(&alphas) != 1);
            p_alphas = mccIPR(&alphas);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_alphas+=I_alphas)
                  {
                     *p_alphas = ((int)(1 + i_+j_*m_));
                  }
               }
            }
         }
         {
            int i_, j_;
            int m_=1, n_=1, cx_ = 0;
            int *p_ts;
            int I_ts=1;
            m_ = mcmCalcResultSize(m_, &n_, 1, ((int)(tmax - 1) + 1));
            mccAllocateMatrix(&ts, m_, n_);
            I_ts = (mccM(&ts) != 1 || mccN(&ts) != 1);
            p_ts = mccIPR(&ts);
            if (m_ != 0)
            {
               for (j_=0; j_<n_; ++j_)
               {
                  for (i_=0; i_<m_; ++i_, p_ts+=I_ts)
                  {
                     *p_ts = ((int)(1 + i_+j_*m_));
                  }
               }
            }
         }
         /* M=zeros(nalpha,tmax); */
         mccZerosMN(&M, nalpha, tmax);
         
         /* % result is 2-column matrix analysis,first column alpha values, */
         /* % second columnt t values */
         /* ID=rhos; */
         mccCopy(&ID_1, &rhos);
         /* l=length(ID); */
         l = mccGetLength(&ID_1);
         /* analysis=zeros(l,2); */
         mccZerosMN(&analysis, l, 2);
         /* u=ones(1,nalpha);v=ones(1,tmax); */
         mccOnesMN(&u_1, 1, nalpha);
         mccOnesMN(&v_2, 1, tmax);
         /* % cost starting at t with (alphavalue,tvalue) */
         /* cost=zeros(l,nalpha,tmax); */
         mccError(&S18_);
         
         /* % initialize dynamic programming (backwards) */
         /* for alpha =1:nalpha */
         for (I1_ = 1; I1_ <= nalpha; I1_ = I1_ + 1)
         {
            alpha = I1_;
            /* for t=1:tmax */
            for (I0_ = 1; I0_ <= tmax; I0_ = I0_ + 1)
            {
               t = I0_;
               /* diffs=mod(t-1-ts,templs(alpha)); */
               I2_ = (t - 1);
               {
                  int i_, j_;
                  int m_=1, n_=1, cx_ = 0;
                  int *p_IM2_;
                  int I_IM2_=1;
                  int *p_ts;
                  int I_ts=1;
                  m_ = mcmCalcResultSize(m_, &n_, mccM(&ts), mccN(&ts));
                  mccAllocateMatrix(&IM2_, m_, n_);
                  I_IM2_ = (mccM(&IM2_) != 1 || mccN(&IM2_) != 1);
                  p_IM2_ = mccIPR(&IM2_);
                  I_ts = (mccM(&ts) != 1 || mccN(&ts) != 1);
                  p_ts = mccIPR(&ts);
                  if (m_ != 0)
                  {
                     for (j_=0; j_<n_; ++j_)
                     {
                        for (i_=0; i_<m_; ++i_, p_IM2_+=I_IM2_, p_ts+=I_ts)
                        {
                           *p_IM2_ = (I2_ - ((int)*p_ts));
                        }
                     }
                  }
               }
               if(mccNOTSET(&templs))
               {
                  mexErrMsgTxt( "variable templs undefined, line 233" );
               }
               Mprhs_[0] = &IM2_;
               Mprhs_[1] = mccTempVectorElement(&templs, alpha);
               Mplhs_[0] = &diffs;
               mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "mod", 233);
               
               /* M=(ID(l)-Itemps(alpha,t))^2+jumpc.*kron((alphas~=alpha)',v)+... */
               mccFindIndex(&IM2_, mccTempMatrix(l, 0., mccINT, l_flags_ ));
               C3__r = (mccGetRealMatrixElement(&Itemps, mccRint(alpha), mccRint(t)));
               C3__i = mccGetImagMatrixElement(&Itemps, mccRint(alpha), mccRint(t));
               {
                  int i_, j_;
                  int m_=1, n_=1, cx_ = 0;
                  double *p_CM4_;
                  int I_CM4_=1;
                  double *q_CM4_;
                  double *p_ID_1;
                  double *q_ID_1;
                  int *p_IM2_;
                  int I_IM2_=1;
                  m_ = mccCalcSubscriptDimensions(m_, &n_, mccM(&IM2_), mccN(&IM2_), &ID_1);
                  mccAllocateMatrix(&CM4_, m_, n_);
                  mccCheckVectorSize(&ID_1, mccGetMaxIndex(&IM2_ ,mccM(&ID_1)*mccN(&ID_1)));
                  I_CM4_ = (mccM(&CM4_) != 1 || mccN(&CM4_) != 1);
                  p_CM4_ = mccPR(&CM4_);
                  q_CM4_ = mccPI(&CM4_);
                  I_IM2_ = (mccM(&IM2_) != 1 || mccN(&IM2_) != 1);
                  p_IM2_ = mccIPR(&IM2_);
                  if (m_ != 0)
                  {
                     for (j_=0; j_<n_; ++j_)
                     {
                        for (i_=0; i_<m_; ++i_, p_CM4_+=I_CM4_, q_CM4_+=I_CM4_, p_IM2_+=I_IM2_)
                        {
                           *p_CM4_ = (mccPR(&ID_1)[((int)(*p_IM2_ - .5))] - C3__r);
                           *q_CM4_ = (mccPI(&ID_1)[((int)(*p_IM2_ - .5))] - C3__i);
                        }
                     }
                  }
               }
               mccPower(&CM3_, &CM4_, mccTempMatrix(2, 0., mccINT, 0 ));
               {
                  int i_, j_;
                  int m_=1, n_=1, cx_ = 0;
                  unsigned short *p_BM0_;
                  int I_BM0_=1;
                  int *p_alphas;
                  int I_alphas=1;
                  m_ = mcmCalcResultSize(m_, &n_, mccM(&alphas), mccN(&alphas));
                  mccAllocateMatrix(&BM0_, m_, n_);
                  I_BM0_ = (mccM(&BM0_) != 1 || mccN(&BM0_) != 1);
                  p_BM0_ = mccSPR(&BM0_);
                  I_alphas = (mccM(&alphas) != 1 || mccN(&alphas) != 1);
                  p_alphas = mccIPR(&alphas);
                  if (m_ != 0)
                  {
                     for (j_=0; j_<n_; ++j_)
                     {
                        for (i_=0; i_<m_; ++i_, p_BM0_+=I_BM0_, p_alphas+=I_alphas)
                        {
                           *p_BM0_ = (((int)*p_alphas) != alpha);
                        }
                     }
                  }
               }
               mccConjTrans(&BM1_, &BM0_);
               Mprhs_[0] = &BM1_;
               Mprhs_[1] = &v_2;
               Mplhs_[0] = &IM0_;
               mccCallMATLAB(1, Mplhs_, 2, Mprhs_, "kron", 235);
               {
                  int i_, j_;
                  int m_=1, n_=1, cx_ = 0;
                  unsigned short *p_BM2_;
                  int I_BM2_=1;
                  int *p_alphas;
                  int I_alphas=1;
                  m_ = mcmCalcResultSize(m_, &n_, mccM(&alphas), mccN(&alphas));
                  mccAllocateMatrix(&BM2_, m_, n_);
                  I_BM2_ = (mccM(&BM2_) != 1 || mccN(&BM2_) != 1);
                  p_BM2_ = mccSPR(&BM2_);
                  I_alphas = (mccM(&alphas) != 1 || mccN(&alphas) != 1);
                  p_alphas = mccIPR(&alphas);
                  if (m_ != 0)
                  {
                     for (j_=0; j_<n_; ++j_)
                     {
                        for (i_=0; i_<m_; ++i_, p_BM2_+=I_BM2_, p_alphas+=I_alphas)
                        {
                           *p_BM2_ = (((int)*p_alphas) == alpha);
                        }
                     }
                  }
               }
               mccConjTrans(&BM3_, &BM2_);
               C2__r = ((mccGetRealVectorElement(&templs, mccRint(alpha))) - 1);
               C2__i = (mccGetImagVectorElement(&templs, mccRint(alpha)) - 0.);
               {
                  int i_, j_;
                  int m_=1, n_=1, cx_ = 0;
                  double *p_CM5_;
                  int I_CM5_=1;
                  double *q_CM5_;
                  double *p_diffs;
                  int I_diffs=1;
                  double *p_1diffs;
                  int I_1diffs=1;
                  m_ = mcmCalcResultSize(m_, &n_, mccM(&diffs), mccN(&diffs));
                  m_ = mcmCalcResultSize(m_, &n_, mccM(&diffs), mccN(&diffs));
                  mccAllocateMatrix(&CM5_, m_, n_);
                  I_CM5_ = (mccM(&CM5_) != 1 || mccN(&CM5_) != 1);
                  p_CM5_ = mccPR(&CM5_);
                  q_CM5_ = mccPI(&CM5_);
                  I_diffs = (mccM(&diffs) != 1 || mccN(&diffs) != 1);
                  p_diffs = mccPR(&diffs);