1/*******************************************************************************
2*
3* McStas, neutron ray-tracing package
4*         Copyright (C) 1997-2009, All rights reserved
5*         Risoe National Laboratory, Roskilde, Denmark
6*         Institut Laue Langevin, Grenoble, France
7*
8* Instrument: mcformat, McStas format converter
9*
10* %Identification
11* Written by: <a href="mailto:farhi@ill.fr">Emmanuel Farhi</a>
12* Date: 1st Feb 2001.
13* Origin: <a href="http://www.ill.fr">ILL (France)</a>
14* Release: McXtrace 1.2 - Jul. 02, 2015
15* Version: $Revision$
16*
17* A McCode format converter to merge/convert data files.
18*
19* %Description
20*
21* A McCode format converter to merge/convert data files.
22* Parameters are the files/directories to process, and option flags.
23* May be used to:
24* 1- continue an interupted simulation, and then add (--merge option) the similar
25*    data files.
26* 2- convert original simple text format into NeXus format (--format=NeXus)
27* 3- catenate list files
28* 4- reconstruct parameter scan sumulations (--scan)
29*
30* %Parameters
31* INPUT PARAMETERS:
32* files [string]:  list of files/directories to convert, separated by spaces
33*
34* %Link
35*
36* %End
37*******************************************************************************/
38
39#ifndef MCFORMAT
40#define MCFORMAT  "$Revision$" /* avoid memory.c to define Pool functions */
41#endif
42
43#ifdef USE_MPI
44#undef USE_MPI
45#endif
46
47/* Instead of including mccode.h file, which would then require to link most
48   of the functions of the mcstas executable, we just put there some parts of
49   the mccode.h file.
50   Build: cc -o mcformat mcformat.c -lm
51*/
52
53#define fatal_error printf  /* remove debug.c dependency */
54#define debug(msg)
55
56#define MCCODE_H                  /* avoids memory.c to import mccode.h */
57typedef struct Pool_header *Pool; /* allows memory to be included */
58#include "memory.c"
59
60#include "../lib/share/mccode-r.h" /* with decl of MC_PATHSEP */
61#include "../lib/share/mccode-r.c"
62
63/* end of parts copied from mccode.h */
64#include "../lib/share/read_table-lib.h" /* independent library */
65#include "../lib/share/read_table-lib.c"
66
67#include <dirent.h>
68#include <errno.h>
69
70/* Functions defined in memory.c */
71
72void *mem(size_t);    /* Allocate memory. */
73void  memfree(void *);    /* Free memory. */
74char *str_dup(char *);    /* Allocate new copy of string. */
75char *str_dup_n(char *string, int n); /* Copies only first N chars. */
76char *str_cat(char *first, ...);/* Concatenate strings to allocated string. */
77char *str_quote(char *string);  /* Quote string for inclusion in C code */
78void  str_free(char *);   /* Free memory for string. */
79
80/* default global variables required by mccode-r (usually generated in cogen) */
81int  mcdefaultmain         = 1;
82int  mctraceenabled        = 0;
83char mcinstrument_name[CHAR_BUF_LENGTH];
84char mcinstrument_source[CHAR_BUF_LENGTH];
85char *mcinstrument_exe = NULL;
86int  mcnumipar             = 0;
87struct mcinputtable_struct mcinputtable[CHAR_BUF_LENGTH];
88mcstatic FILE *mcsiminfo_file        = NULL;
89
90/* default global variables for mcformat converter */
91int  files_to_convert_NB    = 0;          /* nb of files to convert */
92int  files_to_convert_Array[CHAR_BUF_LENGTH];  /* index of argv[] for these files */
93char mcforcemode =0;
94char mcverbose   =0;
95char mcmergemode =0;
96char mcscanmode  =0;
97char mcmergesamedir=0;
98int  mcdircount  =0; /* number of directories scanned */
99int  mcnbconvert =0; /* index of current file */
100FILE **mcsimfiles;
101char **mcdirnames;
102char **mcinstrnames;
103char **mcsources;
104char *mcoutputdir=NULL;
105int  ipar_var    =0;  /* column index of scan variable (used in mcformat_scan_compare) */
106
107struct fileparts_struct {
108  char *FullName;
109  char *Path;
110  char *Name;
111  char *Extension;
112};
113
114/* list of functions ******************************************************** */
115/*
116  char                     *str_dup_numeric(char *orig)
117  char                     *str_dup_name(char *orig, int length)
118  char                     *str_dup_label(char *orig)
119  char                     *str_last_word(char *orig)
120  struct                    fileparts_struct fileparts_init(void)
121  struct                    fileparts_struct fileparts(char *name)
122  void                      fileparts_free(struct fileparts_struct parts)
123  struct McStas_file_format mcformat_init_mcstas_struct(void)
124  void                      mcformat_print_mcstas_struct(struct McStas_file_format McStasStruct)
125  struct McStas_file_format mcformat_free_mcstas_struct(struct McStas_file_format McStasStruct)
126  struct McStas_file_format mcformat_read_mcstas(char *filename)
127  int                       mcformat_dirwalk(char *dir, int (*fcn)(char *))
128  static void               mcformat_usedir(char *dir)
129  int                       mcformat_output(struct McStas_file_format McStasStruct)
130  int                       mcformat_convert(char *name)
131  int                       mcformat_count(char *name)
132  int                       mcformat_merge_compare(int nb)
133  void                      mcformat_scan_compare(int nb)
134  int                       mcformat_merge_output(int nb)
135  void                      mcformat_usage(char *pgmname)
136  void                      mcformat_parseoptions(int argc, char *argv[])
137  int                       main(int argc, char *argv[])
138*/
139
140/*******************************************************************************
141* str_dup_numeric: makes a clean copy of a string and allocate as numeric
142*******************************************************************************/
143char *str_dup_numeric(char *orig)
144{
145  long i;
146  char *valid;
147
148  if (!orig || !strlen(orig)) return(NULL);
149
150  for (i=0; i < strlen(orig); i++)
151  {
152    if ( (orig[i] > 122)
153      || (orig[i] < 32)
154      || (strchr("!\"#$%&'()*,:;<=>?@[\\]^`/ ", orig[i]) != NULL) )
155    {
156      orig[i] = ' ';
157    }
158  }
159  orig[i] = '\0';
160  /* now skip spaces */
161  for (i=0; i < strlen(orig); i++) {
162    if (*orig == ' ') orig++;
163    else break;
164  }
165
166  valid=str_dup(orig);
167
168  return(valid);
169} /* str_dup_numeric */
170
171/*******************************************************************************
172* str_dup_name: makes a clean copy of a string and allocate as name
173*******************************************************************************/
174char *str_dup_name(char *orig, int length)
175{
176  char *valid;
177
178  if (!orig || !strlen(orig)) return(NULL);
179
180  valid=str_dup_n(orig, length > 0 ? length : strlen(orig));
181  strcpy_valid(valid, orig);
182
183  return(valid);
184} /* str_dup_name */
185
186/*******************************************************************************
187* str_dup_label: makes a clean copy of a string and allocate as label/title
188*******************************************************************************/
189char *str_dup_label(char *orig)
190{
191  long i;
192  char *valid;
193  int  length=0;
194
195  if (!orig || !strlen(orig)) return(NULL);
196
197  for (i=0; i < strlen(orig); i++)
198  {
199    if ( (orig[i] > 122)
200      || (orig[i] < 32)
201      || (strchr("\"'=:;", orig[i]) != NULL) )
202    {
203      orig[i] = ' ';
204    }
205  }
206  orig[i] = '\0';
207  /* now skip spaces at begining */
208  for (i=0; i < strlen(orig); i++) {
209    if (*orig == ' ') orig++;
210    else break;
211  }
212  length = strlen(orig);
213  /* skip trailing blanks */
214  for (i=strlen(orig)-1; i > 0; i--) {
215    if (orig[i] == ' ') length--;
216    else break;
217  }
218
219  valid=str_dup_n(orig, length);
220
221  return(valid);
222} /* str_dup_label */
223
224/*******************************************************************************
225* str_last_word: points to last word of string
226*******************************************************************************/
227char *str_last_word(char *orig)
228{
229  char *pos_end  =NULL;
230  char *pos_begin=NULL;
231  char separators[]= MC_PATHSEP_S " !\"#$%&'()*,:;<=>?@[\\]^`/.";
232
233  /* first skip trailing separators */
234  pos_end = orig+strlen(orig);
235  while (pos_end > orig) {
236    if (strchr(separators, *pos_end)) pos_end--; /* pass separators at the end */
237    else break;
238  }
239  pos_begin = pos_end-1;
240  /* search for non separators (pass word) */
241  while (pos_begin >= orig) {
242    if (!strchr(separators, *pos_begin)) pos_begin--; /* pass non separators */
243    else break;
244  }
245  pos_begin++;
246
247  if (pos_begin < orig) pos_begin=orig;
248
249  return(pos_begin);
250} /* str_last_word */
251
252
253/*******************************************************************************
254* fileparts_init: Initialize a zero fileparts structure
255*******************************************************************************/
256struct fileparts_struct fileparts_init(void)
257{
258  struct fileparts_struct parts;
259
260  parts.FullName  = NULL;
261  parts.Path      = NULL;
262  parts.Name      = NULL;
263  parts.Extension = NULL;
264
265  return(parts);
266} /* fileparts_init */
267
268/*******************************************************************************
269* fileparts: Split a fully qualified file name/path into pieces
270* Returns a zero structure if called with NULL argument.
271* Returns: fields are non NULL if they exist
272*    Path is NULL if no Path
273*    Name is NULL if just a Path
274*    Extension is "" if just a dot
275*******************************************************************************/
276struct fileparts_struct fileparts(char *name)
277{
278  struct fileparts_struct parts;
279
280  parts = fileparts_init();
281
282  if (name) {
283    char *dot_pos    = NULL;
284    char *path_pos   = NULL;
285    char *end_pos    = NULL;
286    char *name_pos   = NULL;
287    long  dot_length = 0;
288    long  path_length= 0;
289    long  name_length= 0;
290
291    parts.FullName  = str_dup(name);
292    /* extract path+filename+extension from full filename */
293
294    if (strlen(name) == 0) return(parts);
295
296    end_pos = name+strlen(name);  /* end of file name */
297
298    /* extract path: searches for last file separator */
299    path_pos= strrchr(name, MC_PATHSEP_C);  /* last PATHSEP */
300
301    if (!path_pos) {
302      path_pos   =name;
303      path_length=0;
304      name_pos   =name;
305      parts.Path = str_dup("");
306    } else {
307      name_pos    = path_pos+1;
308      path_length = name_pos - name;  /* from start to path+sep */
309      if (path_length) {
310        parts.Path = str_cat(name, MC_PATHSEP_S, NULL);
311        strncpy(parts.Path, name, path_length);
312        parts.Path[path_length]='\0';
313      } else parts.Path = str_dup("");
314    }
315
316    /* extract ext: now looks for the 'dot' */
317    dot_pos = strrchr(name_pos, '.');           /* last dot */
318    if (dot_pos > name_pos) {
319      dot_length = end_pos - dot_pos;
320      if (dot_length > 0) {
321        parts.Extension = str_dup(name);
322        strncpy(parts.Extension, dot_pos+1, dot_length);  /* skip the dot */
323        parts.Extension[dot_length]='\0';
324      }
325    } else dot_pos = end_pos;
326
327    /* extract Name (without extension) */
328    name_length = dot_pos - name_pos; /* from path to dot */
329    if (name_length) {
330      parts.Name = str_dup(name);
331      strncpy(parts.Name, name_pos, name_length);
332      parts.Name[name_length]='\0';
333    }
334  } /* if (name) */
335  return (parts);
336} /* fileparts */
337
338/*******************************************************************************
339* fileparts_free: Free a fileparts_struct fields
340*******************************************************************************/
341void fileparts_free(struct fileparts_struct parts)
342{
343  str_free(parts.FullName);
344  str_free(parts.Path);
345  str_free(parts.Name);
346  str_free(parts.Extension);
347} /* fileparts_free */
348
349struct McStas_file_format {
350  char   *Source;
351  char   *Creator;
352  char   *Editor;
353  char   *InstrName;
354  char   *filename;
355  char   *Format;
356  char   *Date;
357  char   *EndDate;
358  double  Ncount;
359  char   *component;
360  char   *type;
361  char   *title;
362  char   *xlabel;
363  char   *ylabel;
364  char   *zlabel;
365  char   *xvar;
366  char   *yvar;
367  char   *zvar;
368  char   *xylimits;
369  char   *position;
370  char   *gravitation;
371  char   *outputname;
372  char   *Ncount_str;
373
374  t_Table *Data;
375
376  long    m, n, p;
377  double  x1, x2, y1, y2, z1, z2;
378  char   *ratio;
379  double  RunNum;
380  double *p0, *p1, *p2;
381  double Nsum, Psum, P2sum;
382  Coords  POSITION;
383
384  int  mcnumipar;
385  struct mcinputtable_struct *mcinputtable;
386  char *mcdirname;
387  double Scan_ipar_value;
388  int    Scan_ipar_index;
389  double Scan_mon_distance;
390  int    Scan_mon_index;
391}; /* McStas_file_format */
392
393struct McStas_file_format *Files_to_Merge=NULL;
394int                       *Scans_to_merge=NULL;
395
396/*******************************************************************************
397* mcformat_init_mcstas_struct: Init an empty McStas data structure
398*******************************************************************************/
399struct McStas_file_format mcformat_init_mcstas_struct(void)
400{
401  struct McStas_file_format McStasStruct;
402
403  McStasStruct.Format     = NULL;
404  McStasStruct.Creator    = NULL;
405  McStasStruct.Editor     = NULL;
406  McStasStruct.Date       = NULL;
407  McStasStruct.EndDate    = NULL;
408  McStasStruct.Ncount     = 0;
409  McStasStruct.type       = NULL;
410  McStasStruct.Source     = NULL;
411  McStasStruct.InstrName  = NULL;
412  McStasStruct.component  = NULL;
413  McStasStruct.title      = NULL;
414  McStasStruct.ratio      = NULL;
415  McStasStruct.filename   = NULL;
416  McStasStruct.xvar       = NULL;
417  McStasStruct.yvar       = NULL;
418  McStasStruct.xlabel     = NULL;
419  McStasStruct.ylabel     = NULL;
420  McStasStruct.zvar       = NULL;
421  McStasStruct.zlabel     = NULL;
422  McStasStruct.xylimits   = NULL;
423  McStasStruct.position   = NULL;
424  McStasStruct.outputname = NULL;
425  McStasStruct.Ncount_str = NULL;
426
427  McStasStruct.Data       = NULL;
428  McStasStruct.p0         = NULL;
429  McStasStruct.p1         = NULL;
430  McStasStruct.p2         = NULL;
431  McStasStruct.Nsum       = 0;
432  McStasStruct.Psum       = 0;
433  McStasStruct.P2sum      = 0;
434  McStasStruct.POSITION.x=McStasStruct.POSITION.y=McStasStruct.POSITION.z=0;
435
436  McStasStruct.RunNum     = 0;
437  McStasStruct.m          = 0;
438  McStasStruct.n          = 0;
439  McStasStruct.p          = 0;
440  McStasStruct.x1         = 0;
441  McStasStruct.x2         = 0;
442  McStasStruct.y1         = 0;
443  McStasStruct.y2         = 0;
444  McStasStruct.z1         = 0;
445  McStasStruct.z2         = 0;
446
447  McStasStruct.gravitation= NULL;
448
449  McStasStruct.mcnumipar          =0;
450  McStasStruct.mcinputtable       =NULL;
451  McStasStruct.mcdirname          =NULL;
452
453  McStasStruct.Scan_ipar_value  =0;
454  McStasStruct.Scan_ipar_index  =0;
455  McStasStruct.Scan_mon_distance=0;
456  McStasStruct.Scan_mon_index   =0;
457
458  return(McStasStruct);
459} /* mcformat_init_mcstas_struct */
460
461
462void mcformat_print_mcstas_struct(struct McStas_file_format McStasStruct)
463{
464  printf("Structure from file %s\n", McStasStruct.filename);
465  printf("  Format     = %s\n", McStasStruct.Format     );
466  printf("  Creator    = %s\n", McStasStruct.Creator    );
467  printf("  Editor     = %s\n", McStasStruct.Editor     );
468  printf("  Date       = %s\n", McStasStruct.Date       );
469  printf("  EndDate    = %s\n", McStasStruct.EndDate    );
470  printf("  type       = %s\n", McStasStruct.type       );
471  printf("  Source     = %s\n", McStasStruct.Source     );
472  printf("  InstrName  = %s\n", McStasStruct.InstrName  );
473  printf("  component  = %s\n", McStasStruct.component  );
474  printf("  title      = %s\n", McStasStruct.title      );
475  printf("  ratio      = %s\n", McStasStruct.ratio      );
476  printf("  xvar       = %s\n", McStasStruct.xvar       );
477  printf("  yvar       = %s\n", McStasStruct.yvar       );
478  printf("  xlabel     = %s\n", McStasStruct.xlabel     );
479  printf("  ylabel     = %s\n", McStasStruct.ylabel     );
480  printf("  zvar       = %s\n", McStasStruct.zvar       );
481  printf("  zlabel     = %s\n", McStasStruct.zlabel     );
482  printf("  xylimits   = %s\n", McStasStruct.xylimits   );
483  printf("  position   = %s\n", McStasStruct.position   );
484  printf("  outputname = %s\n", McStasStruct.outputname );
485  printf("  Ncount_str = %s\n", McStasStruct.Ncount_str );
486  printf("  mcdirname  = %s\n", McStasStruct.mcdirname );
487
488  printf("  gravitation= %s\n", McStasStruct.gravitation );
489  printf("  Ncount     = %g\n", McStasStruct.Ncount);
490  printf("  RunNum     = %g\n", McStasStruct.RunNum);
491
492  printf("  x1         = %g\n", McStasStruct.x1);
493  printf("  x2         = %g\n", McStasStruct.x2);
494  printf("  y1         = %g\n", McStasStruct.y1);
495  printf("  y2         = %g\n", McStasStruct.y2);
496  printf("  z1         = %g\n", McStasStruct.z1);
497  printf("  z2         = %g\n", McStasStruct.z2);
498
499  printf("  Data       = %s\n", McStasStruct.Data ? "OK" : "NULL");
500
501  printf("  m          = %ld\n", McStasStruct.m);
502  printf("  n          = %ld\n", McStasStruct.n);
503  printf("  p          = %ld\n", McStasStruct.p);
504
505  printf("  p0         = %s\n", McStasStruct.p0 ? "OK": "NULL");
506  printf("  p1         = %s\n", McStasStruct.p1 ? "OK": "NULL");
507  printf("  p2         = %s\n", McStasStruct.p2 ? "OK": "NULL");
508  printf("  POSITION   = %g %g %g\n", McStasStruct.POSITION.x, McStasStruct.POSITION.y, McStasStruct.POSITION.z);
509  printf("  mcnumipar  = %d\n", McStasStruct.mcnumipar);
510  printf("  inputtable = %s\n", McStasStruct.mcinputtable ? "OK": "NULL");
511}
512
513
514/*******************************************************************************
515* mcformat_free_mcstas_struct: Free a McStas data structure
516*******************************************************************************/
517struct McStas_file_format mcformat_free_mcstas_struct(struct McStas_file_format McStasStruct)
518{
519  memfree(McStasStruct.Format   );
520  memfree(McStasStruct.Date     );
521  memfree(McStasStruct.EndDate  );
522  memfree(McStasStruct.Creator  );
523  memfree(McStasStruct.Editor   );
524  memfree(McStasStruct.type     );
525  memfree(McStasStruct.Source   );
526  memfree(McStasStruct.InstrName);
527  memfree(McStasStruct.component);
528  memfree(McStasStruct.gravitation);
529  memfree(McStasStruct.title    );
530  memfree(McStasStruct.ratio    );
531  memfree(McStasStruct.filename );
532  memfree(McStasStruct.xvar     );
533  memfree(McStasStruct.yvar     );
534  memfree(McStasStruct.xlabel   );
535  memfree(McStasStruct.ylabel   );
536  memfree(McStasStruct.zvar     );
537  memfree(McStasStruct.zlabel   );
538  memfree(McStasStruct.xylimits );
539  memfree(McStasStruct.position );
540  memfree(McStasStruct.outputname );
541  memfree(McStasStruct.Ncount_str );
542
543  if (McStasStruct.Data) Table_Free_Array(McStasStruct.Data);
544  if (McStasStruct.p0) memfree(McStasStruct.p0 );
545  if (McStasStruct.p1) memfree(McStasStruct.p1 );
546  if (McStasStruct.p2) memfree(McStasStruct.p2 );
547
548  memfree(McStasStruct.mcdirname          );
549
550  /* free mcinputtable */
551  int i;
552  for (i=0; i<McStasStruct.mcnumipar; i++) {
553    memfree(McStasStruct.mcinputtable[i].name);
554    memfree(McStasStruct.mcinputtable[i].val);
555  }
556  memfree(McStasStruct.mcinputtable       );
557
558  return(McStasStruct);
559} /* mcformat_free_mcstas_struct */
560
561/*******************************************************************************
562* mcformat_read_mcstas:  Reads filename and checks for McStas format
563*                        extracts header and data and return structure
564*                        structure.p1 is NULL in case of failure
565*******************************************************************************/
566struct McStas_file_format mcformat_read_mcstas(char *filename)
567{
568  struct McStas_file_format McStasStruct;
569  int    m=1,n=1,p=1, i,j;
570  long   array_length=0;
571  char   flag_abort=0;
572  char   flag_pgplot1d=0;
573
574  /* init default empty return value */
575  McStasStruct = mcformat_init_mcstas_struct();
576
577  /* opens the file as a t_Table */
578  t_Table *rTable=NULL;
579
580  /* gets header and block 1 */
581  rTable = Table_Read_Array(filename, &array_length); /* reads all data */
582  if (!rTable) { if (mcverbose) fprintf(stderr, "mcformat: file %s not found\n", filename);
583    return(McStasStruct);
584  }
585
586  /* returns if this is not a McStas format file, or invalid */
587  if (!rTable[0].data || (array_length != 1 && array_length != 3))   {
588    if (mcverbose) fprintf(stderr, "mcformat: file %s is invalid (%ld blocks)\n",
589      filename, array_length);
590    flag_abort=1;
591  } else if (!rTable[0].header) {
592    if (mcverbose) fprintf(stderr, "mcformat: Can not read header from file %s\n"
593      "          Data is [%ld x %ld]\n", filename, rTable[0].rows, rTable[0].columns);
594    flag_abort=1;
595  }
596
597  if (array_length == 3)
598      for (i=1; i<=2; i++)
599      if (rTable[0].rows != rTable[i].rows || rTable[0].columns != rTable[i].columns) {
600        if (mcverbose) fprintf(stderr, "mcformat: file %s Data dimensions are not consistent\n"
601          "          Data[%i] is [%ld x %ld]\n",
602          filename, i, rTable[i].rows, rTable[i].columns);
603          flag_abort=1;
604      }
605  if (flag_abort) {
606    Table_Free_Array(rTable);
607    return(McStasStruct);
608  }
609
610  /* analyze content of file : use read-table:Table_ParseHeader ********** */
611  char **parsing;
612  parsing = Table_ParseHeader(rTable[0].header,
613      "Format", /* original Format */
614      "Editor",
615      "Creator",
616      "Date",
617      "EndDate",
618      "File",   /* original filename */
619      "Source","Instrument-source",
620      "instrument",
621      "Ncount",
622      "Gravitation",
623      "type",
624      "component","parent",
625      "position",
626      "title",
627      "ratio",
628      "xlabel",
629      "ylabel",
630      "zlabel",
631      "xvar",
632      "yvar",
633      "zvar",
634      "xlimits", "xylimits",
635      NULL);
636
637    if (parsing) {
638      if (parsing[0])   McStasStruct.Format     = str_dup_label(parsing[0]);
639      if (parsing[1])   McStasStruct.Editor     = str_dup_label(parsing[1]);
640      if (parsing[2])   McStasStruct.Creator    = str_dup_label(parsing[2]);
641      if (parsing[3])   McStasStruct.Date       = str_dup_numeric(parsing[3]);  /* mcstartdate */
642      if (parsing[4])   McStasStruct.EndDate    = str_dup_numeric(parsing[4]);  /* DATL */
643/*   if (parsing[5])   McStasStruct.filename   = str_dup(parsing[5]); */
644      if (parsing[6])   McStasStruct.Source     = str_dup_label(parsing[6]);  /* mcinstrument_source */
645      if (parsing[7] && !McStasStruct.Source)
646                        McStasStruct.Source     = str_dup_label(parsing[7]);
647      if (parsing[8])   McStasStruct.InstrName  = str_dup_label(parsing[8]);  /* mcinstrument_name */
648      if (parsing[9])   McStasStruct.Ncount_str = str_dup_numeric(parsing[9]);
649      if (parsing[10])  McStasStruct.gravitation= str_dup_numeric(parsing[10]);/* mcgravitation */
650      if (parsing[11])  McStasStruct.type       = str_dup_label(parsing[11]); /* array_Xd */
651      if (parsing[12])  McStasStruct.component  = str_dup_label(parsing[12]); /* NAME_CURRENT_COMP */
652      if (parsing[13] && !McStasStruct.component)
653                        McStasStruct.component  = str_dup_label(parsing[13]);
654      if (parsing[14])  McStasStruct.position   = str_dup_numeric(parsing[14]); /* string(POS_A_CURRENT_COMP) */
655      if (parsing[15])  McStasStruct.title      = str_dup_label(parsing[15]);
656      if (parsing[16])  McStasStruct.ratio      = str_dup_numeric(parsing[16]); /* -> RunNum */
657      if (parsing[17])  McStasStruct.xlabel     = str_dup_label(parsing[17]);
658      if (parsing[18])  McStasStruct.ylabel     = str_dup_label(parsing[18]);
659      if (parsing[19])  McStasStruct.zlabel     = str_dup_label(parsing[19]);
660      if (parsing[20])  McStasStruct.xvar       = str_dup_label(parsing[20]);
661      if (parsing[21])  McStasStruct.yvar       = str_dup_label(parsing[21]);
662      if (parsing[22])  McStasStruct.zvar       = str_dup_label(parsing[22]);
663      if (parsing[23])  McStasStruct.xylimits   = str_dup_numeric(parsing[23]); /* x/y min/max */
664      if (parsing[24] && !McStasStruct.xylimits)
665                        McStasStruct.xylimits   = str_dup_numeric(parsing[24]);
666      for (i=0; i<=24; i++) {
667        if (parsing[i]) free(parsing[i]);
668      }
669      free(parsing);
670    }
671
672/* will be determined over original values when writting  new data file:
673        filename
674        format
675        statistics
676        signal
677        values
678     */
679  if (!McStasStruct.InstrName && !McStasStruct.Source)
680    McStasStruct.Source = str_dup("McStas_Instrument");
681  else if (!McStasStruct.InstrName && McStasStruct.Source) {
682    char *ext=NULL; /* if not found, use instr file name without extension */
683    ext = strstr(McStasStruct.Source, ".ins");
684    if (!ext) ext = strstr(McStasStruct.Source, " ins");
685    if (ext) McStasStruct.InstrName = str_dup_n(McStasStruct.Source, (ext-McStasStruct.Source));
686    else McStasStruct.InstrName = str_dup(McStasStruct.Source);
687  } else if (McStasStruct.InstrName && !McStasStruct.Source) {
688    char *ext=NULL; /* if not found, use instr name without extension */
689    ext = strstr(McStasStruct.InstrName, ".ins");
690    if (!ext) ext = strstr(McStasStruct.InstrName, " ins");
691    if (ext) McStasStruct.Source = str_dup_n(McStasStruct.InstrName, (ext-McStasStruct.InstrName));
692    else McStasStruct.Source = str_dup(McStasStruct.InstrName);
693  }
694
695  if (McStasStruct.Ncount_str) McStasStruct.Ncount=atof(McStasStruct.Ncount_str);
696
697  /* header analysis: general keywords */
698  if (McStasStruct.ratio) {  /* extracts RunNum and Ncount values separated by '/' */
699    double runnum=0, ncount=0;
700    int   i;
701    i = sscanf(McStasStruct.ratio, "%lg %lg", &runnum, &ncount);
702    if (i) McStasStruct.RunNum = runnum;
703    if (i>1) {
704      if (!McStasStruct.Ncount && ncount) McStasStruct.Ncount = ncount;
705      else if (ncount && McStasStruct.Ncount != ncount)
706        fprintf(stderr, "Warning: %s: conflicting Ncount %f value with 'ratio' one %f\n",
707          filename, McStasStruct.Ncount, ncount);
708    } else if (!McStasStruct.Ncount) /* no Ncount and ratio is a single value: Ncount == runnum */
709      McStasStruct.Ncount = McStasStruct.RunNum;
710  }
711  if (!McStasStruct.Ncount) {
712    McStasStruct.Ncount = 1e6;
713    fprintf(stderr, "Warning: %s: can not extract Ncount. using default (%g)\n",
714      filename, McStasStruct.Ncount);
715  }
716  if (!McStasStruct.RunNum) McStasStruct.RunNum = McStasStruct.Ncount;
717  if (McStasStruct.RunNum < McStasStruct.Ncount*0.99)
718    fprintf(stderr, "Warning: %s: Temporary results with ratio %g/%g. Simulation results are not completed.\n",
719      filename, McStasStruct.RunNum, McStasStruct.Ncount);
720
721  /* analyse type, and check that Data dims are m,n,p */
722  if (McStasStruct.type) {
723    if (!strcmp(McStasStruct.type, "array_0d")) m=n=p=1;
724    else if (sscanf(McStasStruct.type, "array_1d(%d)",&m) == 1) n=p=1;
725    else if (sscanf(McStasStruct.type, "array_2d(%d, %d)",&n,&m) == 2) p=1;
726    else if (sscanf(McStasStruct.type, "array_3d(%d, %d, %d)",&n,&m,&p) == 3) { /* void */ }
727    else if (sscanf(McStasStruct.type, "multiarray_1d(%d)", &m)) {
728      if (!mcscanmode) fprintf(stderr, "Warning: %s: use --scan flag for multiarray/scans (skipped)\n", filename);
729      return(McStasStruct);
730    }
731  } else {
732    fprintf(stderr, "Warning: %s: invalid data type '%s' (not 1d/2d/3d/multiarray)\n", filename, McStasStruct.type);
733    return(McStasStruct);
734  }
735
736/* Scans: just test existence. If yes, then will skip the SIM file */
737  /*
738  "# Numpoints: "
739  "# xvars: "
740  "# yvars: "
741*/
742  /* check if data block is transposed */
743  if (m != rTable[0].rows    && n != rTable[0].columns
744   && m == rTable[0].columns && n == rTable[0].rows) {
745    if (mcverbose) fprintf(stderr, "Warning: %s: Data block is transposed. Fixing.\n", filename);
746    int tmp=m;
747    m=n; n=tmp;
748  }
749  if (m != rTable[0].rows) {
750    fprintf(stderr, "Warning: %s: conflicting Data row numbers\n"
751                    "         expected=%d found=%ld. Fixing. Check first line of Data block.\n", filename, m, rTable[0].rows);
752    m = rTable[0].rows;
753  }
754  if (McStasStruct.Format && (strstr(McStasStruct.Format, "binary") || strstr(McStasStruct.Format, "float") || strstr(McStasStruct.Format, "double")))
755    fprintf(stderr, "WARNING: %s: Format of data file indicates binary blocks."
756                    "         Not supported. Might crash.\n", filename);
757  if (McStasStruct.Format && McStasStruct.type
758    && strstr(McStasStruct.Format, "PGPLOT") && array_length == 1
759    && strstr(McStasStruct.type, "array_1d") && rTable[0].columns == 4) flag_pgplot1d=1;
760  if (flag_pgplot1d) n = 4;
761  if (n != rTable[0].columns) {
762    fprintf(stderr, "Warning: %s: conflicting Data column numbers\n"
763                    "         expected=%d found=%ld. Fixing. Check first line of Data block.\n", filename, n, rTable[0].columns);
764    n = rTable[0].columns;
765  }
766  if (flag_pgplot1d) n = 1;
767  McStasStruct.m = rTable[0].rows; /* dimensions from the Table */
768  McStasStruct.n = flag_pgplot1d ? 1 : rTable[0].columns;
769  McStasStruct.p = p; /* extracted from type */
770
771  m = McStasStruct.m;
772  n = McStasStruct.n;
773
774
775  /* extract limits values */
776  McStasStruct.x1=1; McStasStruct.x2=McStasStruct.m;
777  McStasStruct.y1=1; McStasStruct.y2=McStasStruct.n;
778  McStasStruct.z1=1; McStasStruct.z2=McStasStruct.p;
779  if (McStasStruct.xylimits) {
780    i = sscanf(McStasStruct.xylimits, "%lg %lg %lg %lg %lg %lg",
781      &(McStasStruct.x1), &(McStasStruct.x2),
782      &(McStasStruct.y1), &(McStasStruct.y2),
783      &(McStasStruct.z1), &(McStasStruct.z2));
784    if (i != 2 && i != 4 && i != 6)
785      fprintf(stderr, "Warning: %s: invalid xylimits '%s'. extracted %i values\n",
786        filename, McStasStruct.xylimits, i);
787  }
788
789  /* special case of 2D binary files: axes must be exchanged */
790  if (strstr(mcformat, "binary") || strstr(mcformat, "float") || strstr(mcformat, "double"))
791  if (p == 1 && m>1 && n>1) {
792    double tmp; char* c;
793    tmp=McStasStruct.x1; McStasStruct.x1=McStasStruct.y1; McStasStruct.y1=tmp;
794    tmp=McStasStruct.x2; McStasStruct.x2=McStasStruct.y2; McStasStruct.y2=tmp;
795    c=McStasStruct.xlabel; McStasStruct.xlabel=McStasStruct.ylabel; McStasStruct.ylabel=c;
796    c=McStasStruct.xvar;   McStasStruct.xvar=McStasStruct.yvar;     McStasStruct.yvar=c;
797  }
798
799  McStasStruct.filename   = str_dup(filename);
800  if (!McStasStruct.component) McStasStruct.component=str_dup(filename);
801  struct fileparts_struct file_parts=fileparts(filename);
802  char *tmp=str_cat(file_parts.Name, ".", file_parts.Extension, NULL);
803  McStasStruct.outputname = str_dup_name(tmp, 0);
804  memfree(tmp);
805  fileparts_free(file_parts);
806  McStasStruct.mcdirname  = str_dup(mcdirname);
807  if (McStasStruct.position) sscanf(McStasStruct.position, "%lg %lg %lg",
808    &McStasStruct.POSITION.x,&McStasStruct.POSITION.y,&McStasStruct.POSITION.z);
809
810  /* header analysis: get instrument parameters and fills numipar and inputtable */
811  char *s = rTable[0].header;
812  char *tok=s;
813  char *equal_sign=NULL;
814  char *name_start=NULL;
815  mcnumipar = 0;
816  while (tok) { /* extract parameter=value as clean names */
817    tok = strstr(s, "Parameters");
818    if (!tok) tok = strstr(s, "parameters");
819    if (!tok) tok = strstr(s, "Param");
820    if (!tok) tok = strstr(s, "param");
821    if (!tok) break;
822    parsing = Table_ParseHeader(tok, "Parameters", NULL); /* get line */
823    if (!parsing[0]) parsing = Table_ParseHeader(tok, "parameters", NULL); /* get line */
824    if (!parsing[0]) parsing = Table_ParseHeader(tok, "Param", NULL); /* get line */
825    if (!parsing[0]) parsing = Table_ParseHeader(tok, "param", NULL); /* get line */
826    name_start = (parsing[0] ? str_dup(parsing[0]) : NULL);
827    memfree(parsing[0]); free(parsing);
828    if (!name_start) break;
829    equal_sign = strchr(name_start+1, '=');
830    if (equal_sign > name_start && strlen(name_start)) {
831      char *name_to_equal=str_dup_n(name_start, equal_sign-name_start);
832      char *name=str_last_word(name_to_equal);
833      char *value      = str_dup(equal_sign+1);
834      if (name && value && strlen(name) && strlen(value)) {
835        char *name_label = str_dup_label(name);
836        /*printf("name_to_equal='%s' name='%s' value='%s'\n", name_to_equal, name, value); */
837        mcinputtable[mcnumipar].name = name_label;
838        mcinputtable[mcnumipar].type = instr_type_string;
839        mcinputtable[mcnumipar].val  = value;
840        mcinputtable[mcnumipar].par  = NULL;
841        mcnumipar++;
842      }
843      memfree(name_to_equal);
844    }
845    memfree(name_start);
846    s = tok+strlen("Param");
847  } /* end while tok */
848
849  /* now transfer mcinputtable into McStasStruct */
850  McStasStruct.mcnumipar = mcnumipar;
851  McStasStruct.mcinputtable = (struct mcinputtable_struct *)mem(mcnumipar*sizeof(struct mcinputtable_struct));
852  for (i=0; i<mcnumipar; i++) {
853    McStasStruct.mcinputtable[i] = mcinputtable[i];
854  }
855
856  McStasStruct.Data = rTable;
857
858  /* transfer Data into p0, p1, p2 */
859  McStasStruct.p1   = (double*)mem(McStasStruct.m*McStasStruct.n*McStasStruct.p*sizeof(double));
860  if ((array_length == 3 || flag_pgplot1d) && !strstr(McStasStruct.Format, " list ")) {
861    McStasStruct.p0 = (double*)mem(McStasStruct.m*McStasStruct.n*McStasStruct.p*sizeof(double));
862    McStasStruct.p2 = (double*)mem(McStasStruct.m*McStasStruct.n*McStasStruct.p*sizeof(double));
863  }
864
865  if (flag_pgplot1d)
866  { /* this is a 1D PGPLOT data file.
867       columns [ xvar, p, err, N ] */
868    memfree(McStasStruct.yvar); McStasStruct.yvar=str_dup("I");
869    for (i=0; i<m; i++) {
870      McStasStruct.p1[i] = Table_Index(McStasStruct.Data[0], i, 1);
871      McStasStruct.p2[i] = Table_Index(McStasStruct.Data[0], i, 2);
872      McStasStruct.p0[i] = Table_Index(McStasStruct.Data[0], i, 3);
873    }
874  } else if (array_length == 1 || strstr(McStasStruct.Format, " list ")) {
875    /* lists are stored as is */
876    for (i=0; i<m; i++)
877      for (j=0; j <n*p; j++) {
878        int tmp=i*n*p+j;
879        McStasStruct.p1[tmp] =
880          Table_Index(McStasStruct.Data[0], i, j);
881    }
882  } else if (array_length == 3) { /* order in files is [p, err, N] */
883    for (i=0; i<m; i++)
884      for (j=0; j <n*p; j++) {
885        int tmp=i*n*p+j;
886        McStasStruct.p1[tmp] = Table_Index(McStasStruct.Data[0], i, j);
887        McStasStruct.p2[tmp] = Table_Index(McStasStruct.Data[1], i, j);
888        McStasStruct.p0[tmp] = Table_Index(McStasStruct.Data[2], i, j);
889      }
890  } else
891    fprintf(stderr, "Warning: %s: can not store data blocks (total %ld).\n",
892                    filename, array_length);
893
894  /* compute back the original sigma->p2 array (see mcestimate_error) */
895  if (McStasStruct.p0 && McStasStruct.p1 && McStasStruct.p2) {
896    double *p0 = McStasStruct.p0;
897    double *p1 = McStasStruct.p1;
898    double *p2 = McStasStruct.p2;
899    for (i=0; i < m*n*p; i++) {
900      if (p2[i] > 0) p2[i] = (p0[i] > 1 ? ((p0[i]-1)*p2[i]*p2[i] + p1[i]*p1[i]/p0[i])/p0[i]
901                         : p1[i]);
902    }
903  }
904  McStasStruct.m *= -1; /* always transposed in files w/r to memory */
905
906  /* free McStasStruct.Data for better memory management */
907  Table_Free_Array(McStasStruct.Data); McStasStruct.Data=NULL;
908
909  return(McStasStruct);
910} /* end mcformat_read_mcstas */
911
912/*******************************************************************************
913* mcformat_dirwalk:  apply fcn to all files in dir (see K & R 'C language')
914*                    returns number of processed files
915*******************************************************************************/
916int mcformat_dirwalk(char *dir, int (*fcn)(char *))
917{
918  char name[CHAR_BUF_LENGTH];
919  int  ret=0;
920  struct dirent *dp;
921  DIR *dfd;
922
923  if ((dfd = opendir(dir)) == NULL) {
924    fprintf(stderr, "mcformat: can't open %s\n", dir);
925    return 0;
926  }
927  while ((dp = readdir(dfd)) != NULL) {
928    if (strcmp(dp->d_name, ".") == 0 || strcmp(dp->d_name, "..") == 0)
929      continue;    /* skip self and parent */
930    if (strlen(dir)+strlen(dp->d_name)+2 > sizeof(name))
931      fprintf(stderr, "mcformat: name %s%c%s too long\n",
932        dir, MC_PATHSEP_C, dp->d_name);
933    else {
934      sprintf(name, "%s%c%s", dir, MC_PATHSEP_C, dp->d_name);
935      ret += (*fcn)(name);
936    }
937  }
938  closedir(dfd);
939  return(ret);
940} /* end mcformat_dirwalk */
941
942/*******************************************************************************
943* mcformat_usedir: set directory to use. Same as mcuse_dir with force option.
944*******************************************************************************/
945static void mcformat_usedir(char *dir)
946{
947#ifdef MC_PORTABLE
948  fprintf(stderr, "Error: "
949          "Directory output cannot be used with portable simulation (mcformat_usedir)\n");
950  return;
951#else  /* !MC_PORTABLE */
952#ifndef WIN32
953  if(!mcdisable_output_files && mkdir(dir, 0777))
954#else
955  if(!mcdisable_output_files && mkdir(dir))
956#endif
957  {
958    int errno_mkdir = errno;
959    if (errno_mkdir == ENOENT) {
960      if (mcverbose) fprintf(stderr, "mkdir: ENOENT A directory component in pathname '%s' does not exist or is a dangling symbolic link.\n", dir);
961      /* we build the required elements in the path */
962      struct fileparts_struct dir_parts = fileparts(dir);
963      if (strlen(dir_parts.Path)) {
964        char *path_pos= strrchr(dir_parts.Path, MC_PATHSEP_C);  /* last PATHSEP */
965        if (path_pos == dir_parts.Path+strlen(dir_parts.Path)-1) dir_parts.Path[strlen(dir_parts.Path)-1] = '\0';
966        if (mcverbose) fprintf(stderr, "mcformat: Warning: building output directory '%s' from '%s'.\n", dir_parts.Path, dir);
967        mcformat_usedir(dir_parts.Path);
968        fileparts_free(dir_parts);
969        mcformat_usedir(dir);
970      }
971    } else if (errno_mkdir == EEXIST) {
972      fprintf(stderr, "mkdir: EEXIST pathname %s already exists (not necessarily as a directory).\n", dir);
973      if (!mcforcemode && mcscanmode != 2) {
974        fprintf(stderr, "Error: unable to create directory '%s' (mcformat_usedir)\n", dir);
975        fprintf(stderr, "(Maybe the directory already exists? Use --force, --scan-only or --test before -d %s to override)\n", dir);
976        exit(1);
977      }
978      fprintf(stderr, "mcformat: Warning: re-using output directory '%s'.\n", dir);
979    } else {
980      switch (errno_mkdir) {
981#ifdef EACCES
982      case EACCES:
983        fprintf(stderr, "mkdir: EACCES The parent directory does not allow write permission\n"
984        "       or one of the directories in pathname did not allow search permission.\n"); break;
985#endif
986#ifdef EFAULT
987      case EFAULT:
988        fprintf(stderr, "mkdir: EFAULT pathname points outside your accessible address space.\n"); break;
989#endif
990#ifdef ELOOP
991      case ELOOP:
992        fprintf(stderr, "mkdir: ELOOP Too many symbolic links were encountered in resolving pathname.\n"); break;
993#endif
994#ifdef ENAMETOOLONG
995      case ENAMETOOLONG:
996        fprintf(stderr, "mkdir: ENAMETOOLONG pathname was too long.\n"); break;
997#endif
998#ifdef ENOMEM
999      case ENOMEM:
1000        fprintf(stderr, "mkdir: ENOMEM Insufficient kernel memory was available.\n"); break;
1001#endif
1002#ifdef ENOSPC
1003      case ENOSPC:
1004        fprintf(stderr, "mkdir: ENOSPC The device containing pathname has no room for the new directory.\n"); break;
1005#endif
1006#ifdef ENOTDIR
1007      case ENOTDIR:
1008        fprintf(stderr, "mkdir: ENOTDIR A component used as a directory in pathname is not, in fact, a directory.\n"); break;
1009#endif
1010#ifdef EPERM
1011      case EPERM:
1012        fprintf(stderr, "mkdir: EPERM The filesystem containing pathname does not support the creation of directories.\n"); break;
1013#endif
1014#ifdef EROFS
1015      case EROFS:
1016        fprintf(stderr, "mkdir: EROFS pathname refers to a file on a read-only filesystem.\n"); break;
1017#endif
1018      default:
1019        fprintf(stderr, "mkdir: ERROR %i using mkdir.\n", errno_mkdir); break;
1020      }
1021      fprintf(stderr, "mcformat: Fatal error accessing %s. Aborting.\n", dir);
1022      exit(-1);
1023    }
1024  }
1025  if (mcverbose) printf("mcformat: Creating directory %s.\n", dir);
1026  mcdirname = dir;
1027
1028#endif /* !MC_PORTABLE */
1029} /* mcformat_usedir */
1030
1031/*******************************************************************************
1032* mcformat_output: write monitor file and optionally free memory
1033*******************************************************************************/
1034int mcformat_output(struct McStas_file_format McStasStruct)
1035{
1036  char *currentdir= mcdirname; /* save current dir */
1037  int i;
1038
1039
1040  if (mcdisable_output_files) return(1);
1041  if (!McStasStruct.p1) return(0); /* empty data */
1042  /* determine in which directory we are and set SIM file */
1043
1044  if (!mcdircount) {
1045    if (!strlen(mcinstrument_name)) strcpy(mcinstrument_name, McStasStruct.InstrName);
1046    if (!strlen(mcinstrument_source)) strcpy(mcinstrument_source, McStasStruct.Source);
1047    if (!mcsiminfo_file) mcsiminfo_init(NULL); /* open SIM file once */
1048  } else {
1049    if (!mcmergemode) i=mcdircount-1;
1050    else {
1051      for (i=0; i<mcdircount; i++)
1052        if (!strcmp(mcdirnames[i], McStasStruct.mcdirname)) break;
1053    }
1054    if (i >= mcdircount) {
1055        fprintf(stderr, "ERROR: unable to find directory '%s' in scanned list\n", McStasStruct.mcdirname);
1056        return(0);
1057    }
1058    if (!mcsimfiles[i]) {
1059      mcdirname      = McStasStruct.mcdirname;
1060      mcinstrnames[i]= str_dup(McStasStruct.InstrName);
1061      mcsources[i]   = str_dup(McStasStruct.Source);
1062      strncpy(mcinstrument_name,     str_last_word(mcinstrnames[i]), CHAR_BUF_LENGTH);
1063      strncpy(mcinstrument_source  , str_dup(mcsources[i]), CHAR_BUF_LENGTH);
1064      mcsiminfo_init(NULL); /* open new SIM file in this dir for the first time */
1065      mcsimfiles[i]  = mcsiminfo_file;
1066    } else mcsiminfo_file = mcsimfiles[i];
1067    mcdirname = mcdirnames[i];
1068  }
1069
1070/* transfer to global variables used in output functions */
1071  if (!McStasStruct.Date) mcstartdate = 0;
1072  else {
1073    mcstartdate         = atol(McStasStruct.Date);
1074    if (!mcstartdate) sscanf(McStasStruct.Date, "Simulation started %ld", &mcstartdate);
1075  }
1076  mcgravitation       = (McStasStruct.gravitation && strstr(McStasStruct.gravitation, "yes") ? 1 : 0);
1077  mcrun_num = McStasStruct.RunNum;
1078  mcncount  = McStasStruct.Ncount;
1079  strncpy(mcinstrument_source, str_dup(McStasStruct.Source), CHAR_BUF_LENGTH);
1080  strncpy(mcinstrument_name  , str_last_word(McStasStruct.InstrName), CHAR_BUF_LENGTH);
1081  /* transfer mcnumipar */
1082  mcnumipar = McStasStruct.mcnumipar;
1083  for (i=0; i<mcnumipar; i++) {
1084    mcinputtable[i] = McStasStruct.mcinputtable[i];
1085  }
1086  mcdirname = McStasStruct.mcdirname;
1087
1088  /* write data file from McStasStruct */
1089  if (mcverbose) printf("mcformat: Writing %s from %s\n", McStasStruct.outputname, McStasStruct.filename);
1090
1091  if (strstr(McStasStruct.Format, "PGPLOT") && strstr(McStasStruct.type, "array_1d") && !mcdisable_output_files) {
1092    mcdetector_out_1D(McStasStruct.title, McStasStruct.xlabel, McStasStruct.ylabel,
1093                  McStasStruct.xvar, McStasStruct.x1, McStasStruct.x2, McStasStruct.m,
1094                  McStasStruct.p0, McStasStruct.p1, McStasStruct.p2, McStasStruct.outputname,
1095                  McStasStruct.component, McStasStruct.POSITION);
1096  } else if (!mcdisable_output_files)
1097  mcdetector_out_2D(McStasStruct.title,
1098    McStasStruct.xlabel, McStasStruct.ylabel,
1099    McStasStruct.x1, McStasStruct.x2, McStasStruct.y1, McStasStruct.y2,
1100    McStasStruct.m, McStasStruct.n,
1101    McStasStruct.p0,
1102    McStasStruct.p1,
1103    McStasStruct.p2,
1104    McStasStruct.outputname,
1105    McStasStruct.component,
1106    McStasStruct.POSITION);
1107
1108  mcdirname = currentdir;
1109  return(1);
1110} /* end mcformat_output */
1111
1112/*******************************************************************************
1113* mcformat_convert: recursive wrapper to the conversion routine (see K & R)
1114*                   returns number of processed files
1115*******************************************************************************/
1116int mcformat_convert(char *name)
1117{
1118  struct stat               stbuf;
1119  struct McStas_file_format McStasFile;
1120  int    ret=0;
1121
1122  if (!name || !strlen(name)) return(0);
1123
1124  if (stat(name, &stbuf) == -1) {
1125    fprintf(stderr, "mcformat: can't access file %s\n", name);
1126    return 0;
1127  }
1128  if ((stbuf.st_mode & S_IFMT) == S_IFDIR) { /* dir: recursive call */
1129    char *upper_dir= mcdirname;
1130    char *name_win=name;
1131#ifdef WIN32
1132    if (isalpha(name[0]) && name[1] == ':') name_win += 2;  // skip Windows disk label before cat
1133#endif
1134    char *lower_dir= mcoutputdir ? str_cat(mcoutputdir, MC_PATHSEP_S, name_win, NULL)
1135                               : str_dup(name);
1136    /* entering new directory */
1137    if (!mcdisable_output_files) mcformat_usedir(lower_dir);
1138    else mcdirname = lower_dir;
1139    mcdirnames[mcdircount] = lower_dir;
1140    mcsimfiles[mcdircount] = NULL;
1141    mcinstrnames[mcdircount]=NULL;
1142    mcsources[mcdircount]   =NULL;
1143    mcdircount++; /* number of directories scanned */
1144    ret += mcformat_dirwalk(name, mcformat_convert);
1145    if (mcverbose) printf("mcformat: coming back to upper directory %s\n", upper_dir);
1146    mcdirname = upper_dir;
1147  } else {
1148    /* process the current file */
1149    if (mcverbose) printf("mcformat: reading %s (%ld bytes)\n", name, (long int)stbuf.st_size);
1150    McStasFile = mcformat_read_mcstas(name);
1151
1152    if (!McStasFile.p1) {
1153      if (mcverbose) fprintf(stderr, "mcformat: warning: skipping %s (invalid format).\n", name);
1154      return 0;
1155    }
1156
1157    if (!mcmergemode && !mcscanmode) { /* direct conversion */
1158
1159      /* now calls the output routines: sim and data files (0d, 1d, 2d) and free */
1160      if (mcverbose) mcformat_print_mcstas_struct(McStasFile); fflush(stdout);
1161      if (mcformat_output(McStasFile))
1162        if (mcverbose) {
1163          printf("mcformat: converting %s (%ld bytes) ", name, (long int)stbuf.st_size);
1164          printf("into %s%s ", mcdirname ? mcdirname : ".", mcdirname ? MC_PATHSEP_S : "");
1165          if (mcdisable_output_files) printf(" (--test mode)\n"); else printf("\n");
1166        }
1167
1168      mcformat_free_mcstas_struct(McStasFile);
1169
1170    } else { /* store: merging of all files from either a scan or similar files */
1171      Files_to_Merge[mcnbconvert] = McStasFile; /* copy */
1172      mcnbconvert++; /* global index */
1173    }
1174
1175    ret++; /* number of valid data files loaded from DIR */
1176  } /* end else single file */
1177  return (ret);
1178} /* end mcformat_convert */
1179
1180/*******************************************************************************
1181* mcformat_count: count the number of files to store for merge mode
1182*******************************************************************************/
1183int mcformat_count(char *name)
1184{
1185  struct stat               stbuf;
1186  int    ret=0;
1187
1188  if (!name || !strlen(name)) return(0);
1189
1190  if (stat(name, &stbuf) == -1) {
1191    fprintf(stderr, "mcformat: can't access file %s\n", name);
1192    return 0;
1193  }
1194  if ((stbuf.st_mode & S_IFMT) == S_IFDIR) { /* recursive call */
1195    mcdircount++; /* number of directories scanned */
1196    ret += mcformat_dirwalk(name, mcformat_count);
1197  } else {
1198    ret++;
1199  }
1200  return (ret);
1201} /* end mcformat_count */
1202
1203/*******************************************************************************
1204* mcformat_merge_compare: merge files and free the redundant items. Find scans.
1205*******************************************************************************/
1206int mcformat_merge_compare(int nb)
1207{
1208  int i,j,k;
1209
1210  if (!Files_to_Merge) return(0);
1211  /* loop on Files_to_Merge */
1212  for (i=0; i<nb; i++) {
1213    struct McStas_file_format McStasStruct=Files_to_Merge[i];
1214
1215    /* skip empty elements */
1216    if (!McStasStruct.p1) continue;
1217
1218    /* we multiply the p1 and p2 by the Ncount so that the operations take
1219       into account the relative Ncount weight of each simulation. We will
1220       divide by the sum(Ncount) at the end of operation */
1221    if (!strstr(McStasStruct.Format, " list ")) /* NOT FOR LISTs (content non additive) */
1222    for (j=0; j<abs(McStasStruct.m*McStasStruct.n*McStasStruct.p); j++) {
1223      Files_to_Merge[i].p1[j] *= Files_to_Merge[i].Ncount;
1224      if (Files_to_Merge[i].p2) Files_to_Merge[i].p2[j] *= Files_to_Merge[i].Ncount;
1225    }
1226
1227    /* second loop to search for similar items to add to Files_to_Merge[i] */
1228    for (j=i+1; j<nb; j++) {
1229      struct McStas_file_format ThisStruct=Files_to_Merge[j];
1230      if (!ThisStruct.p1) continue;
1231
1232      /* Parameters check: mcinputtable.name and mcinputtable.val */
1233      /* flag_equiv_parname must be true  to go on */
1234      /* flag_equiv_parval  is      false if this is a SCAN item */
1235      char flag_equiv_parname = (ThisStruct.mcnumipar == McStasStruct.mcnumipar);
1236      char flag_equiv_parval  = flag_equiv_parname;
1237      if (flag_equiv_parname) {
1238        for (k=0; k<McStasStruct.mcnumipar; k++) {
1239          /* test if parameter names are the same */
1240          if (strcmp(ThisStruct.mcinputtable[k].name, McStasStruct.mcinputtable[k].name))
1241            { flag_equiv_parname=0; break; }
1242          if (strcmp(ThisStruct.mcinputtable[k].val, McStasStruct.mcinputtable[k].val))
1243            flag_equiv_parval=0;
1244        }
1245      }
1246      if (!flag_equiv_parname) continue;
1247
1248      /* Data/List: must be of same kind */
1249      char flag_list = 0;
1250        if ( strstr(ThisStruct.Format, " list ")      &&  strstr(McStasStruct.Format, " list "))
1251          flag_list=1; /* two lists */
1252        else if (!strstr(ThisStruct.Format, " list ") && !strstr(McStasStruct.Format, " list "))
1253          flag_list=2; /* two monitors */
1254      if (!flag_list) continue;
1255
1256      /* Data type check: dimension, limits, xvar, yvar, zvar, xlabel, ylabel, zlabel */
1257      /* limits may be forced (to avoid rounding errors) */
1258      /* number of rows may be different for lists */
1259      char flag_data = (
1260        (flag_list==1 || ThisStruct.m == McStasStruct.m) &&
1261        ThisStruct.n == McStasStruct.n &&
1262        ThisStruct.p == McStasStruct.p &&
1263        (!ThisStruct.xvar || !McStasStruct.xvar || !strcmp(ThisStruct.xvar, McStasStruct.xvar)) &&
1264        (!ThisStruct.yvar || !McStasStruct.yvar || !strcmp(ThisStruct.yvar, McStasStruct.yvar)) &&
1265        (!ThisStruct.zvar || !McStasStruct.zvar || !strcmp(ThisStruct.zvar, McStasStruct.zvar)) &&
1266        (!ThisStruct.xlabel || !McStasStruct.xlabel || !strcmp(ThisStruct.xlabel, McStasStruct.xlabel)) &&
1267        (!ThisStruct.ylabel || !McStasStruct.ylabel || !strcmp(ThisStruct.ylabel, McStasStruct.ylabel)) &&
1268        (!ThisStruct.zlabel || !McStasStruct.zlabel || !strcmp(ThisStruct.zlabel, McStasStruct.zlabel)) );
1269      if (!flag_data) continue;
1270
1271      char flag_limits=(mcforcemode || mcscanmode ||
1272          (ThisStruct.n == 1 && ThisStruct.x1 == McStasStruct.x1) ||
1273          (ThisStruct.x1 == McStasStruct.x1 && ThisStruct.y1 == McStasStruct.y1)
1274          );
1275      if (!flag_limits)  {
1276        fprintf(stderr, "Warning: Axes limits are not identical for %s and %s. Skipping (use --force to override).\n",
1277          McStasStruct.filename, ThisStruct.filename);
1278        continue;
1279      }
1280
1281      /* Warning if gravitation not constant (may be forced) */
1282      char flag_gravitation = mcforcemode || (
1283           ThisStruct.gravitation && McStasStruct.gravitation
1284        && !strcmp(ThisStruct.gravitation, McStasStruct.gravitation));
1285      if (!flag_gravitation && (ThisStruct.gravitation || McStasStruct.gravitation))
1286        fprintf(stderr, "Warning: Gravitation is not constant for %s and %s.\n",
1287          McStasStruct.filename, ThisStruct.filename);
1288
1289      /* Names check: InstrName, Source and component */
1290      /* directories must be different except if mcmergesamedir */
1291      char flag_names = (
1292        (!ThisStruct.InstrName || !McStasStruct.InstrName
1293            || !strcmp(ThisStruct.InstrName, McStasStruct.InstrName)) &&
1294        (!ThisStruct.Source || !McStasStruct.Source
1295            || !strcmp(ThisStruct.Source, McStasStruct.Source)) &&
1296        (!ThisStruct.component || !McStasStruct.component
1297            || !strcmp(ThisStruct.component, McStasStruct.component)) &&
1298        (mcmergesamedir || !ThisStruct.mcdirname || !McStasStruct.mcdirname
1299            || strcmp(ThisStruct.mcdirname, McStasStruct.mcdirname)) );
1300      if (!flag_names) continue;
1301
1302      /* handle scan index */
1303      if (!flag_equiv_parval) {
1304        /* attach index j to scan column origin monitor 'i' */
1305        if (Scans_to_merge[i] < 0) Scans_to_merge[i] = i;
1306        if (mcverbose && Scans_to_merge[j] < 0 && Scans_to_merge[i] >= 0) printf("  Gathering Scan step %s/%s (%d) with %s/%s (%d)\n",
1307          McStasStruct.mcdirname, McStasStruct.outputname, i,
1308          ThisStruct.mcdirname,  ThisStruct.outputname,    j);
1309        if (Scans_to_merge[j] < 0) Scans_to_merge[j] = i;
1310        /* next i if this is a scan (no add/cat) */
1311        continue; /* for j */
1312      }
1313
1314      if (!mcmergemode) continue;
1315
1316      if (mcverbose) fprintf(stderr,"  %s %s/%s (%d) with %s/%s (%d) total %ld elements\n",
1317        flag_list==2 ? "Adding" : "Appending",
1318        McStasStruct.mcdirname, McStasStruct.outputname, i,
1319        ThisStruct.mcdirname,  ThisStruct.outputname,    j,
1320        (long)abs(McStasStruct.m*McStasStruct.n*McStasStruct.p));
1321
1322      if (flag_list==1) {  /* if list: catenate data j to end of i */
1323        /* allocate new array of size rows(i+j), same n,p */
1324        double *p1=mem(abs((ThisStruct.m+McStasStruct.m)*McStasStruct.n*McStasStruct.p)*sizeof(double));
1325        /* copy data from i */
1326        int index_i, index_j, index; /* data index is index_i*columns+index_j */
1327        for (index_i=0; index_i<abs(McStasStruct.m); index_i++)
1328          for (index_j=0; index_j<abs(McStasStruct.n*McStasStruct.p); index_j++) {
1329            index=index_i*abs(McStasStruct.n*McStasStruct.p)+index_j;
1330            p1[index] = McStasStruct.p1[index];
1331        }
1332        /* catenate data from j */
1333        for (index_i=0; index_i<abs(ThisStruct.m); index_i++)
1334          for (index_j=0; index_j<abs(McStasStruct.n*McStasStruct.p); index_j++) {
1335            int index_shifted=(index_i+abs(McStasStruct.m))*abs(McStasStruct.n*McStasStruct.p)+index_j;
1336            index=index_i*abs(McStasStruct.n*McStasStruct.p)+index_j;
1337            p1[index_shifted] = ThisStruct.p1[index];
1338        }
1339        Files_to_Merge[i].m += ThisStruct.m;
1340        /* free original arrays in i and j */
1341        memfree(McStasStruct.p1); memfree(ThisStruct.p1);
1342        Files_to_Merge[j].p1 = NULL;
1343        /* attach new array to i and leave j in empty state */
1344        Files_to_Merge[i].p1 = p1;
1345      } else { /* else add data i+j */
1346        int index_i;
1347        /* add data from j to i. p1 and p2 are weighted by their Ncount */
1348        for (index_i=0; index_i<abs(McStasStruct.m*McStasStruct.n*McStasStruct.p); index_i++) {
1349          if (McStasStruct.p0) Files_to_Merge[i].p0[index_i] += ThisStruct.p0[index_i];
1350                               Files_to_Merge[i].p1[index_i] += ThisStruct.Ncount*ThisStruct.p1[index_i];
1351          if (McStasStruct.p2) Files_to_Merge[i].p2[index_i] += ThisStruct.Ncount*ThisStruct.p2[index_i];
1352
1353        }
1354        /* free and leave j in empty state */
1355        memfree(ThisStruct.p0); memfree(ThisStruct.p1); memfree(ThisStruct.p2);
1356        Files_to_Merge[j].p0=Files_to_Merge[j].p1=Files_to_Merge[j].p2=NULL;
1357      }
1358      Files_to_Merge[i].Ncount += ThisStruct.Ncount;
1359      Files_to_Merge[i].RunNum += ThisStruct.RunNum;
1360    } /* end for j */
1361
1362    /* divide p1 and p2 by total Ncount to account for the weighted sum */
1363    if (!strstr(McStasStruct.Format, " list ")) /* NOT FOR LISTs (content non additive) */
1364    for (j=0; j<abs(McStasStruct.m*McStasStruct.n*McStasStruct.p); j++) {
1365      Files_to_Merge[i].p1[j] /= Files_to_Merge[i].Ncount;
1366      if (Files_to_Merge[i].p2) Files_to_Merge[i].p2[j] /= Files_to_Merge[i].Ncount;
1367    }
1368
1369    /* now count integral values for item Files_to_Merge[i] */
1370    double Nsum=0, Psum=0, P2sum=0;
1371    for (j=0; j<abs(Files_to_Merge[i].m*Files_to_Merge[i].n*Files_to_Merge[i].p); j++) {
1372      double N=1,I=0,E=0;
1373      if (Files_to_Merge[i].p0) N = Files_to_Merge[i].p0[j];
1374      if (Files_to_Merge[i].p1) I = Files_to_Merge[i].p1[j];
1375      if (Files_to_Merge[i].p2) E = Files_to_Merge[i].p2[j];
1376      Nsum += N;
1377      Psum += I;
1378      P2sum += Files_to_Merge[i].p2 ? E : I*I;
1379    }
1380    Files_to_Merge[i].Nsum = Nsum;
1381    Files_to_Merge[i].Psum = Psum;
1382    Files_to_Merge[i].P2sum= P2sum;
1383  } /* end for i */
1384  return(nb);
1385} /* mcformat_merge_compare */
1386
1387/* sorting functions for qsort: sort monitor files with distance */
1388#ifdef HAVE_QSORT
1389int sort_distances (const void *a, const void *b)
1390{
1391  const int *pa = (const int *) a;
1392  const int *pb = (const int *) b;
1393  int ia=*pa;
1394  int ib=*pb;
1395  double da=Files_to_Merge[ia].POSITION.x*Files_to_Merge[ia].POSITION.x
1396           +Files_to_Merge[ia].POSITION.y*Files_to_Merge[ia].POSITION.y
1397           +Files_to_Merge[ia].POSITION.z*Files_to_Merge[ia].POSITION.z;
1398  double db=Files_to_Merge[ib].POSITION.x*Files_to_Merge[ib].POSITION.x
1399           +Files_to_Merge[ib].POSITION.y*Files_to_Merge[ib].POSITION.y
1400           +Files_to_Merge[ib].POSITION.z*Files_to_Merge[ib].POSITION.z;
1401  if       (da > db) return 1;
1402  else if (da < db) return -1;
1403  else /* same distance, sort on names */
1404  return strcmp(Files_to_Merge[ia].filename, Files_to_Merge[ib].filename);
1405}
1406
1407/* sorting functions for qsort: sort monitor column with ipar.val */
1408int sort_ipar_mon (const void *a, const void *b)
1409{
1410  const int *pa = (const int *) a;
1411  const int *pb = (const int *) b;
1412  int ia=*pa;
1413  int ib=*pb;
1414  double da=Files_to_Merge[ia].mcinputtable[ipar_var].type == instr_type_string ?
1415    0 : atof(Files_to_Merge[ia].mcinputtable[ipar_var].val);
1416  double db=Files_to_Merge[ib].mcinputtable[ipar_var].type == instr_type_string ?
1417    0 : atof(Files_to_Merge[ib].mcinputtable[ipar_var].val);
1418  if       (da > db) return 1;
1419  else if (da < db) return -1;
1420  else {/* same distance, sort on ipar value then filenames */
1421   int tmp=strcmp(Files_to_Merge[ia].mcinputtable[ipar_var].val,
1422                  Files_to_Merge[ib].mcinputtable[ipar_var].val);
1423   if (tmp) return(tmp);
1424   else return strcmp(Files_to_Merge[ia].filename, Files_to_Merge[ib].filename);
1425  }
1426}
1427#endif
1428
1429/*********************************************************************
1430* mcformat_scan_compare: assemble scanned monitors into multiarray sets
1431*                        warning: this function is a real headache
1432*********************************************************************/
1433void mcformat_scan_compare(int nb)
1434{
1435  /* definition of a scan
1436    Must be 'mergeable' data but with different ipar values
1437    Each Scans_to_merge[] index set is a column of multiarray_1d (this is a scanned monitor)
1438    Must analyse scan columns (Scans_to_merge[]) to gather them with
1439      -common ipar names (set of scanned monitors)
1440      -same set (monitor number) length and monitor names
1441      -same ipar values define rows
1442   */
1443  int scan_index1=0;
1444  int i=0,j;
1445
1446  /* test if there is scan data to be processed */
1447  for (scan_index1=0; scan_index1<nb; scan_index1++) {
1448    if (Scans_to_merge[scan_index1] < 0 || Scans_to_merge[scan_index1] < scan_index1) continue;
1449    else { i=1; break; }
1450  }
1451
1452  if (i == 0)
1453    fprintf(stderr, "Warning: Could not find scanned parameter within 'Param' lines in data file headers.\n"
1454          "Ignoring [mcformat:mcformat_scan_compare].\n");
1455
1456  /* loop on Scans_to_merge: index > -1  */
1457  for (scan_index1=0; scan_index1<nb; scan_index1++) {
1458    if (Scans_to_merge[scan_index1] < 0 || Scans_to_merge[scan_index1] < scan_index1) continue;
1459    /* get all sets of given index: Scans_to_merge[j] = index */
1460    int scan_length1= 1;
1461    int next_in_scan=-1;
1462    ipar_var        = 0; /* global variable, as it is used in sorting function 'sort_ipar_mon' */
1463    int scan_index2;
1464    /* compute length of monitor column (length of scan): scan_length1 */
1465    for (scan_index2=scan_index1+1; scan_index2<nb; scan_index2++)
1466      if (Scans_to_merge[scan_index1] == Scans_to_merge[scan_index2]) {
1467        scan_length1++;
1468        if (next_in_scan < 0) next_in_scan = scan_index2;
1469      }
1470    if (scan_length1 <= 1) continue; /* not a scan, single data set: go to next data file */
1471
1472    /* count other monitors in Scans_to_merge with same ipar names and values: mon_count */
1473    int mon_count=0;
1474    int Scan_columns[nb];
1475    for (scan_index2=0; scan_index2<nb; Scan_columns[scan_index2++]=-1);
1476
1477    for (scan_index2=0; scan_index2<nb; scan_index2++) {
1478      char flag_matchpar=0;
1479      if (scan_index1 == scan_index2) {
1480        mon_count++; Scan_columns[scan_index2] = scan_index1; continue;
1481      }
1482      /* must be a scan */
1483      if (Scans_to_merge[scan_index2] < 0) continue; /* scan_index2 */
1484      /* must have same numipar */
1485      if (Files_to_Merge[scan_index1].mcnumipar != Files_to_Merge[scan_index2].mcnumipar) continue; /* scan_index2 */
1486      /* must have same scan length as first monitor */
1487      int scan_length2= 0;
1488      for (j=scan_index2; j<nb; j++)
1489        if (Scans_to_merge[scan_index2] == Scans_to_merge[j]) {
1490          scan_length2++;
1491        }
1492      if (scan_length1 != scan_length2) continue; /* scan_index2 */
1493      /* must have same ipar.name and values along column */
1494      int ipar2;
1495      for (ipar2=0; ipar2<Files_to_Merge[scan_index2].mcnumipar; ipar2++) {
1496        char *this_parname = Files_to_Merge[scan_index2].mcinputtable[ipar2].name;
1497        char *this_parval  = Files_to_Merge[scan_index2].mcinputtable[ipar2].val;
1498        /* look back into 'scan_index1' series for similar parameters */
1499        for(j=scan_index1; j<nb; j++) {
1500          if (Scans_to_merge[j] == Scans_to_merge[scan_index1]) {
1501            /* in original column, look if we can find ipar */
1502            int ipar1;
1503            for (ipar1=0; ipar1<Files_to_Merge[j].mcnumipar; ipar1++) {
1504              if (strcmp(Files_to_Merge[j].mcinputtable[ipar1].name, this_parname)
1505               || strcmp(Files_to_Merge[j].mcinputtable[ipar1].val,  this_parval) ) continue; /* ipar1 */
1506              flag_matchpar=1;
1507              break;
1508            } /* for ipar1 */
1509            if (flag_matchpar) break; /* j: we found matching ipar both in scan_index1 and scan_index2 */
1510          } /* if */
1511        } /* for scan_index1b */
1512        if (!flag_matchpar) break; /* ipar2: ipar of scan_index2 not found in scan_index1 */
1513      } /* for ipar2 */
1514      if (!flag_matchpar) continue; /* scan_index2: an ipar of scan_index2 was not found in scan_index1. try other scan_index2 */
1515      else {
1516        mon_count++;
1517        /* attach that file to first monitor of dir/column */
1518        Scan_columns[scan_index2] = scan_index2;
1519      }
1520    } /* for scan_index2 */
1521
1522    /* get first varying ipar: ipar_var */
1523    for (j=0; j<Files_to_Merge[scan_index1].mcnumipar; j++) {
1524      /* compare ipar value of scan_index1 and next_in_scan: stored in j */
1525      if (strcmp(Files_to_Merge[scan_index1].mcinputtable[j].val, Files_to_Merge[next_in_scan].mcinputtable[j].val)) break;
1526    }
1527    /* test if ipar has changed in scan else go to next data file */
1528    if (j < Files_to_Merge[scan_index1].mcnumipar) ipar_var=j;
1529    /* now we know how many scanned monitor there are */
1530
1531    if (mcverbose)
1532      printf("scan from %s has %d rows (scan steps) of %d columns (monitors) varying %s\n",
1533        Files_to_Merge[scan_index1].filename, scan_length1, mon_count,
1534        Files_to_Merge[scan_index1].mcinputtable[ipar_var].name);
1535
1536    /* single Scan_columns[] is an unsorted row
1537      and have same ipar name/value in same dir/row */
1538
1539    /* single Scans_to_merge[i] value is an unsorted column
1540      and have same name in different directories/columns
1541      but different ipar values */
1542
1543    /* allocate array: scan_length1(rows)*(mcnumipar+2*mon_count) */
1544    t_Table Scan;
1545    if (!Table_Init(&Scan,
1546      scan_length1,
1547      Files_to_Merge[scan_index1].mcnumipar+2*mon_count)) {
1548      fprintf(stderr, "Warning: %s: Can not allocate Scan %ix(%i+2*%i) multiarray_1d. Skipping.\n",
1549          Files_to_Merge[scan_index1].filename, scan_length1,
1550          Files_to_Merge[scan_index1].mcnumipar,
1551          mon_count);
1552      continue;
1553    }
1554    char *header=(char*)mem(64*CHAR_BUF_LENGTH); strcpy(header, "");
1555    char *youts =(char*)mem(64*CHAR_BUF_LENGTH); strcpy(youts,  "");
1556    double ipar_min=FLT_MAX, ipar_max=0;
1557
1558    /* we first sort Scan_columns(monitors) with distance */
1559    int Scan_distances[mon_count];  /* this is column index */
1560    j = 0;
1561
1562    for (scan_index2=0; scan_index2<nb; scan_index2++)
1563      if (Scan_columns[scan_index2]>=0) Scan_distances[j++] = Scan_columns[scan_index2];
1564#ifdef HAVE_QSORT
1565    qsort(Scan_distances, mon_count, sizeof(int), sort_distances);
1566#endif
1567    for (j=0; j<mon_count; j++) { /* loop for each column */
1568      int Monitor_column[scan_length1];
1569      int k=0;
1570      for (i=0; i<scan_length1; Monitor_column[i]=i) { i++; }
1571      for (scan_index2=0; scan_index2<nb; scan_index2++) {
1572        /* find Scan_distances[j] in Scans_to_merge */
1573        if (Scans_to_merge[scan_index2] >= 0 && Scans_to_merge[scan_index2] == Scan_distances[j]) {
1574          Monitor_column[k++] = scan_index2;
1575          Scans_to_merge[scan_index2]=-1; /*unactivate that scan element */
1576        }
1577      }
1578      /* sort each monitor column with ipar value */
1579#ifdef HAVE_QSORT
1580      qsort(Monitor_column, scan_length1, sizeof(int), sort_ipar_mon);
1581#endif
1582      /* extract sorted column and set Scan */
1583      for (i=0; i<scan_length1; i++) { /* loop on column elements(row) */
1584        Table_SetElement(&Scan,
1585          i,
1586          Files_to_Merge[scan_index1].mcnumipar+2*j,
1587          Files_to_Merge[Monitor_column[i]].Psum);
1588        Table_SetElement(&Scan,
1589          i, Files_to_Merge[scan_index1].mcnumipar+2*j+1,
1590          mcestimate_error(
1591            Files_to_Merge[Monitor_column[i]].Ncount,
1592            Files_to_Merge[Monitor_column[i]].Psum,
1593            Files_to_Merge[Monitor_column[i]].P2sum));
1594        if (j==0) { /* first monitor in row also sets ipar */
1595          for (k=0; k<Files_to_Merge[Monitor_column[i]].mcnumipar; k++) {
1596            Table_SetElement(&Scan,
1597              i, k,
1598              atof(Files_to_Merge[Monitor_column[i]].mcinputtable[k].val));
1599            /* set name of ipar */
1600            if (i==0) {
1601              strcat(header, Files_to_Merge[Monitor_column[i]].mcinputtable[k].name);
1602              strcat(header, " ");
1603            }
1604          }
1605        } /* if j==0 */
1606        double this=atof(Files_to_Merge[Monitor_column[i]].mcinputtable[ipar_var].val);
1607        if (ipar_min > this) ipar_min=this;
1608        if (ipar_max < this) ipar_max=this;
1609      } /* for i */
1610      strcat(header, Files_to_Merge[Scan_distances[j]].outputname); strcat(header, "_I ");
1611      strcat(header, Files_to_Merge[Scan_distances[j]].outputname); strcat(header, "_ERR ");
1612      strcat(youts, "("); strcat(youts, Files_to_Merge[Scan_distances[j]].outputname); strcat(youts, "_I,");
1613      strcat(youts, Files_to_Merge[Scan_distances[j]].outputname); strcat(youts, "_ERR) ");
1614    } /* for j */
1615    Scan.header=header;
1616
1617    Scans_to_merge[scan_index1]=-1; /* unactivate that scan element */
1618
1619    char *title=str_cat("Scan of ",
1620      Files_to_Merge[scan_index1].mcinputtable[ipar_var].name, NULL);
1621    /* output scan multiarray */
1622    /* for PGPLOT: open mcstas.sim */
1623    if (strcasestr(mcformat, "McCode")) mcsiminfo_init(NULL); /* open new SIM file in this dir for the first time */
1624    else mcsiminfo_name=NULL;
1625    char *datfile=NULL;
1626    if (!datfile) datfile = str_cat("mcstas.", "dat", NULL);
1627    strcat(mcformat, " scan ");
1628    if (mcverbose)
1629      printf("Writing scan file=%s (%s) into directory %s: %s=%g:%g\n",
1630        datfile, mcsiminfo_name, mcdirname,
1631        Files_to_Merge[scan_index1].mcinputtable[ipar_var].name, ipar_min, ipar_max);
1632    strcpy(mcinstrument_name,   Files_to_Merge[scan_index1].InstrName);
1633    strcpy(mcinstrument_source, Files_to_Merge[scan_index1].Source);
1634    if (!mcdisable_output_files)
1635    mcdetector_out_2D(title,
1636      Files_to_Merge[scan_index1].mcinputtable[ipar_var].name,
1637      Files_to_Merge[scan_index1].ylabel,
1638      ipar_min, ipar_max, 0, 0,
1639      -scan_length1, Files_to_Merge[scan_index1].mcnumipar+2*mon_count,
1640      NULL, Scan.data, NULL, datfile,
1641      Files_to_Merge[scan_index1].component,
1642      Files_to_Merge[scan_index1].POSITION);
1643
1644    /* for PGPLOT: close mcstas.sim */
1645    if (strstr(mcformat, "McCode")) mcsiminfo_close();
1646    Table_Free(&Scan);
1647    memfree(datfile);
1648    memfree(youts);
1649    memfree(title);
1650    /* go to end of scan and continue to search for scans */
1651  } /* for scan_index1 */
1652} /*  mcformat_scan_compare */
1653
1654/*******************************************************************************
1655* mcformat_merge_output: output non empty files
1656*                        special handling for multiarray/scans
1657*******************************************************************************/
1658int mcformat_merge_output(int nb)
1659{
1660  int i;
1661  char mcdisable_output_files_sav=mcdisable_output_files;
1662  if (mcscanmode == 2) mcdisable_output_files=1; /* scan only will skip writing for non scan files */
1663  /* output files for non empty elements */
1664  for (i=0; i<nb; i++) {
1665    if (mcformat_output(Files_to_Merge[i])) {
1666      if (mcverbose) {
1667        printf("mcformat: merging/scanning %s ", Files_to_Merge[i].outputname);
1668        printf("into %s%s ",
1669          Files_to_Merge[i].mcdirname ? Files_to_Merge[i].mcdirname : ".",
1670          Files_to_Merge[i].mcdirname ? MC_PATHSEP_S : "");
1671        if (mcdisable_output_files) printf(" (--test or --scan-only mode)\n"); else printf("\n");
1672      }
1673    }
1674  }
1675  mcdisable_output_files = mcdisable_output_files_sav;
1676
1677  if (mcscanmode) {
1678  /* build and output scans (if any) for non empty elements of same index */
1679  mcformat_scan_compare(nb);
1680  /* allocate nb_scan rows of nb_vars+(nb_monitors*3) columns */
1681  /* fill in scan variables extracted from mcinputtable */
1682  /* fill in monitor values */
1683  /* call outout routine for scans (can not be binary format): from mcrun.pl */
1684  }
1685
1686  /* free merged/scan elements */
1687  for (i=0; i<nb; i++) mcformat_free_mcstas_struct(Files_to_Merge[i]);
1688  return(nb);
1689}
1690
1691/*******************************************************************************
1692* mcformat_usage: print mcformat usage/help
1693*******************************************************************************/
1694void mcformat_usage(char *pgmname)
1695{
1696  int i;
1697
1698  fprintf(stderr, "mcformat version %s format conversion tool (" MCCODE_STRING ")\n", MCFORMAT);
1699  fprintf(stderr, "Usage: %s [options] file1|dir1 file2|dir2 ...\n", pgmname);
1700  fprintf(stderr,
1701"Convert/merge files and directories from McStas format to an other specified format\n"
1702"Options are:\n"
1703"  -d DIR    --dir=DIR        Put all data files in directory DIR.\n"
1704"  -f FILE   --file=FILE      Put all data in a single file.\n"
1705"  --no-output-files          Do not write any data files (test data).\n"
1706"  -h        --help           Show this help message.\n"
1707"  -p FORMAT --format=FORMAT  Output data files using format FORMAT\n"
1708"  -c        --force          Force writting in existing directories\n"
1709"  -t        --test           Test mode, does not write files\n"
1710"  -m        --merge          Add/Append equivalent data files and lists\n"
1711"            --merge-samedir  Merges inside same directories (dangerous)\n"
1712"  -s        --scan           Gather simulations per scan series\n"
1713"  -so       --scan-only      Create scan series but does not merge data\n"
1714"            --verbose        Verbose mode\n"
1715"\n"
1716"Examples:\n"
1717"mcformat -d target_dir original_dir\n"
1718"  # using default target format %s\n"
1719"mcformat -d target_dir original_dir --format=target_format\n"
1720"mcformat -d target_dir original_dir1 original_dir2 --merge\n"
1721"  # merge simulation data sets (grids)\n"
1722"mcformat -d target_dir dir1 dir2 ... --scan\n"
1723"  # merge scan data sets\n\n",
1724getenv(FLAVOR_UPPER "_FORMAT") ? getenv(FLAVOR_UPPER "_FORMAT") : "MCCODE");
1725
1726  fprintf(stderr, "Available output formats are (default is %s):\n  ", mcformat);
1727  fprintf(stderr, "\n  Format modifiers: FORMAT may be followed by 'binary float' or \n");
1728  fprintf(stderr, "  'binary double' to save data blocks as binary. This removes text headers.\n");
1729  fprintf(stderr, "  The " FLAVOR_UPPER "_FORMAT environment variable may set the default FORMAT to use.\n");
1730} /* mcformat_usage */
1731
1732/*******************************************************************************
1733* mcformat_parseoptions: parse command line parameters
1734*******************************************************************************/
1735void
1736mcformat_parseoptions(int argc, char *argv[])
1737{
1738  int i;
1739  char cwd[CHAR_BUF_LENGTH];
1740  mcdirname = NULL;
1741  for(i = 1; i < argc; i++)
1742  {
1743    if(!strcmp("-d", argv[i]) && (i + 1) < argc)
1744      mcformat_usedir(argv[++i]);
1745    else if(!strncmp("-d", argv[i], 2))
1746      mcformat_usedir(&argv[i][2]);
1747    else if(!strcmp("--dir", argv[i]) && (i + 1) < argc)
1748      mcformat_usedir(argv[++i]);
1749    else if(!strncmp("--dir=", argv[i], 6))
1750      mcformat_usedir(&argv[i][6]);
1751    else if(!strcmp("-h", argv[i]) || !strcmp("--help", argv[i]))
1752    {  mcformat_usage(argv[0]); exit(-1); }
1753    else if(!strcmp("-v", argv[i]) || !strcmp("--version", argv[i]))
1754    {  fprintf(stderr, "mcformat version %s format conversion tool (" MCCODE_STRING ")\n", MCFORMAT); exit(-1); }
1755    else if(!strcmp("-c", argv[i]))
1756      mcforcemode = 1;
1757    else if(!strcmp("--force", argv[i]))
1758      mcforcemode = 1;
1759    else if(!strncmp("--format=", argv[i], 9)) {
1760      mcformat=&argv[i][9];
1761    }
1762    else if(!strcmp("--format", argv[i]) && (i + 1) < argc) {
1763      mcformat=argv[++i];
1764    }
1765    else if(!strcmp("--no-output-files", argv[i]))
1766      mcdisable_output_files = 1;
1767    else if(!strcmp("--verbose", argv[i]))
1768      mcverbose  = 1;
1769    else if(!strcmp("-t", argv[i]) || !strcmp("--test", argv[i]))
1770      mcdisable_output_files = 1;
1771    else if(!strcmp("-m", argv[i]) || !strcmp("--merge", argv[i]))
1772      mcmergemode= 1;
1773    else if(!strcmp("-s", argv[i]) || !strcmp("--scan", argv[i]))
1774      mcscanmode = 1;
1775    else if(!strcmp("-so", argv[i]) || !strcmp("--scan-only", argv[i]))
1776      mcscanmode = 2;
1777    else if(!strcmp("--merge-samedir", argv[i]))
1778      mcmergesamedir=mcmergemode=1;
1779    else {
1780      /* convert argv[i]: store index of argument */
1781      if (files_to_convert_NB <CHAR_BUF_LENGTH)
1782        files_to_convert_Array[files_to_convert_NB++] = i;
1783      else
1784        fprintf(stderr, "Warning: Exceeding maximum number of files to process (%d).\n"
1785          "Ignoring %s [mcformat:mcformat_parseoptions].\n",
1786          CHAR_BUF_LENGTH, argv[i]);
1787    }
1788  }
1789
1790  if (!mcdirname) {
1791    getcwd(cwd, CHAR_BUF_LENGTH);
1792    mcdirname = str_dup(cwd); /* default is to export to PWD */
1793  }
1794} /* mcformat_parseoptions */
1795
1796/*******************************************************************************
1797* main: program entry point (start):calls parseoptions, convert and merge
1798*******************************************************************************/
1799int main(int argc, char *argv[])
1800{
1801  time_t t;
1802  int j;
1803
1804#ifdef MAC
1805  argc = ccommand(&argv);
1806#endif
1807
1808  t = (time_t)mcstartdate;
1809  time(&t);
1810  mcstartdate = t;
1811  mcformat=getenv(FLAVOR_UPPER "_FORMAT") ?
1812           getenv(FLAVOR_UPPER "_FORMAT") : "MCCODE";
1813  /* default is to output as McCode format */
1814
1815  /* parse parameters from the command line and get files to convert */
1816  mcformat_parseoptions(argc, argv);
1817
1818  if (!files_to_convert_NB) {  mcformat_usage(argv[0]); exit(-1); }
1819
1820  /* check format */
1821  if (!mcformat || !strlen(mcformat)
1822   || !strcasecmp(mcformat, "McStas") || !strcasecmp(mcformat, "McXtrace")
1823   || !strcasecmp(mcformat, "PGPLOT"))
1824    mcformat = "MCCODE";
1825
1826  mcnbconvert = mcdircount = 0;
1827  mcoutputdir = mcdirname; /* base output dir */
1828
1829  /* count the number of files to store */
1830  for(j = 0; j < files_to_convert_NB; j++) {
1831    int    i;
1832    i = files_to_convert_Array[j]; /* does the job for each file/dir */
1833    mcnbconvert += mcformat_count(argv[i]);
1834  }
1835  if (mcdircount) { /* allocate list of directories amd FILE handles */
1836    mcdirnames  = (char**)mem(mcdircount*sizeof(char*));
1837    mcsimfiles  = (FILE**)mem(mcdircount*sizeof(FILE*));
1838    mcinstrnames= (char**)mem(mcdircount*sizeof(char*));
1839    mcsources   = (char**)mem(mcdircount*sizeof(char*));
1840  }
1841  if (mcmergemode || mcscanmode) {
1842    printf("mcformat: Will process and merge/scan %d data file%s in %d director%s.\n",
1843      mcnbconvert, mcnbconvert > 1 ? "s" : "", mcdircount, mcdircount > 1 ? "ies" : "y");
1844    /* allocate array of Original file structures */
1845    Files_to_Merge = mem(mcnbconvert*sizeof(struct McStas_file_format));
1846    Scans_to_merge = mem(mcnbconvert*sizeof(int));
1847    for (j=0; j<mcnbconvert; j++) {
1848      Files_to_Merge[j] = mcformat_init_mcstas_struct();
1849      Scans_to_merge[j] = -1;
1850    }
1851
1852  }
1853  mcnbconvert = mcdircount = 0;
1854
1855  strcpy(mcinstrument_source, "");
1856  strcpy(mcinstrument_name  , "");
1857
1858  /* calls the conversion routine. SIM file opened at first conversion in mcformat_convert */
1859  for(j = 0; j < files_to_convert_NB; j++) {
1860    int    i;
1861    i = files_to_convert_Array[j]; /* does the job for each file/dir */
1862    mcformat_convert(argv[i]);
1863  }
1864
1865  if (mcmergemode || mcscanmode) {
1866    /* call to merging routine (free some Files_to_Merge elements) */
1867    mcformat_merge_compare(mcnbconvert);
1868    /* iterative call to output routine for remaining elements of Files_to_Merge */
1869    mcformat_merge_output(mcnbconvert);
1870    if (Files_to_Merge) memfree(Files_to_Merge);
1871    if (Scans_to_merge) memfree(Scans_to_merge);
1872  }
1873
1874  /* clean up: close all SIM files */
1875  if (!mcdisable_output_files) {
1876    if (!mcdircount) mcsiminfo_close();
1877    else {
1878      for (j=0; j<mcdircount; j++) {
1879        if (mcsimfiles[j]) {
1880          mcdirname = mcdirnames[j];
1881          mcsiminfo_file     = mcsimfiles[j];
1882          strncpy(mcinstrument_source, str_last_word(mcinstrnames[j]), CHAR_BUF_LENGTH);
1883          strncpy(mcinstrument_name  , str_last_word(mcsources[j]), CHAR_BUF_LENGTH);
1884          mcsiminfo_close();
1885          mcsimfiles[j] = NULL;
1886          memfree(mcinstrnames[j]);
1887          memfree(mcsources[j]);
1888          memfree(mcdirnames[j]);
1889        }
1890      }
1891      memfree(mcdirnames);
1892      memfree(mcsimfiles);
1893      memfree(mcinstrnames);
1894      memfree(mcsources);
1895    }
1896  }
1897
1898  return mcnbconvert; /* number of converted files */
1899
1900} /* main */
1901