Index: /trunk/Ohana/src/dvomerge/src/dvoverify_catalogs.c =================================================================== --- /trunk/Ohana/src/dvomerge/src/dvoverify_catalogs.c (revision 36083) +++ /trunk/Ohana/src/dvomerge/src/dvoverify_catalogs.c (revision 36084) @@ -143,5 +143,10 @@ continue; } - fscanf (results, "%*s %d %*s %d", &NbadHost, &NNotSortedHost); + int Nresults = fscanf (results, "%*s %d %*s %d", &NbadHost, &NNotSortedHost); + if (Nresults != 2) { + fprintf (stderr, "warning: failure scanning results\n"); + fclose (results); + continue; + } *Nbad += NbadHost; NNotSorted += NNotSortedHost; Index: /trunk/Ohana/src/getstar/src/dvoImagesAtCoords.c =================================================================== --- /trunk/Ohana/src/getstar/src/dvoImagesAtCoords.c (revision 36083) +++ /trunk/Ohana/src/getstar/src/dvoImagesAtCoords.c (revision 36084) @@ -16,5 +16,4 @@ int status; Image *dbImages; - int Npoints; FITS_DB db; @@ -25,4 +24,5 @@ Point *points; int readStdin = 0; + int Npoints = 0; if (coordsFile) { if (strcmp(coordsFile, "-") == 0) { Index: /trunk/Ohana/src/kapa2/include/constants.h =================================================================== --- /trunk/Ohana/src/kapa2/include/constants.h (revision 36083) +++ /trunk/Ohana/src/kapa2/include/constants.h (revision 36084) @@ -15,5 +15,5 @@ // XXX for the moment, this is set to match the values in SetColorScale3D_CC # define NPIXELS_STATIC 128 -// # define NPIXELS_STATIC 1000 +// # define NPIXELS_STATIC 4600 -- if we want fullcolor to work, need this # define PAD1 3 Index: /trunk/Ohana/src/libdvo/include/dvo.h =================================================================== --- /trunk/Ohana/src/libdvo/include/dvo.h (revision 36083) +++ /trunk/Ohana/src/libdvo/include/dvo.h (revision 36084) @@ -696,4 +696,5 @@ void sort_coords_index (double *X, double *Y, off_t *S, off_t N); void sort_coords_indexonly (double *X, double *Y, off_t *S, off_t N); +void sort_IDs_indexonly (opihi_int *X, off_t *S, off_t N); void sort_regions (SkyRegion *region, off_t N); Index: /trunk/Ohana/src/libdvo/include/dvodb.h =================================================================== --- /trunk/Ohana/src/libdvo/include/dvodb.h (revision 36083) +++ /trunk/Ohana/src/libdvo/include/dvodb.h (revision 36084) @@ -174,4 +174,7 @@ AVE_TMEAN, AVE_TRANGE, + AVE_PSF_QF, + AVE_PSF_QF_PERF, + AVE_STARGAL, AVE_Xp, AVE_NMEAS, Index: /trunk/Ohana/src/libdvo/src/dbExtractAverages.c =================================================================== --- /trunk/Ohana/src/libdvo/src/dbExtractAverages.c (revision 36083) +++ /trunk/Ohana/src/libdvo/src/dbExtractAverages.c (revision 36084) @@ -151,4 +151,14 @@ case AVE_TRANGE: value.Flt = TimeValue (average[0].Trange, 0, TimeFormat); + break; + + case AVE_PSF_QF: + value.Flt = average[0].psfQF; + break; + case AVE_PSF_QF_PERF: + value.Flt = average[0].psfQFperf; + break; + case AVE_STARGAL: + value.Flt = average[0].stargal; break; Index: /trunk/Ohana/src/libdvo/src/dbFields.c =================================================================== --- /trunk/Ohana/src/libdvo/src/dbFields.c (revision 36083) +++ /trunk/Ohana/src/libdvo/src/dbFields.c (revision 36084) @@ -314,34 +314,39 @@ } - if (!strcasecmp (fieldName, "RA")) ESCAPE (AVE_RA, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "DEC")) ESCAPE (AVE_DEC, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "GLON")) ESCAPE (AVE_GLON, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "GLAT")) ESCAPE (AVE_GLAT, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "ELON")) ESCAPE (AVE_ELON, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "ELAT")) ESCAPE (AVE_ELAT, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "dRA")) ESCAPE (AVE_RA_ERR, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "dDEC")) ESCAPE (AVE_DEC_ERR, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "uRA")) ESCAPE (AVE_U_RA, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "uDEC")) ESCAPE (AVE_U_DEC, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "duRA")) ESCAPE (AVE_U_RA_ERR, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "duDEC")) ESCAPE (AVE_U_DEC_ERR, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "PAR")) ESCAPE (AVE_PAR, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "dPAR")) ESCAPE (AVE_PAR_ERR, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "ChiSqPos")) ESCAPE (AVE_CHISQ_POS, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "ChiSqPM")) ESCAPE (AVE_CHISQ_PM, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "ChiSqPar")) ESCAPE (AVE_CHISQ_PAR, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "Tmean")) ESCAPE (AVE_TMEAN, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "Trange")) ESCAPE (AVE_TRANGE, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "NMEAS")) ESCAPE (AVE_NMEAS, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "NMISS")) ESCAPE (AVE_NMISS, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "NPOS")) ESCAPE (AVE_NPOS, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "NASTROM")) ESCAPE (AVE_NPOS, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "FLAGS")) ESCAPE (AVE_OBJ_FLAGS, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "OBJ_FLAGS")) ESCAPE (AVE_OBJ_FLAGS, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "OBJFLAGS")) ESCAPE (AVE_OBJ_FLAGS, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "OBJID")) ESCAPE (AVE_OBJID, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "CATID")) ESCAPE (AVE_CATID, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "EXTID_HI")) ESCAPE (AVE_EXTID_HI, MAG_NONE, OPIHI_INT); - if (!strcasecmp (fieldName, "EXTID_LO")) ESCAPE (AVE_EXTID_LO, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "RA")) ESCAPE (AVE_RA, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "DEC")) ESCAPE (AVE_DEC, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "GLON")) ESCAPE (AVE_GLON, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "GLAT")) ESCAPE (AVE_GLAT, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "ELON")) ESCAPE (AVE_ELON, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "ELAT")) ESCAPE (AVE_ELAT, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "dRA")) ESCAPE (AVE_RA_ERR, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "dDEC")) ESCAPE (AVE_DEC_ERR, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "uRA")) ESCAPE (AVE_U_RA, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "uDEC")) ESCAPE (AVE_U_DEC, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "duRA")) ESCAPE (AVE_U_RA_ERR, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "duDEC")) ESCAPE (AVE_U_DEC_ERR, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "PAR")) ESCAPE (AVE_PAR, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "dPAR")) ESCAPE (AVE_PAR_ERR, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "ChiSqPos")) ESCAPE (AVE_CHISQ_POS, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "ChiSqPM")) ESCAPE (AVE_CHISQ_PM, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "ChiSqPar")) ESCAPE (AVE_CHISQ_PAR, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "Tmean")) ESCAPE (AVE_TMEAN, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "Trange")) ESCAPE (AVE_TRANGE, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "psfqf")) ESCAPE (AVE_PSF_QF, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "psfqfperf")) ESCAPE (AVE_PSF_QF_PERF, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "psf_qf")) ESCAPE (AVE_PSF_QF, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "psf_qf_perf")) ESCAPE (AVE_PSF_QF_PERF, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "stargal")) ESCAPE (AVE_STARGAL, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "NMEAS")) ESCAPE (AVE_NMEAS, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "NMISS")) ESCAPE (AVE_NMISS, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "NPOS")) ESCAPE (AVE_NPOS, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "NASTROM")) ESCAPE (AVE_NPOS, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "FLAGS")) ESCAPE (AVE_OBJ_FLAGS, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "OBJ_FLAGS")) ESCAPE (AVE_OBJ_FLAGS, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "OBJFLAGS")) ESCAPE (AVE_OBJ_FLAGS, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "OBJID")) ESCAPE (AVE_OBJID, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "CATID")) ESCAPE (AVE_CATID, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "EXTID_HI")) ESCAPE (AVE_EXTID_HI, MAG_NONE, OPIHI_INT); + if (!strcasecmp (fieldName, "EXTID_LO")) ESCAPE (AVE_EXTID_LO, MAG_NONE, OPIHI_INT); // check for code:mode in photcode name Index: /trunk/Ohana/src/libdvo/src/dvosorts.c =================================================================== --- /trunk/Ohana/src/libdvo/src/dvosorts.c (revision 36083) +++ /trunk/Ohana/src/libdvo/src/dvosorts.c (revision 36084) @@ -60,4 +60,19 @@ /* sort a coordinate pair (X,Y) and the associated index (S) */ +void sort_IDs_indexonly (opihi_int *X, off_t *S, off_t N) { + +# define SWAPFUNC(A,B){ off_t itmp; \ + itmp = S[A]; S[A] = S[B]; S[B] = itmp; \ +} +# define COMPARE(A,B)(X[S[A]] < X[S[B]]) + + OHANA_SORT (N, COMPARE, SWAPFUNC); + +# undef SWAPFUNC +# undef COMPARE + +} + +/* sort a coordinate pair (X,Y) and the associated index (S) */ void sort_coords_indexonly (double *X, double *Y, off_t *S, off_t N) { Index: /trunk/Ohana/src/libkapa/src/IOfuncs.c =================================================================== --- /trunk/Ohana/src/libkapa/src/IOfuncs.c (revision 36083) +++ /trunk/Ohana/src/libkapa/src/IOfuncs.c (revision 36084) @@ -54,5 +54,8 @@ va_end (argp); - if (!Nbyte) return (FALSE); + if (!Nbyte) { + KiiSendCommand (device, 16, "LEN: %11d", 0); + return (FALSE); + } /* the message may contain up to 99,999,999,999 bytes (100MB) */ @@ -84,6 +87,11 @@ /* find the message length, allocate space */ sscanf (buffer, "%*s %d", &Nbytes); + if (Nbytes == 0) { + return TRUE; + } + ALLOCATE (message, char, Nbytes + 1); memset (message, 0, Nbytes + 1); + /* read Nbytes from the device */ Index: /trunk/Ohana/src/libohana/src/gaussj.c =================================================================== --- /trunk/Ohana/src/libohana/src/gaussj.c (revision 36083) +++ /trunk/Ohana/src/libohana/src/gaussj.c (revision 36084) @@ -1,5 +1,5 @@ # include # define GROWTHTEST 0 -# define MAX_RANGE 1.0e7 +# define MAX_RANGE 1.0e9 // Gauss-Jordan elimination using full pivots based on Press et al's description. Substantially @@ -91,4 +91,5 @@ /* rescale by pivot reciprocal */ double tmpval = 1.0 / A[maxcol][maxcol]; + // XXX why is this here (don't I double count this element?) A[maxcol][maxcol] = 1.0; A[maxcol][maxcol] = 1.0; for (col = 0; col < N; col++) A[maxcol][col] *= tmpval; @@ -117,4 +118,15 @@ } } + +# if (GROWTHTEST) + fprintf (stderr, "final A matrix:\n"); + for (row = 0; row < N; row++) { + for (col = 0; col < N; col++) { + fprintf (stderr, "%10.3e ", A[row][col]); + } + fprintf (stderr, "\n"); + } + fprintf (stderr, "\n"); +# endif // swap back the inverse matrix based on the row swaps above Index: /trunk/Ohana/src/opihi/cmd.astro/Makefile =================================================================== --- /trunk/Ohana/src/opihi/cmd.astro/Makefile (revision 36083) +++ /trunk/Ohana/src/opihi/cmd.astro/Makefile (revision 36084) @@ -60,4 +60,5 @@ $(SRC)/scale.$(ARCH).o \ $(SRC)/sexigesimal.$(ARCH).o \ +$(SRC)/sersic.$(ARCH).o \ $(SRC)/spec.$(ARCH).o \ $(SRC)/specpairfit.$(ARCH).o \ Index: /trunk/Ohana/src/opihi/cmd.astro/init.c =================================================================== --- /trunk/Ohana/src/opihi/cmd.astro/init.c (revision 36083) +++ /trunk/Ohana/src/opihi/cmd.astro/init.c (revision 36084) @@ -49,4 +49,5 @@ int scale PROTO((int, char **)); int sexigesimal PROTO((int, char **)); +int sersic PROTO((int, char **)); int spec PROTO((int, char **)); int specpairfit PROTO((int, char **)); @@ -107,4 +108,5 @@ {1, "scale", scale, "get / set real bzero / bscale values"}, {1, "sexigesimal", sexigesimal, "convert to/from sexigesimal/decimal"}, + {1, "sersic", sersic, "generate sub-pixel resolved sersic model"}, {1, "spec", spec, "extract a spectrum"}, {1, "specpairfit", specpairfit, "fit spectrum to another spectrum"}, Index: /trunk/Ohana/src/opihi/cmd.astro/sersic.c =================================================================== --- /trunk/Ohana/src/opihi/cmd.astro/sersic.c (revision 36084) +++ /trunk/Ohana/src/opihi/cmd.astro/sersic.c (revision 36084) @@ -0,0 +1,115 @@ +# include "astro.h" + +float sersic_norm (float index); + +// the sersic model has an integral of 2 \pi Rmin Rmax Io norm * Nsubpix^2 + +int sersic (int argc, char **argv) { + + int i, j; + Buffer *buf; + + // generate a sersic profile in the given buffer using N + + if (argc != 10) { + gprint (GP_ERR, "USAGE: sersic (buffer) (Io) (Xo) (Yo) (Rmajor) (Rminor) (theta) (index) (Nsubpix)\n"); + gprint (GP_ERR, " generate a sersic profile for the given parameters\n"); + return (FALSE); + } + + /* select input / output buffers */ + if ((buf = SelectBuffer (argv[1], OLDBUFFER, TRUE)) == NULL) return (FALSE); + int Nx = buf[0].header.Naxis[0]; + int Ny = buf[0].header.Naxis[1]; + + // Io is the central pixel normalization + float Io = atof(argv[2]); + + // Xo, Yo are the central pixel coordinates in real pixels + float Xo = atof(argv[3]); + float Yo = atof(argv[4]); + + float Rmaj = atof(argv[5]); + float Rmin = atof(argv[6]); + + float theta = atof(argv[7]); + float index = atof(argv[8]); + float rindex = 0.5 / index; + + float Nsubpix = atof(argv[9]); + + float f1 = 1.0 / SQ(Rmin) + 1.0 / SQ(Rmaj); + float f2 = 1.0 / SQ(Rmin) - 1.0 / SQ(Rmaj); + + float sxr = 0.5*f1 - 0.5*f2*cos(2.0*theta*RAD_DEG); + float syr = 0.5*f1 + 0.5*f2*cos(2.0*theta*RAD_DEG); + + // sxr, syr cannot be < 0 (f1 >= f2) + + float Rxx = +1.0 / sqrt(sxr); + float Ryy = +1.0 / sqrt(syr); + float Rxy = -f2*sin(2.0*theta*RAD_DEG); + + float kappa = -0.275552 + 1.972625*index + 0.003487 * SQ(index); + + // dereference pointer + float *in = (float *) buf[0].matrix.buffer; + + for (j = 0; j < Ny; j++) { + for (i = 0; i < Nx; i++) { + + float x = ((i / Nsubpix) - Xo); + float y = ((j / Nsubpix) - Yo); + + // z = r^2 + float z = SQ(x/Rxx) + SQ(y/Ryy) + Rxy*x*y; + + float q = pow (z, rindex); + float f = Io * exp(-kappa*q); + + in[i + j*Nx] = f; + } + } + + float norm = sersic_norm(index); + set_variable ("sersic_norm", norm); + set_variable ("sersic_kappa", kappa); + + return (TRUE); +} + +float sersic_norm (float index) { + + float lnorm, norm; + + if ((index >= 0.0) && (index < 1.0)) { + lnorm = 0.201545 - 0.950965 * index - 0.315248 * SQ(index); + norm = exp(lnorm); + return norm; + } + + if ((index >= 1.0) && (index < 2.0)) { + lnorm = 0.402084 - 1.357775 * index - 0.105102 * SQ(index); + norm = exp(lnorm); + return norm; + } + + if ((index >= 2.0) && (index < 3.0)) { + lnorm = 0.619093 - 1.591674 * index - 0.041576 * SQ(index); + norm = exp(lnorm); + return norm; + } + + if ((index >= 3.0) && (index < 4.0)) { + lnorm = 0.770263 - 1.696421 * index - 0.023363 * SQ(index); + norm = exp(lnorm); + return norm; + } + + if ((index >= 4.0) && (index < 5.5)) { + lnorm = 0.885891 - 1.755684 * index - 0.015753 * SQ(index); + norm = exp(lnorm); + return norm; + } + return NAN; +} Index: /trunk/Ohana/src/opihi/cmd.basic/list.c =================================================================== --- /trunk/Ohana/src/opihi/cmd.basic/list.c (revision 36083) +++ /trunk/Ohana/src/opihi/cmd.basic/list.c (revision 36084) @@ -2,4 +2,6 @@ # define D_NLINES 100 static char prompt[] = ">> "; + +static int set_list_varname (char *line, char *base, int N, int excelStyle); int list (int argc, char **argv) { @@ -21,4 +23,14 @@ } + int EXCEL_STYLE = FALSE; + if ((N = get_argument (argc, argv, "-excel-style"))) { + remove_argument (N, &argc, argv); + EXCEL_STYLE = TRUE; + } + if ((N = get_argument (argc, argv, "-excel"))) { + remove_argument (N, &argc, argv); + EXCEL_STYLE = TRUE; + } + if ((N = get_argument (argc, argv, "-vectors"))) { remove_argument (N, &argc, argv); @@ -49,5 +61,5 @@ for (i = 0; i < argc - 3; i++) { - sprintf (line, "%s:%d", argv[1], i); + set_list_varname (line, argv[1], i, EXCEL_STYLE); set_str_variable (line, argv[i+3]); } @@ -83,5 +95,7 @@ if (!word) break; - sprintf (line, "%s:%d", argv[1], nWords); + // sprintf (line, "%s:%d", argv[1], nWords); + set_list_varname (line, argv[1], nWords, EXCEL_STYLE); + set_str_variable (line, word); FREE (word); @@ -120,5 +134,6 @@ sprintf (line, "%s:%d", argv[3], i); value = get_variable (line); - sprintf (line, "%s:%d", argv[1], i); + // sprintf (line, "%s:%d", argv[1], i); + set_list_varname (line, argv[1], i, EXCEL_STYLE); set_str_variable (line, value); } @@ -136,5 +151,6 @@ N = get_int_variable (line, &found); for (i = 0; i < argc - 3; i++) { - sprintf (line, "%s:%d", argv[1], N + i); + // sprintf (line, "%s:%d", argv[1], N + i); + set_list_varname (line, argv[1], N + i, EXCEL_STYLE); set_str_variable (line, argv[i+3]); } @@ -159,5 +175,6 @@ N = get_int_variable (line, &found); for (i = 0; i < N; i++) { - sprintf (line, "%s:%d", argv[1], i); + // sprintf (line, "%s:%d", argv[1], i); + set_list_varname (line, argv[1], i, EXCEL_STYLE); value = get_variable (line); if (value == NULL) continue; @@ -165,5 +182,6 @@ free (value); for (j = i + 1; j < N; j++) { - sprintf (line2, "%s:%d", argv[1], j); + // sprintf (line2, "%s:%d", argv[1], j); + set_list_varname (line2, argv[1], j, EXCEL_STYLE); next_value = get_variable (line2); set_str_variable (line, next_value); @@ -235,5 +253,6 @@ if (B != (char *) NULL) { *B = 0; } if (*A != 0) { - sprintf (line, "%s:%d", argv[1], i); + // sprintf (line, "%s:%d", argv[1], i); + set_list_varname (line, argv[1], i, EXCEL_STYLE); set_str_variable (line, A); A = B + 1; @@ -283,10 +302,41 @@ if (*input) { - sprintf (line, "%s:%d", argv[1], i); + // sprintf (line, "%s:%d", argv[1], i); + set_list_varname (line, argv[1], i, EXCEL_STYLE); set_str_variable (line, input); free (input); i++; - } + } } return (TRUE); } + +static int set_list_varname (char *line, char *base, int N, int excelStyle) { + + int i; + + // A-Z correspond to 0 - 25 + + if (excelStyle) { + float f = log(26.0); + float g = (N == 0) ? 0.0 : log(1.0*N); + int Ndigit = (int) (g / f) + 1; + if (Ndigit > 10) { + sprintf (line, "%s:ZZZZZZZZZZ", base); + return FALSE; + } + char name[12]; + memset (name, 0, 12); + for (i = 0; i < Ndigit; i++) { + float Npow = Ndigit - i - 1; + float g = pow(26.0, Npow); + int V = (int) (N / g); + name[i] = (Npow == 0.0) ? 'A' + V : 'A' + V - 1; + N -= V * g; + } + sprintf (line, "%s:%s", base, name); + } else { + sprintf (line, "%s:%d", base, N); + } + return TRUE; +} Index: /trunk/Ohana/src/opihi/cmd.data/Makefile =================================================================== --- /trunk/Ohana/src/opihi/cmd.data/Makefile (revision 36083) +++ /trunk/Ohana/src/opihi/cmd.data/Makefile (revision 36084) @@ -66,4 +66,5 @@ $(SRC)/integrate.$(ARCH).o \ $(SRC)/interpolate.$(ARCH).o \ +$(SRC)/join.$(ARCH).o \ $(SRC)/jpeg.$(ARCH).o \ $(SRC)/kern.$(ARCH).o \ @@ -143,4 +144,5 @@ $(SRC)/vmaxwell.$(ARCH).o \ $(SRC)/vgrid.$(ARCH).o \ +$(SRC)/vlist.$(ARCH).o \ $(SRC)/vload.$(ARCH).o \ $(SRC)/vzload.$(ARCH).o \ Index: /trunk/Ohana/src/opihi/cmd.data/init.c =================================================================== --- /trunk/Ohana/src/opihi/cmd.data/init.c (revision 36083) +++ /trunk/Ohana/src/opihi/cmd.data/init.c (revision 36084) @@ -55,4 +55,5 @@ int integrate PROTO((int, char **)); int interpolate PROTO((int, char **)); +int join PROTO((int, char **)); int jpeg PROTO((int, char **)); int kern PROTO((int, char **)); @@ -133,4 +134,5 @@ int vmaxwell PROTO((int, char **)); int vload PROTO((int, char **)); +int vlist PROTO((int, char **)); int vzload PROTO((int, char **)); int vstats PROTO((int, char **)); @@ -209,4 +211,5 @@ {1, "integrate", integrate, "integrate a vector"}, {1, "interpolate", interpolate, "interpolate between vector pairs"}, + {1, "join", join, "find the join of two ID vectors"}, {1, "jpeg", jpeg, "convert display image to JPEG"}, {1, "kern", kern, "convolve with 3x3 kernel"}, @@ -290,4 +293,5 @@ {1, "vgrid", vgrid, "generate an image from a triplet of vectors"}, {1, "vhistogram", histogram, "generate histogram from vector"}, + {1, "vlist", vlist, "append values to a vector from command line"}, {1, "vload", vload, "load vectors as overlay on image display"}, {1, "vmaxwell", vmaxwell, "fit a Maxwellian to a vector"}, Index: /trunk/Ohana/src/opihi/cmd.data/join.c =================================================================== --- /trunk/Ohana/src/opihi/cmd.data/join.c (revision 36084) +++ /trunk/Ohana/src/opihi/cmd.data/join.c (revision 36084) @@ -0,0 +1,194 @@ +# include "data.h" + +int join_IDs_inner (Vector *ID1, Vector *ID2, Vector *index1, Vector *index2); +int join_IDs_outer (Vector *ID1, Vector *ID2, Vector *index); + +// join (ID1) (ID2) [-index1 (index1)] [-index2 (index2)] +// generate indexes for ID1 and ID2 to identify the matches of ID1 and ID2 (must be int?) +int join (int argc, char **argv) { + + int N, OUTER; + Vector *ID1vec, *ID2vec; + Vector *index1, *index2; + + if ((N = get_argument (argc, argv, "-h"))) goto usage; + if ((N = get_argument (argc, argv, "--help"))) goto usage; + + OUTER = FALSE; + if ((N = get_argument (argc, argv, "-outer"))) { + remove_argument (N, &argc, argv); + OUTER = TRUE; + } + + if ((N = get_argument (argc, argv, "-index1"))) { + remove_argument (N, &argc, argv); + if ((index1 = SelectVector (argv[N], ANYVECTOR, TRUE)) == NULL) return (FALSE); + remove_argument (N, &argc, argv); + } else { + if ((index1 = SelectVector ("index1", ANYVECTOR, TRUE)) == NULL) return (FALSE); + } + + if ((N = get_argument (argc, argv, "-index2"))) { + remove_argument (N, &argc, argv); + if ((index2 = SelectVector (argv[N], ANYVECTOR, TRUE)) == NULL) return (FALSE); + remove_argument (N, &argc, argv); + } else { + if ((index2 = SelectVector ("index2", ANYVECTOR, TRUE)) == NULL) return (FALSE); + } + + if (argc != 3) { + gprint (GP_ERR, "USAGE: join ID1 ID2 [-index1 (index1)] [-index2 (index2)] [-outer]\n"); + gprint (GP_ERR, " use -h or --help for more detail\n"); + } + + if ((ID1vec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE); + if ((ID2vec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE); + + REQUIRE_VECTOR_INT (ID1vec, FALSE); + REQUIRE_VECTOR_INT (ID2vec, FALSE); + + if (OUTER) { + join_IDs_outer (ID1vec, ID2vec, index1); + join_IDs_outer (ID2vec, ID1vec, index2); + } else { + join_IDs_inner (ID1vec, ID2vec, index1, index2); + } + + return (TRUE); + +usage: + gprint (GP_ERR, "we have two modes of operation:\n\n"); + + gprint (GP_ERR, "without -outer, we are finding all matched ID pairs. in this\n"); + gprint (GP_ERR, "case, the two index vectors have the same length, one entry per matched pair.\n"); + gprint (GP_ERR, "ID1[index1] matches to ID2[index2].\n\n"); + + gprint (GP_ERR, "with -outer selected, we are finding the matched element of set 1 for each of set 2\n"); + gprint (GP_ERR, "and vice versa. in this case, index1 is always the same length as ID1, while index2\n"); + gprint (GP_ERR, "is the same length as ID2. ID2[index1] matches ID1 while x1[ID2] matches ID1\n"); + + return FALSE; +} + +// find the entries in ID1 which match ID2 (duplicates allows) +int join_IDs_inner (Vector *ID1, Vector *ID2, Vector *index1, Vector *index2) { + + off_t i, j, first_j, I, J, *N1, *N2, Nmatch, NMATCH, DMATCH; + opihi_int dID; + + NMATCH = MAX(MAX(0.01*ID1->Nelements, 0.01*ID2->Nelements), 100); + DMATCH = NMATCH; + + ResetVector (index1, OPIHI_INT, NMATCH); + ResetVector (index2, OPIHI_INT, NMATCH); + + ALLOCATE (N1, off_t, ID1->Nelements); + ALLOCATE (N2, off_t, ID2->Nelements); + + for (i = 0; i < ID1->Nelements; i++) { N1[i] = i; } + for (i = 0; i < ID2->Nelements; i++) { N2[i] = i; } + + sort_IDs_indexonly (ID1->elements.Int, N1, ID1->Nelements); + sort_IDs_indexonly (ID2->elements.Int, N2, ID2->Nelements); + + Nmatch = 0; + for (i = j = 0; (i < ID1->Nelements) && (j < ID2->Nelements);) { + I = N1[i]; + J = N2[j]; + + dID = ID1->elements.Int[I] - ID2->elements.Int[J]; + + if (dID < 0) { i++; continue; } + if (dID > 0) { j++; continue; } + + // look for all matches of list2() to list1(i) + // this allows for multiple values of ID1 or ID2 + first_j = j; + for (j = first_j; (dID == 0) && (j < ID2->Nelements); j++) { + J = N2[j]; + dID = ID1->elements.Int[I] - ID2->elements.Int[J]; + if (dID == 0) { + index1->elements.Int[Nmatch] = I; + index2->elements.Int[Nmatch] = J; + + Nmatch ++; + if (Nmatch >= NMATCH) { + NMATCH += DMATCH; + REALLOCATE (index1->elements.Int, opihi_int, NMATCH); + REALLOCATE (index2->elements.Int, opihi_int, NMATCH); + } + } + } + j = first_j; + i++; + } + index1->Nelements = Nmatch; + index2->Nelements = Nmatch; + + free (N1); + free (N2); + + return (TRUE); +} + +// find the elements of ID1 which match ID2 (-1 if no match) +int join_IDs_outer (Vector *ID1, Vector *ID2, Vector *index) { + + off_t i, j, Jfirst, Ji, I, J, *N1, *N2, NMATCH; + opihi_int dID; + + NMATCH = ID1->Nelements; + ResetVector (index, OPIHI_INT, NMATCH); + + for (i = 0; i < index->Nelements; i++) { index->elements.Int[i] = -1; } + + ALLOCATE (N1, off_t, ID1->Nelements); + ALLOCATE (N2, off_t, ID2->Nelements); + + for (i = 0; i < ID1->Nelements; i++) { N1[i] = i; } + for (i = 0; i < ID2->Nelements; i++) { N2[i] = i; } + + sort_IDs_indexonly (ID1->elements.Int, N1, ID1->Nelements); + sort_IDs_indexonly (ID2->elements.Int, N2, ID2->Nelements); + + // find the closest entry in list 2 to the current entry in list 1: + for (i = j = 0; (i < ID1->Nelements) && (j < ID2->Nelements);) { + I = N1[i]; + J = N2[j]; + + dID = ID1->elements.Int[I] - ID2->elements.Int[J]; + + if (dID < 0) { + // no match in list 2 to this entry + index->elements.Int[I] = -1; + i++; + continue; + } + if (dID > 0) { j++; continue; } + + // look for closest matches of list2() to list1(i) + Jfirst = -1; + for (Ji = j; (Jfirst == -1) && (dID == 0) && (Ji < ID2->Nelements); Ji++) { + J = N2[Ji]; + dID = ID1->elements.Int[I] - ID2->elements.Int[J]; + if (dID == 0) { + Jfirst = J; + } + } + + // no match in list 2 to this entry + if (Jfirst == -1) { + index->elements.Int[I] = -1; + i++; + continue; + } + index->elements.Int[I] = Jfirst; + i++; + } + index->Nelements = NMATCH; + + free (N1); + free (N2); + + return (TRUE); +} Index: /trunk/Ohana/src/opihi/cmd.data/match2d.c =================================================================== --- /trunk/Ohana/src/opihi/cmd.data/match2d.c (revision 36083) +++ /trunk/Ohana/src/opihi/cmd.data/match2d.c (revision 36084) @@ -14,4 +14,7 @@ Vector *index1, *index2; + if ((N = get_argument (argc, argv, "-h"))) goto usage; + if ((N = get_argument (argc, argv, "--help"))) goto usage; + CLOSEST = FALSE; if ((N = get_argument (argc, argv, "-closest"))) { @@ -38,6 +41,5 @@ if (argc != 6) { gprint (GP_ERR, "USAGE: match2d X1 Y1 X2 Y2 Radius [-index1 (index1)] [-index2 (index2)] [-closest]\n"); - gprint (GP_ERR, " if -closest is provided, index1 & index2 will have the same length as X1 and X2 (respectively)\n"); - gprint (GP_ERR, " with either the index of the match or a value of -1 for non-matches\n"); + gprint (GP_ERR, " use -h or --help for more detail\n"); return (FALSE); } @@ -90,4 +92,21 @@ return (TRUE); + +usage: + gprint (GP_ERR, "we have two modes of operation:\n\n"); + + gprint (GP_ERR, "without -closest, we are finding all matched pairs within the match radius. in this\n"); + gprint (GP_ERR, "case, the two index vectors have the same length, one entry per matched pair.\n"); + gprint (GP_ERR, "x1[index1],y1[index1] matches to x2[index2],y2[index2].\n\n"); + + gprint (GP_ERR, "with -closest selected, we are finding the closest element of set 1 to each of set 2\n"); + gprint (GP_ERR, "and vice versa. in this case, index1 is always the same length as x1,y1, while index2\n"); + gprint (GP_ERR, "is the same lengths as x2,y2. x2[index1],y2[index1] matches x1,y1 while\n"); + gprint (GP_ERR, "x1[index2],y1[index2] matches x2,y2\n\n"); + + gprint (GP_ERR, "if -index1 or -index2 is not supplied, the vectors are created with names index1 or index2\n"); + gprint (GP_ERR, "use 'reindex' to generate new vectors based on these index vectors\n"); + + return FALSE; } Index: /trunk/Ohana/src/opihi/cmd.data/test/join.sh =================================================================== --- /trunk/Ohana/src/opihi/cmd.data/test/join.sh (revision 36084) +++ /trunk/Ohana/src/opihi/cmd.data/test/join.sh (revision 36084) @@ -0,0 +1,110 @@ + +list tests + test1 + memtest1 +end + +# Test if imhist works +macro test1 + + $PASS = 1 + + local i + + delete ID1 ID2 + + vlist -int ID1 1 2 3 4 5 6 + vlist val1 0.1 4.3 2.1 5.5 2.1 4.2 + + vlist -int ID2 5 3 1 6 + vlist val 0.3 2.3 1.1 2.5 + + join ID1 ID2 + +mana: for i 0 index1[] +>> echo index1[$i] index2[$i] +>> end +0 2 +2 1 +4 0 +5 3 +mana: join ID1 ID2 -outer +mana: vectors + N name size + 0 ID1 6 (INT) + 1 val1 6 (FLT) + 2 ID2 4 (INT) + 3 val 4 (FLT) + 4 index1 6 (INT) + 5 index2 4 (INT) +mana: for i 0 index1[] +>> echo index1[$i] +>> end +2 +-1 +1 +-1 +0 +3 +mana: echo ID1[2] +3 +mana: echo ID1[0] ID2[2] +1 1 +mana: for i 0 index1[] +>> if (index1[$i] == -1) continue +>> $N = index1[$i] +>> echo $i index1[$i] ID1[$i] ID2[$N] +>> end +0 2 1 1 +2 1 3 3 +4 0 5 5 +5 3 6 6 +mana: for i 0 index2[] +>> if (index2[$i] == -1) continue +>> $N = index2[$i] +>> echo $i index2[$i] ID2[$i] ID1[$N] +>> end +0 4 5 5 +1 2 3 3 +2 0 1 1 +3 5 6 6 +mana: + + # if ((xvec[1024] != 17) || (yvec[1024] != 100)) + # $PASS = 0 + # echo "Value mismatch: xvec[1024] yvec[1024] (should be 17,100)" + # end + # + # imhist -q buff xvec yvec -region 40 0 25 10 -range 0 10 + # + # if ((xvec[1024] != 10) || (yvec[1024] != 100)) + # $PASS = 0 + # echo "Value mismatch: xvec[1024] yvec[1024] (should be 10,100)" + # end + +end + + +# Memory test +macro memtest1 + + local i + + list word -x "ps -p $PID -o rss" + $startmem = $word:1 + + for i 0 1000 + imhist -q buff xvec yvec -region 40 0 25 10 -range 0 10 + end + + list word -x "ps -p $PID -o rss" + $endmem = $word:1 + + $PASS = 1 + + if ($endmem - $startmem > 10) + $PASS = 0 + echo "growth: {$endmem-$startmem}" + echo "kB/loop: {($endmem-$startmem)/1000}" + end +end Index: /trunk/Ohana/src/opihi/cmd.data/vlist.c =================================================================== --- /trunk/Ohana/src/opihi/cmd.data/vlist.c (revision 36084) +++ /trunk/Ohana/src/opihi/cmd.data/vlist.c (revision 36084) @@ -0,0 +1,63 @@ +# include "data.h" + +int vlist (int argc, char **argv) { + + int i, N; + Vector *vec; + + int APPEND = FALSE; + if ((N = get_argument (argc, argv, "-append"))) { + APPEND = TRUE; + remove_argument (N, &argc, argv); + } + + int INT = FALSE; + if ((N = get_argument (argc, argv, "-int"))) { + if (APPEND) { + gprint (GP_ERR, " ERROR : cannot mix -append and -int\n"); + return (FALSE); + } + INT = TRUE; + remove_argument (N, &argc, argv); + } + + if (argc < 3) { + gprint (GP_ERR, "USAGE: vlist [-append] [-int] vector value value ...\n"); + gprint (GP_ERR, " -int : resulting vector is integer type (values are truncated to integer)\n"); + gprint (GP_ERR, " -append : append value to end of vector (otherwise deleted)\n"); + return (FALSE); + } + + if ((vec = SelectVector (argv[1], ANYVECTOR, TRUE)) == NULL) return (FALSE); + + // for append, type comes from input vector + if (APPEND) { + INT = (vec[0].type == OPIHI_INT); + } + + int Nitem = argc - 2; + int Noffset = APPEND ? vec[0].Nelements : 0; + vec[0].Nelements = APPEND ? vec[0].Nelements + Nitem : Nitem; + + if (INT) { + vec[0].type = OPIHI_INT; + REALLOCATE (vec[0].elements.Int, opihi_int, vec[0].Nelements); + for (i = 0; i < Nitem; i++) { + vec[0].elements.Int[i + Noffset] = atoi (argv[i+2]); + } + } else { + vec[0].type = OPIHI_FLT; + REALLOCATE (vec[0].elements.Flt, opihi_flt, vec[0].Nelements); + for (i = 0; i < Nitem; i++) { + vec[0].elements.Flt[i + Noffset] = atof (argv[i+2]); + } + } + + return (TRUE); +} + +// vlist name 0 1 2 3 +// argc = 6 +// Nitem = 4 +// Noffset = 0 +// Int[0,1,2,3] = argv[2,3,4,5] Index: /trunk/Ohana/src/opihi/dvo/avmatch.c =================================================================== --- /trunk/Ohana/src/opihi/dvo/avmatch.c (revision 36083) +++ /trunk/Ohana/src/opihi/dvo/avmatch.c (revision 36084) @@ -70,40 +70,6 @@ dbExtractAveragesInit (); - // this does all the work of re-packaging the command, calling it on the remote machines, then loading in the results RAvec = NULL; DECvec = NULL; - if (PARALLEL && !HOST_ID) { - if (!CoordsFile) { - // get vectors corresponding to coordinates of interest - if ((RAvec = SelectVector (argv[1], ANYVECTOR, TRUE)) == NULL) goto help; - if ((DECvec = SelectVector (argv[2], ANYVECTOR, TRUE)) == NULL) goto help; - - ALLOCATE (vec, Vector *, 2); - vec[0] = RAvec; - vec[1] = DECvec; - - CoordsFile = abspath("coords.fits", 1024); - int status = WriteVectorTableFITS (CoordsFile, "COORDS", vec, 2, FALSE, NULL); - if (!status) goto escape; - } - - char *targv1 = argv[1]; - char *targv2 = argv[2]; - argv[1] = strcreate ("-coords"); - argv[2] = strcreate (CoordsFile); - free (CoordsFile); - - // I need to pass the RA & DEC vectors to the remote clients... - int status = HostTableParallelOps (skylist, argc, argv, RESULT_FILE, TRUE, RAvec->Nelements, VERBOSE); - if (vec) free (vec); - - free (argv[1]); - free (argv[2]); - argv[1] = targv1; - argv[2] = targv2; - - return status; - } - // get vectors corresponding to coordinates of interest if (CoordsFile) { @@ -120,4 +86,55 @@ remove_argument (1, &argc, argv); } + + /* load regions which contain all supplied RA,DEC coordinates */ + if ((skylist = SelectRegionsByCoordVectors (RAvec, DECvec)) == NULL) goto escape; + + // this does all the work of re-packaging the command, calling it on the remote machines, then loading in the results + if (PARALLEL && !HOST_ID) { + + // We need to copy the args to a temp array and modify them so that we send the + // correct set to the remote client. The args list looks like this: + // if (!CoordsFile) : avmatch (RADIUS) field, ... [we removed RA & DEC above] + // if ( CoordsFile) : avmatch (RADIUS) field, ... [because we stripped off the -coords filename elements] + + // allocate the temp array and copy all but (RA) (DEC) + int targc = 0; + char **targv = NULL; + ALLOCATE (targv, char *, argc + 2); + for (i = 0; i < argc; i++) { + targv[targc] = strcreate (argv[i]); + targc ++; + } + + // if not specified, create the coords.fits input file + // NOTE: RAvec, DECvec were set above + if (!CoordsFile) { + ALLOCATE (vec, Vector *, 2); + vec[0] = RAvec; + vec[1] = DECvec; + + CoordsFile = abspath("coords.fits", 1024); + int status = WriteVectorTableFITS (CoordsFile, "COORDS", vec, 2, FALSE, NULL); + if (!status) goto escape; + } + + // add the coords file to the args list + targv[targc+0] = strcreate ("-coords"); + targv[targc+1] = CoordsFile; // this gets freed with targv + targc += 2; + + // I need to pass the RA & DEC vectors to the remote clients... + int status = HostTableParallelOps (skylist, targc, targv, RESULT_FILE, TRUE, RAvec->Nelements, VERBOSE); + if (vec) free (vec); + + // free up targv + for (i = 0; i < targc; i++) { + free (targv[i]); + } + free (targv); + + return status; + } + RADIUS = atof (argv[1]); remove_argument (1, &argc, argv); @@ -140,7 +157,4 @@ if (fields[i].photcode[0].type == PHOT_DEP) needMeasures = TRUE; } - - /* load regions which contain all supplied RA,DEC coordinates */ - if ((skylist = SelectRegionsByCoordVectors (RAvec, DECvec)) == NULL) goto escape; /* create output storage vectors */ Index: /trunk/Ohana/src/opihi/dvo/hosts.c =================================================================== --- /trunk/Ohana/src/opihi/dvo/hosts.c (revision 36083) +++ /trunk/Ohana/src/opihi/dvo/hosts.c (revision 36084) @@ -12,5 +12,5 @@ gprint (GP_ERR, " commands:\n"); gprint (GP_ERR, " purge-temp : delete all tempfiles for this shell\n"); - gprint (GP_ERR, " : [-old-pid] [-all-pid] [-v] [-verbose] [-commit] [-type type]\n"); + gprint (GP_ERR, " : [-old-pid] [-all-pid] [-v] [-verbose] [-commit] [-type type] [-age hours]\n"); gprint (GP_ERR, " get-results : determine name of RESULTS file\n"); return FALSE; @@ -33,4 +33,15 @@ } + struct timeval now; + gettimeofday (&now, NULL); + + float AGE = 0; + int NOW = now.tv_sec; + if ((N = get_argument (argc, argv, "-age"))) { + remove_argument (N, &argc, argv); + AGE = atof(argv[N]); + remove_argument (N, &argc, argv); + } + int VERBOSE = FALSE; if ((N = get_argument (argc, argv, "-v"))) { @@ -121,5 +132,14 @@ glob (name, 0, NULL, &pglob); int j; + struct stat filestats; for (j = 0; j < pglob.gl_pathc; j++) { + if (AGE > 0) { + if (stat(pglob.gl_pathv[j], &filestats)) { + gprint (GP_ERR, "failed to get stats for %s\n", pglob.gl_pathv[j]); + continue; + } + float myAge = (NOW - filestats.st_mtime) / 3600.0; + if (myAge < AGE) continue; + } if (VERBOSE) gprint (GP_ERR, "unlink %s\n", pglob.gl_pathv[j]); if (!DRYRUN) unlink (pglob.gl_pathv[j]); Index: /trunk/Ohana/src/opihi/lib.shell/ListOps.c =================================================================== --- /trunk/Ohana/src/opihi/lib.shell/ListOps.c (revision 36083) +++ /trunk/Ohana/src/opihi/lib.shell/ListOps.c (revision 36084) @@ -189,12 +189,19 @@ int is_list_data (char *line) { - char *comm, *temp; - - temp = thisword (nextword (nextword (line))); - comm = thisword (line); - + char *comm = NULL; + char *temp = NULL; + char *ptr = NULL; + + comm = thisword (line); if (comm == NULL) goto escape; - if (strcmp (comm, "list")) goto escape; + + ptr = nextword (line); + if (!strncmp("-excel", ptr, strlen("-excel"))) ptr = nextword (ptr); + if (!strncmp("-excel-style", ptr, strlen("-excel-style"))) ptr = nextword (ptr); + ptr = nextword (ptr); + if (!strncmp("-excel", ptr, strlen("-excel"))) ptr = nextword (ptr); + if (!strncmp("-excel-style", ptr, strlen("-excel-style"))) ptr = nextword (ptr); + temp = thisword (ptr); /* if (cond) (command) does not define a complete block */ Index: /trunk/Ohana/src/tools/src/medianfilter.c =================================================================== --- /trunk/Ohana/src/tools/src/medianfilter.c (revision 36083) +++ /trunk/Ohana/src/tools/src/medianfilter.c (revision 36084) @@ -96,14 +96,26 @@ for (j = 0; j < Npixrd; j++, O++) { Nvalid = 0; + int debug = FALSE && (n == 5) && (j >= 1776779) && (j <= 1776789); for (k = 0; k < Nfiles; k++) { v = (float *)tmpmatr[k].buffer; + if (debug) { + fprintf (stderr, "%8.2f ", v[j]); + } if (v[j] < MinValid) continue; - if (isnan(v[j])) continue; - if (isinf(v[j])) continue; - list[k] = v[j]; + if (isnan(v[j])) { + continue; + } + if (isinf(v[j])) { + continue; + } + list[Nvalid] = v[j]; Nvalid ++; } + if (debug) fprintf (stderr, "\n"); if (Nvalid == 0) { *O = 0.0; + if (debug) { + fprintf (stderr, "Nvalid is 0\n"); + } continue; } @@ -111,9 +123,17 @@ sum = 0; Nval = 0; + if (debug) fprintf (stderr, "list : "); for (k = fmin*Nvalid; k < fmax*Nvalid; k++) { + if (debug) { + fprintf (stderr, "%8.2f ", list[k]); + } sum += list[k]; Nval ++; } *O = (sum / Nval); + if (debug) { + fprintf (stderr, " = %f / %f = %f\n", sum, Nval, *O); + fprintf (stderr, "\n"); + } }