00001 #include "SUMA_suma.h"
00002
00003 SUMA_SurfaceViewer *SUMAg_cSV = NULL;
00004 SUMA_SurfaceViewer *SUMAg_SVv = NULL;
00005
00006 int SUMAg_N_SVv = 0;
00007 SUMA_DO *SUMAg_DOv = NULL;
00008 int SUMAg_N_DOv = 0;
00009 SUMA_CommonFields *SUMAg_CF = NULL;
00010
00011 void usage_SurfToSurf (SUMA_GENERIC_ARGV_PARSE *ps)
00012 {
00013 static char FuncName[]={"usage_SurfToSurf"};
00014 char * s = NULL, *sio=NULL, *st = NULL, *sts = NULL;
00015 int i;
00016 s = SUMA_help_basics();
00017 sio = SUMA_help_IO_Args(ps);
00018 printf ( "\n"
00019 "Usage: SurfToSurf <-i_TYPE S1> [<-sv SV1>]\n"
00020 " <-i_TYPE S2> [<-sv SV1>]\n"
00021 " [<-prefix PREFIX>]\n"
00022 " [<-output_params PARAM_LIST>]\n"
00023 " [<-node_indices NODE_INDICES>]\n"
00024 " [<-proj_dir PROJ_DIR>]\n"
00025 " [<-data DATA>]\n"
00026 " [<-node_debug NODE>]\n"
00027 " [<-debug DBG_LEVEL>]\n"
00028 " [-make_consistent]\n"
00029 " \n"
00030 " This program is used to interpolate data from one surface (S2)\n"
00031 " to another (S1), assuming the surfaces are quite similar in\n"
00032 " shape but having different meshes (non-isotopic).\n"
00033 " This is done by projecting each node (nj) of S1 along the normal\n"
00034 " at nj and finding the closest triangle t of S2 that is intersected\n"
00035 " by this projection. Projection is actually bidirectional.\n"
00036 " If such a triangle t is found, the nodes (of S2) forming it are \n"
00037 " considered to be the neighbors of nj.\n"
00038 " Values (arbitrary data, or coordinates) at these neighboring nodes\n"
00039 " are then transfered to nj using barycentric interpolation or \n"
00040 " nearest-node interpolation.\n"
00041 " Nodes whose projections fail to intersect triangles in S2 are given\n"
00042 " nonsensical values of -1 and 0.0 in the output.\n"
00043 "\n"
00044 " Mandatory input:\n"
00045 " Two surfaces are required at input. See -i_TYPE options\n"
00046 " below for more information. \n"
00047 "\n"
00048 " Optional input:\n"
00049 " -prefix PREFIX: Specify the prefix of the output file.\n"
00050 " The output file is in 1D format at the moment.\n"
00051 " Default is SurfToSurf\n"
00052 " -output_params PARAM_LIST: Specify the list of mapping\n"
00053 " parameters to include in output\n"
00054 " PARAM_LIST can have any or all of the following:\n"
00055 " NearestTriangleNodes: Use Barycentric interpolation (default)\n"
00056 " and output indices of 3 nodes from S2\n"
00057 " that neighbor nj of S1\n"
00058 " NearestNode: Use only the closest node from S2 (of the three \n"
00059 " closest neighbors) to nj of S1 for interpolation\n"
00060 " and output the index of that closest node.\n"
00061 " NearestTriangle: Output index of triangle t from S2 that\n"
00062 " is the closest to nj along its projection\n"
00063 " direction. \n"
00064 " DistanceToSurf: Output distance (signed) from nj, along \n"
00065 " projection direction to S2.\n"
00066 " This is the parameter output by the precursor\n"
00067 " program CompareSurfaces\n"
00068 " ProjectionOnSurf: Output coordinates of projection of nj onto \n"
00069 " triangle t of S2.\n"
00070 " Data: Output the data from S2, interpolated onto S1\n"
00071 " If no data is specified via the -data option, then\n"
00072 " the XYZ coordinates of SO2's nodes are considered\n"
00073 " the data.\n"
00074 " -data DATA: 1D file containing data to be interpolated.\n"
00075 " Each row i contains data for node i of S2.\n"
00076 " You must have one row for each node making up S2.\n"
00077 " In other terms, if S2 has N nodes, you need N rows\n"
00078 " in DATA. \n"
00079 " Each column of DATA is processed separately (think\n"
00080 " sub-bricks, and spatial interpolation).\n"
00081 " You can use [] selectors to choose a subset \n"
00082 " of columns.\n"
00083 " If -data option is not specified and Data is in PARAM_LIST\n"
00084 " then the XYZ coordinates of SO2's nodes are the data.\n"
00085 " -node_indices NODE_INDICES: 1D file containing the indices of S1\n"
00086 " to consider. The default is all of the\n"
00087 " nodes in S1. Only one column of values is\n"
00088 " allowed here, use [] selectors to choose\n"
00089 " the column of node indices if NODE_INDICES\n"
00090 " has multiple columns in it.\n"
00091 " -proj_dir PROJ_DIR: 1D file containing projection directions to use\n"
00092 " instead of the node normals of S1.\n"
00093 " Each row should contain one direction for each\n"
00094 " of the nodes forming S1.\n"
00095 " -make_consistent: Force a consistency check and correct triangle \n"
00096 " orientation of S1 if needed. Triangles are also\n"
00097 " oriented such that the majority of normals point\n"
00098 " away from center of surface.\n"
00099 " The program might not succeed in repairing some\n"
00100 " meshes with inconsistent orientation.\n"
00101 "\n"
00102 "%s"
00103 "%s"
00104 "\n", sio, s);
00105 SUMA_free(s); s = NULL; SUMA_free(st); st = NULL; SUMA_free(sio); sio = NULL;
00106 s = SUMA_New_Additions(0, 1); printf("%s\n", s);SUMA_free(s); s = NULL;
00107 printf(" Ziad S. Saad SSCC/NIMH/NIH ziad@nih.gov \n"
00108 " Shruti Japee LBC/NIMH/NIH shruti@codon.nih.gov \n");
00109 exit(0);
00110 }
00111
00112 SUMA_GENERIC_PROG_OPTIONS_STRUCT *SUMA_SurfToSurf_ParseInput(char *argv[], int argc, SUMA_GENERIC_ARGV_PARSE *ps)
00113 {
00114 static char FuncName[]={"SUMA_BrainWrap_ParseInput"};
00115 SUMA_GENERIC_PROG_OPTIONS_STRUCT *Opt=NULL;
00116 int kar;
00117 SUMA_Boolean brk, accepting_out;
00118
00119 SUMA_Boolean LocalHead = NOPE;
00120
00121 SUMA_ENTRY;
00122
00123 Opt = SUMA_Alloc_Generic_Prog_Options_Struct();
00124 kar = 1;
00125 brk = NOPE;
00126 Opt->in_1D = NULL;
00127 Opt->NodeDbg = -1;
00128 Opt->debug = 0;
00129 Opt->NearestNode = 3;
00130 Opt->NearestTriangle = 0;
00131 Opt->DistanceToMesh = 0;
00132 Opt->ProjectionOnMesh = 0;
00133 Opt->Data = 0;
00134 Opt->in_name = NULL;
00135 Opt->out_prefix = NULL;
00136 Opt->fix_winding = 0;
00137 accepting_out = NOPE;
00138 while (kar < argc) {
00139
00140
00141 if (!brk && accepting_out) {
00142 if (*(argv[kar]) == '-') accepting_out = NOPE;
00143 }
00144
00145 if (strcmp(argv[kar], "-h") == 0 || strcmp(argv[kar], "-help") == 0) {
00146 usage_SurfToSurf(ps);
00147 exit (0);
00148 }
00149
00150 SUMA_SKIP_COMMON_OPTIONS(brk, kar);
00151
00152 if (!brk && (strcmp(argv[kar], "-debug") == 0))
00153 {
00154 if (kar+1 >= argc)
00155 {
00156 fprintf (SUMA_STDERR, "need a number after -debug \n");
00157 exit (1);
00158 }
00159
00160 Opt->debug = atoi(argv[++kar]);
00161 brk = YUP;
00162 }
00163
00164 if (!brk && (strcmp(argv[kar], "-node_debug") == 0))
00165 {
00166 if (kar+1 >= argc)
00167 {
00168 fprintf (SUMA_STDERR, "need a number after -node_debug \n");
00169 exit (1);
00170 }
00171
00172 Opt->NodeDbg = atoi(argv[++kar]);
00173 brk = YUP;
00174 }
00175
00176 if (!brk && (strcmp(argv[kar], "-node_indices") == 0))
00177 {
00178 if (kar+1 >= argc)
00179 {
00180 fprintf (SUMA_STDERR, "need a parameter after -node_indices \n");
00181 exit (1);
00182 }
00183
00184 Opt->in_nodeindices = argv[++kar];
00185 brk = YUP;
00186 }
00187
00188 if (!brk && (strcmp(argv[kar], "-output_params") == 0))
00189 {
00190 if (kar+1 >= argc)
00191 {
00192 fprintf (SUMA_STDERR, "need at least one parameter after output_params \n");
00193 exit (1);
00194 }
00195
00196 accepting_out = YUP;
00197 brk = YUP;
00198 }
00199
00200 if (!brk && accepting_out && (strcmp(argv[kar], "NearestNode") == 0)) {
00201 if (Opt->NearestNode < 1) Opt->NearestNode = 1;
00202 brk = YUP;
00203 }
00204
00205 if (!brk && accepting_out && (strcmp(argv[kar], "NearestTriangleNodes") == 0)) {
00206 if (Opt->NearestNode < 3) Opt->NearestNode = 3;
00207 brk = YUP;
00208 }
00209
00210 if (!brk && accepting_out && (strcmp(argv[kar], "NearestTriangle") == 0)) {
00211 if (Opt->NearestTriangle < 1) Opt->NearestTriangle = 1;
00212 brk = YUP;
00213 }
00214
00215 if (!brk && accepting_out && (strcmp(argv[kar], "DistanceToSurf") == 0)) {
00216 Opt->DistanceToMesh = 1;
00217 brk = YUP;
00218 }
00219
00220 if (!brk && accepting_out && (strcmp(argv[kar], "ProjectionOnSurf") == 0)) {
00221 Opt->ProjectionOnMesh = 1;
00222 brk = YUP;
00223 }
00224
00225 if (!brk && accepting_out && (strcmp(argv[kar], "Data") == 0)) {
00226 Opt->Data = 1;
00227 brk = YUP;
00228 }
00229
00230 if (!brk && (strcmp(argv[kar], "-data") == 0))
00231 {
00232 if (kar+1 >= argc)
00233 {
00234 fprintf (SUMA_STDERR, "need a name after -data \n");
00235 exit (1);
00236 }
00237 ++kar;
00238 if (strcmp(argv[kar],"_XYZ_") == 0) {
00239
00240 if (Opt->in_name) {
00241 SUMA_SL_Err("Unexpected non null value");
00242 exit (1);
00243 }
00244 } else {
00245 Opt->in_name = SUMA_copy_string(argv[kar]);
00246 }
00247 Opt->Data = 1;
00248 brk = YUP;
00249 }
00250
00251 if (!brk && (strcmp(argv[kar], "-prefix") == 0))
00252 {
00253 if (kar+1 >= argc)
00254 {
00255 fprintf (SUMA_STDERR, "need a name after -prefix \n");
00256 exit (1);
00257 }
00258
00259 Opt->out_prefix = SUMA_Extension(argv[++kar],".1D", YUP);
00260 brk = YUP;
00261 }
00262
00263 if (!brk && (strcmp(argv[kar], "-proj_dir") == 0))
00264 {
00265 if (kar+1 >= argc)
00266 {
00267 fprintf (SUMA_STDERR, "need a name after -proj_dir \n");
00268 exit (1);
00269 }
00270
00271 Opt->in_1D = argv[++kar];
00272 brk = YUP;
00273 }
00274 if (!brk && (strcmp(argv[kar], "-make_consistent") == 0))
00275 {
00276 Opt->fix_winding = 1;
00277 brk = YUP;
00278 }
00279 if (!brk && !ps->arg_checked[kar]) {
00280 fprintf (SUMA_STDERR,"Error %s:\nOption %s not understood. Try -help for usage\n", FuncName, argv[kar]);
00281 exit (1);
00282 } else {
00283 brk = NOPE;
00284 kar ++;
00285 }
00286 }
00287
00288 if (!Opt->out_prefix) Opt->out_prefix = SUMA_copy_string("SurfToSurf");
00289
00290 SUMA_RETURN(Opt);
00291 }
00292
00293 int main (int argc,char *argv[])
00294 {
00295 static char FuncName[]={"SurfToSurf"};
00296 SUMA_GENERIC_PROG_OPTIONS_STRUCT *Opt;
00297 SUMA_GENERIC_ARGV_PARSE *ps=NULL;
00298 SUMA_SurfaceObject *SO1=NULL, *SO2 = NULL;
00299 SUMA_SurfSpecFile *Spec = NULL;
00300 SUMA_M2M_STRUCT *M2M = NULL;
00301 int N_Spec=0, *nodeind = NULL, N_nodeind, icol, i, j;
00302 MRI_IMAGE *im = NULL, *im_data=NULL;
00303 int nvec, ncol=0, nvec_data, ncol_data, Nchar;
00304 float *far = NULL, *far_data=NULL, *dt = NULL, *projdir=NULL;
00305 char *outname = NULL, *s=NULL, sbuf[100];
00306 void *SO_name = NULL;
00307 FILE *outptr=NULL;
00308 SUMA_Boolean exists = NOPE;
00309 SUMA_INDEXING_ORDER d_order;
00310 SUMA_STRING *SS=NULL;
00311 SUMA_Boolean LocalHead = NOPE;
00312
00313 SUMA_mainENTRY;
00314 SUMA_STANDALONE_INIT;
00315
00316
00317 SUMAg_DOv = SUMA_Alloc_DisplayObject_Struct (SUMA_MAX_DISPLAYABLE_OBJECTS);
00318 ps = SUMA_Parse_IO_Args(argc, argv, "-i;-t;-spec;-s;-sv;-o;");
00319
00320 if (argc < 2) {
00321 usage_SurfToSurf(ps);
00322 exit (1);
00323 }
00324
00325 Opt = SUMA_SurfToSurf_ParseInput (argv, argc, ps);
00326
00327
00328 if (ps->o_N_surfnames) {
00329 SO_name = SUMA_Prefix2SurfaceName(ps->o_surfnames[0], NULL, NULL, ps->o_FT[0], &exists);
00330 if (exists) {
00331 fprintf(SUMA_STDERR,"Error %s:\nOutput file(s) %s* on disk.\nWill not overwrite.\n", FuncName, ps->o_surfnames[0]);
00332 exit(1);
00333 }
00334 }
00335
00336 if (Opt->debug > 2) LocalHead = YUP;
00337
00338 outname = SUMA_append_string(Opt->out_prefix,".1D");
00339 if (SUMA_filexists(outname)) {
00340 fprintf(SUMA_STDERR,"Output file %s exists.\n", outname);
00341 exit(1);
00342 }
00343
00344
00345 Spec = SUMA_IO_args_2_spec(ps, &N_Spec);
00346 if (N_Spec != 1) {
00347 SUMA_S_Err( "Multiple spec at input.\n"
00348 "Do not mix surface input types together\n");
00349 exit(1);
00350 }
00351
00352 if (Spec->N_Surfs != 2) {
00353 SUMA_S_Err("2 surfaces expected.");
00354 exit(1);
00355 }
00356
00357 SO1 = SUMA_Load_Spec_Surf(Spec, 0, ps->sv[0], 0);
00358 if (!SO1) {
00359 fprintf (SUMA_STDERR,"Error %s:\n"
00360 "Failed to find surface\n"
00361 "in spec file. \n",
00362 FuncName );
00363 exit(1);
00364 }
00365 if (!SUMA_SurfaceMetrics(SO1, "EdgeList|MemberFace", NULL)) { SUMA_SL_Err("Failed to create edge list for SO1"); exit(1); }
00366 if (Opt->fix_winding) {
00367 int orient, trouble;
00368 if (LocalHead) fprintf(SUMA_STDERR,"%s: Making sure S1 is consistently orientated\n", FuncName);
00369 if (!SUMA_MakeConsistent (SO1->FaceSetList, SO1->N_FaceSet, SO1->EL, Opt->debug, &trouble)) {
00370 SUMA_SL_Err("Failed in SUMA_MakeConsistent");
00371 }
00372 if (trouble && LocalHead) {
00373 fprintf(SUMA_STDERR,"%s: trouble value of %d from SUMA_MakeConsistent.\n"
00374 "Inconsistencies were found and corrected unless \n"
00375 "stderr output messages from SUMA_MakeConsistent\n"
00376 "indicate otherwise.\n", FuncName, trouble);
00377 }
00378 if (LocalHead) fprintf(SUMA_STDERR,"%s: Checking orientation.\n", FuncName);
00379 orient = SUMA_OrientTriangles (SO1->NodeList, SO1->N_Node, SO1->FaceSetList, SO1->N_FaceSet, 1, 0);
00380 if (orient < 0) {
00381
00382 if (SO1->EL) SUMA_free_Edge_List(SO1->EL); SO1->EL = NULL;
00383 if (!SUMA_SurfaceMetrics(SO1, "EdgeList", NULL)) { SUMA_SL_Err("Failed to create edge list for SO1"); exit(1); }
00384
00385 if (SO1->NodeNormList) SUMA_free(SO1->NodeNormList); SO1->NodeNormList = NULL;
00386 if (SO1->FaceNormList) SUMA_free(SO1->FaceNormList); SO1->FaceNormList = NULL;
00387 }
00388 if (!orient) { fprintf(SUMA_STDERR,"Error %s:\nFailed in SUMA_OrientTriangles\n", FuncName); }
00389 if (LocalHead) {
00390 if (orient < 0) { SUMA_SL_Note("S1 was reoriented"); }
00391 else { SUMA_SL_Note("S1 was properly oriented"); }
00392 }
00393 }
00394
00395
00396 if (!SO1->NodeNormList || !SO1->FaceNormList) { SUMA_LH("Node Normals"); SUMA_RECOMPUTE_NORMALS(SO1); }
00397 if (Opt->NodeDbg >= SO1->N_Node) {
00398 SUMA_SL_Warn("node_debug index is larger than number of nodes in surface, ignoring -node_debug.");
00399 Opt->NodeDbg = -1;
00400 }
00401
00402 SO2 = SUMA_Load_Spec_Surf(Spec, 1, ps->sv[1], 0);
00403 if (!SO2) {
00404 fprintf (SUMA_STDERR,"Error %s:\n"
00405 "Failed to find surface\n"
00406 "in spec file. \n",
00407 FuncName );
00408 exit(1);
00409 }
00410 if (!SUMA_SurfaceMetrics(SO2, "EdgeList|MemberFace", NULL)) { SUMA_SL_Err("Failed to create edge list for SO2"); exit(1); }
00411 if (!SO2->NodeNormList || !SO2->FaceNormList) { SUMA_LH("Node Normals"); SUMA_RECOMPUTE_NORMALS(SO2); }
00412
00413 if (LocalHead) {
00414 SUMA_LH("Surf1");
00415 SUMA_Print_Surface_Object(SO1, NULL);
00416 SUMA_LH("Surf2");
00417 SUMA_Print_Surface_Object(SO2, NULL);
00418 }
00419
00420
00421 nodeind = NULL; N_nodeind = 0;
00422 if (Opt->in_nodeindices) {
00423 im = mri_read_1D(Opt->in_nodeindices);
00424 if (!im) { SUMA_SL_Err("Failed to read 1D file of node indices"); exit(1);}
00425 far = MRI_FLOAT_PTR(im);
00426 N_nodeind = nvec = im->nx;
00427 ncol = im->ny;
00428 if (ncol != 1) { SUMA_SL_Err("More than one column in node index input file."); exit(1);}
00429 nodeind = (int *)SUMA_calloc(nvec, sizeof(int));
00430 if (!nodeind) { SUMA_SL_Crit("Failed to allocate"); exit(1); }
00431 for (i=0;i<nvec;++i) {
00432 nodeind[i] = (int)far[i];
00433 if (nodeind[i] < 0 || nodeind[i] >= SO1->N_Node) {
00434 fprintf(SUMA_STDERR, "Error %s: A node index of %d was found in input file %s, entry %d.\n"
00435 "Acceptable indices are positive and less than %d\n",
00436 FuncName, nodeind[i], Opt->in_nodeindices, i, SO1->N_Node);
00437 exit(1);
00438 }
00439 }
00440 mri_free(im); im = NULL;
00441 }
00442
00443
00444 projdir = NULL;
00445 if (Opt->in_1D) {
00446 im = mri_read_1D(Opt->in_1D);
00447 if (!im) { SUMA_SL_Err("Failed to read 1D file of projection directions"); exit(1);}
00448 far = MRI_FLOAT_PTR(im);
00449 if (im->ny != 3) { SUMA_SL_Err("Need three columns in projection directions file."); exit(1); }
00450 if (im->nx != SO1->N_Node) {
00451 fprintf(SUMA_STDERR, "Error %s: You must have a direction for each node in SO1.\n"
00452 "%d directions found but SO1 has %d nodes.\n", FuncName, im->nx, SO1->N_Node);
00453 exit(1);
00454 }
00455
00456
00457 projdir = (float *)SUMA_calloc(SO1->N_Node*3, sizeof(float));
00458 if (!projdir) { SUMA_SL_Crit("Failed to allocate"); exit(1); }
00459 for (i=0; i<SO1->N_Node; ++i) {
00460 projdir[3*i ] = far[i ];
00461 projdir[3*i+1] = far[i+ SO1->N_Node];
00462 projdir[3*i+2] = far[i+2*SO1->N_Node];
00463 }
00464 mri_free(im); im = NULL;
00465
00466 }
00467
00468 if (SO_name) {
00469
00470 if (nodeind) {
00471 fprintf(SUMA_STDERR, "Error %s: You cannot combine option -o_TYPE with -node_indices", FuncName);
00472 exit(1);
00473 }
00474 if (Opt->in_name) {
00475 fprintf(SUMA_STDERR, "Error %s: You cannot combine option -o_TYPE with -data", FuncName);
00476 exit(1);
00477 }
00478 }
00479
00480 if (Opt->in_name) {
00481 im_data = mri_read_1D(Opt->in_name);
00482 if (!im_data) { SUMA_SL_Err("Failed to read 1D file of data"); exit(1);}
00483 far_data = MRI_FLOAT_PTR(im_data);
00484 nvec_data = im_data->nx;
00485 ncol_data = im_data->ny;
00486 if (nvec_data != SO2->N_Node) {
00487 SUMA_SL_Err("Your data file must have one row for each node in surface 2.\n"); exit(1);
00488 }
00489 d_order = SUMA_COLUMN_MAJOR;
00490 } else {
00491 im_data = NULL;
00492 far_data = SO2->NodeList;
00493 nvec_data = SO2->N_Node;
00494 ncol_data = 3;
00495 d_order = SUMA_ROW_MAJOR;
00496 }
00497
00498
00499 SUMA_LH("Going for the mapping of SO1 --> SO2");
00500 M2M = SUMA_GetM2M_NN( SO1, SO2, nodeind, N_nodeind, projdir, 0, Opt->NodeDbg);
00501
00502
00503 if (Opt->NodeDbg >= 0) {
00504 char *s = NULL;
00505 s = SUMA_M2M_node_Info(M2M, Opt->NodeDbg);
00506 fprintf(SUMA_STDERR,"%s: Debug for node %d ([%f, %f, %f])of SO1:\n%s\n\n",
00507 FuncName, Opt->NodeDbg,
00508 SO1->NodeList[3*Opt->NodeDbg], SO1->NodeList[3*Opt->NodeDbg+1], SO1->NodeList[3*Opt->NodeDbg+2],
00509 s);
00510 SUMA_free(s); s = NULL;
00511 }
00512
00513
00514 if (Opt->Data) {
00515 if (Opt->NearestNode > 1) dt = SUMA_M2M_interpolate(M2M, far_data, ncol_data, nvec_data, d_order, 0 );
00516 else if (Opt->NearestNode == 1) dt = SUMA_M2M_interpolate(M2M, far_data, ncol_data, nvec_data, d_order, 1 );
00517 if (!dt) {
00518 SUMA_SL_Err("Failed to interpolate");
00519 exit(1);
00520 }
00521 }
00522
00523 SUMA_LH("Forming the output");
00524 outptr = fopen(outname,"w");
00525 if (!outptr) {
00526 SUMA_SL_Err("Failed to open file for output.\n");
00527 exit(1);
00528 }
00529
00530
00531 SS = SUMA_StringAppend(NULL, NULL);
00532 SS = SUMA_StringAppend_va(SS, "#Mapping from nodes on surf 1 (S1) to nodes on surf 2 (S2)\n"
00533 "# Surf 1 is labeled %s, idcode:%s\n"
00534 "# Surf 2 is labeled %s, idcode:%s\n",
00535 SO1->Label, SO1->idcode_str, SO2->Label, SO2->idcode_str);
00536 icol = 0;
00537 SS = SUMA_StringAppend_va(SS, "#Col. %d:\n"
00538 "# S1n (or nj): Index of node on S1\n"
00539 , icol); ++icol;
00540 if (Opt->NearestNode > 1) {
00541 SS = SUMA_StringAppend_va(SS,
00542 "#Col. %d..%d:\n"
00543 "# S2ne_S1n: Indices of %d nodes on S2 \n"
00544 "# that are closest neighbors of nj.\n"
00545 "# The first index is that of the node on S2 that is closest to nj.\n"
00546 "# If -1 then thes values should be ignored. This happens when nj's projection failed.\n"
00547 , icol, icol+Opt->NearestNode-1, Opt->NearestNode); icol += Opt->NearestNode;
00548 SS = SUMA_StringAppend_va(SS,
00549 "#Col. %d..%d:\n"
00550 "# S2we_S1n: Weights assigned to nodes on surf 2 (S2) \n"
00551 "# that are closest neighbors of nj.\n"
00552 , icol, icol+Opt->NearestNode-1, Opt->NearestNode); icol += Opt->NearestNode;
00553 } else if (Opt->NearestNode == 1) {
00554 SS = SUMA_StringAppend_va(SS,
00555 "#Col. %d:\n"
00556 "# S2ne_S1n: Index of the node on S2 (label:%s idcode:%s)\n"
00557 "# that is the closest neighbor of nj.\n"
00558 "# If -1 then this value should be ignored. This happens when nj's projection failed.\n"
00559 , icol, SO2->Label, SO2->idcode_str); ++icol;
00560 }
00561 if (Opt->NearestTriangle) {
00562 SS = SUMA_StringAppend_va(SS,
00563 "#Col. %d:\n"
00564 "# S2t_S1n: Index of the triangle on S2 that hosts node nj on S1.\n"
00565 "# In other words, nj's closest projection onto S2 falls on \n"
00566 "# triangle S2t_S1n\n"
00567 "# If -1 then this value should be ignored. This happens when nj's projection failed.\n"
00568 , icol); ++icol;
00569 }
00570 if (Opt->ProjectionOnMesh) {
00571 SS = SUMA_StringAppend_va(SS,
00572 "#Col. %d..%d:\n"
00573 "# S2p_S1n: Coordinates of projection of nj onto S2\n"
00574 , icol, icol+2); icol += 3;
00575 }
00576 if (Opt->DistanceToMesh) {
00577 SS = SUMA_StringAppend_va(SS,
00578 "#Col. %d:\n"
00579 "# Closest distance from nj to S2\n"
00580 , icol); ++icol;
00581 }
00582 if (Opt->Data) {
00583 if (!Opt->in_name) {
00584 SS = SUMA_StringAppend_va(SS,
00585 "#Col. %d..%d:\n"
00586 "# Interpolation using XYZ coordinates of nodes on S2 that neighbor nj\n"
00587 "# (same as coordinates of node's projection onto triangle in S2, if using \n"
00588 " barycentric interpolation)\n"
00589 , icol, icol+2); icol += 3;
00590 } else {
00591 SS = SUMA_StringAppend_va(SS,
00592 "#Col. %d..%d:\n"
00593 "# Interpolation of data at nodes on S2 that neighbor nj\n"
00594 "# Data obtained from %s\n"
00595 , icol, icol+ncol_data-1, Opt->in_name); icol += ncol_data;
00596 }
00597 }
00598 s = SUMA_HistString("SurfToSurf", argc, argv, NULL);
00599 SS = SUMA_StringAppend_va(SS,
00600 "#History:\n"
00601 "#%s\n", s); SUMA_free(s); s = NULL;
00602 SUMA_SS2S(SS,s);
00603 fprintf(outptr,"%s\n",s); SUMA_free(s); s = NULL;
00604
00605
00606 Nchar = 6;
00607 for (i=0; i<icol; ++i) { sprintf(sbuf,"#%s", MV_format_fval2(i, Nchar -1)); fprintf(outptr,"%6s ", sbuf); }
00608 fprintf(outptr,"\n");
00609
00610
00611 for (i=0; i<M2M->M1Nn; ++i) {
00612 fprintf(outptr,"%6s ", MV_format_fval2(M2M->M1n[i], Nchar));
00613 if (Opt->NearestNode > 0) {
00614 for (j=0; j<Opt->NearestNode; ++j) {
00615 if (j < M2M->M2Nne_M1n[i]) fprintf(outptr,"%6s ", MV_format_fval2(M2M->M2ne_M1n[i][j], Nchar));
00616 else fprintf(outptr,"%6s ", "-1");
00617 }
00618 }
00619 if (Opt->NearestNode > 1) {
00620 for (j=0; j<Opt->NearestNode; ++j) {
00621 if (j < M2M->M2Nne_M1n[i]) fprintf(outptr,"%6s ", MV_format_fval2(M2M->M2we_M1n[i][j], Nchar));
00622 else fprintf(outptr,"%6s ", "0.0");
00623 }
00624 }
00625 if (Opt->NearestTriangle) {
00626 fprintf(outptr,"%6s ", MV_format_fval2(M2M->M2t_M1n[i], Nchar));
00627 }
00628 if (Opt->ProjectionOnMesh) {
00629 fprintf(outptr,"%6s ", MV_format_fval2(M2M->M2p_M1n[3*i], Nchar));
00630 fprintf(outptr,"%6s ", MV_format_fval2(M2M->M2p_M1n[3*i+1], Nchar));
00631 fprintf(outptr,"%6s ", MV_format_fval2(M2M->M2p_M1n[3*i+2], Nchar));
00632 }
00633 if (Opt->DistanceToMesh) {
00634 fprintf(outptr,"%6s ", MV_format_fval2(M2M->PD[i], Nchar));
00635 }
00636 if (Opt->Data) {
00637 if (!Opt->in_name) {
00638 fprintf(outptr,"%6s ", MV_format_fval2(dt[3*i], Nchar));
00639 fprintf(outptr,"%6s ", MV_format_fval2(dt[3*i+1], Nchar));
00640 fprintf(outptr,"%6s ", MV_format_fval2(dt[3*i+2], Nchar));
00641 } else {
00642 for (j=0; j<ncol_data; ++j) { fprintf(outptr,"%6s ", MV_format_fval2(dt[i+j*M2M->M1Nn], Nchar)); }
00643 }
00644 }
00645 fprintf(outptr,"\n");
00646 }
00647
00648
00649 if (SO_name) {
00650 float *tmpfv = NULL;
00651 SUMA_LH("Writing surface");
00652 tmpfv = SO1->NodeList;
00653 SO1->NodeList = dt;
00654 if (!SUMA_Save_Surface_Object (SO_name, SO1, ps->o_FT[0], ps->o_FF[0], NULL)) {
00655 fprintf (SUMA_STDERR,"Error %s: Failed to write surface object.\n", FuncName);
00656 exit (1);
00657 }
00658 SO1->NodeList = tmpfv; tmpfv = NULL;
00659 }
00660
00661 if (projdir) SUMA_free(projdir); projdir = NULL;
00662 if (SO_name) SUMA_free(SO_name); SO_name = NULL;
00663 if (outptr) fclose(outptr); outptr = NULL;
00664 if (dt) SUMA_free(dt); dt = NULL;
00665 if (s) SUMA_free(s); s = NULL;
00666 if (im_data) mri_free(im_data); im_data = NULL;
00667 if (nodeind) SUMA_free(nodeind); nodeind = NULL;
00668 if (M2M) M2M = SUMA_FreeM2M(M2M);
00669 if (SO1) SUMA_Free_Surface_Object(SO1); SO1 = NULL;
00670 if (SO2) SUMA_Free_Surface_Object(SO2); SO2 = NULL;
00671 if (Spec) SUMA_free(Spec); Spec = NULL;
00672 if (ps) SUMA_FreeGenericArgParse(ps); ps = NULL;
00673 if (Opt) Opt = SUMA_Free_Generic_Prog_Options_Struct(Opt);
00674 if (!SUMA_Free_CommonFields(SUMAg_CF)) SUMA_error_message(FuncName,"SUMAg_CF Cleanup Failed!",1);
00675 exit(0);
00676
00677 }