add bayesModel

utils/Makefile.am

-bin_PROGRAMS = bayesChires bayesImage bayesSource bayesBest bayesDlsDs bayesGrad bayesMap bayesAmp bayesM200 bayesCl bayesCleanlens bayesClump
+bin_PROGRAMS = bayesChires bayesImage bayesSource bayesBest bayesDlsDs bayesGrad bayesMap bayesAmp bayesM200 bayesCl bayesCleanlens bayesClump bayesModel
 if PGPLOT
 bin_PROGRAMS += Histogram Histogram2D 
 endif
@@ -18,6 +18,7 @@ bayesM200_SOURCES = bayesM200.c
 bayesCl_SOURCES = bayesCl.c
 bayesCleanlens_SOURCES = bayesCleanlens.c
 bayesClump_SOURCES = bayesClump.c
+bayesModel_SOURCES = bayesModel.c
 
 
 bayesChires_LDFLAGS = @CFITSIO_LIB@ @WCS_LIB@ @GSL_LIB@
@@ -68,6 +69,10 @@ bayesClump_LDFLAGS = @CFITSIO_LIB@ @WCS_LIB@ @GSL_LIB@
 bayesClump_LDADD = @top_builddir@/src/liblenstool.a \
 		    @top_builddir@/liblt/liblt.a @CFITSIO_STATIC@ @WCS_STATIC@ -lm @GSL_STATIC@
 
+bayesModel_LDFLAGS = @CFITSIO_LIB@ @WCS_LIB@ @GSL_LIB@
+bayesModel_LDADD = @top_builddir@/src/liblenstool.a \
+		    @top_builddir@/liblt/liblt.a @CFITSIO_STATIC@ @WCS_STATIC@ -lm @GSL_STATIC@
+
 if PGPLOT
 Histogram_SOURCES = histogram.f 
 

utils/bayesModel.c

+#include<stdio.h>
+#include<signal.h>
+#include<stdlib.h>
+#include<string.h>
+#include<time.h>
+#include "dimension.h"
+#include "structure.h"
+#include "constant.h" 
+#include "fonction.h"
+#include "bayesChires.h"
+
+/******************************************************************************
+ * Name: bayesModel.c
+ * Authors: JR
+ * Date: 20/10/19
+ * 
+ * Analyse a bayes.dat file and compute the best.par and bestopt.par files
+ * from the Nth model 
+ *
+ * syntax : bayesModel <.par> 
+ *
+ * The bayes.dat file in the current directory is used.
+ * Output : modelXX.par
+ *
+ ******************************************************************************/
+typedef void (*sighandler_t)(int);
+int optInterrupt;
+
+void help_msg()
+{
+		fprintf(stderr, "Syntax : bayesModel [OPTION] <.par> <modelnumber>\n");
+		exit(-1);
+}
+
+void writeSource()
+{
+
+}
+
+
+int main( int argc, char** argv )
+{
+	int    method;
+	int    i;
+	int    prior;
+	double limit;
+	char   *pch;
+	FILE   *IN;
+    extern int sblock[NFMAX][NPAMAX];
+    extern struct g_source  S;
+    extern struct galaxie  sources[NFMAX];
+    int ipx;
+    double evidence;
+
+    double **array; // contains the bayes.dat data
+    int nParam;
+    long int iVal, nVal;  // size of array
+    int flagswap;
+    char cmd[500];
+
+        flagswap=0;
+	// Check the arguments
+	if( argc < 3 )
+		help_msg();
+
+	if( strstr(argv[1], ".par") == NULL )
+		help_msg();
+
+        sscanf(argv[2],"%ld",&iVal);
+
+	// Read the .par file
+	init_grille( argv[1], 1);
+	
+	// Read constraints
+	readConstraints();
+
+	if( G.nmsgrid != G.nlens ) 
+	{
+		prep_non_param();
+	}
+	
+        // Read the bayes.dat file
+        array = readBayesModels(&nParam, &nVal);
+        if( array == NULL )
+        {
+                fprintf(stderr, "ERROR: bayes.dat file not found\n");
+                return -1;
+        }
+        if( iVal>=nVal)
+        {
+                fprintf(stderr, "ERROR: model index %ld larger than number of bayes.dat lines (%ld)\n",iVal,nVal);
+                return -1;
+        }
+
+	// Set the model in memory
+        setBayesModel( iVal, nVal, array );
+
+	// Reset the lens parameters
+	set_res_par();
+
+	// Print the best.par and bestopt.par files
+	// save the original best.par and bestopt.par
+	if( (IN =  fopen( "best.par" ,"r") ))
+	{
+		fclose(IN);
+		sprintf(cmd,"mv best.par tmp.par\n");
+		system(cmd);
+		sprintf(cmd,"mv bestopt.par tmpopt.par\n");
+		system(cmd);
+                flagswap=1;
+	}
+
+        o_print_res(o_chi(), evidence);  // maybe we know the evidence 
+	sprintf(cmd,"cp best.par model%ld.par",iVal);
+	system(cmd);
+	sprintf(cmd,"cp bestopt.par modelopt%ld.par",iVal);
+	system(cmd);
+        
+
+	if( flagswap==1)
+	{
+		sprintf(cmd,"mv tmp.par best.par");
+		system(cmd);
+		sprintf(cmd,"mv tmpopt.par bestopt.par");
+		system(cmd);
+                flagswap=1;
+	}
+
+    o_global_free();
+
+	return 0;
+}
+