Merge branch 'RoyBranch' into teake/merge-roy (b1205e75) · Commits · csst-pipeline / msc / sextractor

src/analyse.c

+4 −5

Original line number	Diff line number	Diff line
		@@ -7,7 +7,7 @@
		*
		* This file part of: SExtractor
		*
		* Copyright: (C) 1993-2016 IAP/CNRS/UPMC
		* Copyright: (C) 1993-2016 Emmanuel Bertin -- IAP/CNRS/UPMC
		*
		* License: GNU General Public License
		*
		@@ -280,8 +280,7 @@ void examineiso(picstruct field, picstruct dfield, objstruct *obj,
		emy2 /= flux2; /* variance of ym */
		emxy /= flux2; /* covariance */

		/*-- Handle fully correlated profile
		s (which cause a singularity...) */
		/-- Handle fully correlated profiles (which cause a singularity...) /
		esum *= 0.08333/flux2;
		if (obj->singuflag && (emx2emy2-emxyemxy) < esum*esum)
		{

src/back.c

+44 −15

Original line number	Diff line number	Diff line
		@@ -52,8 +52,10 @@ void makeback(picstruct field, picstruct wfield, int wscale_flag)

		{
		backstruct backmesh,wbackmesh, bm,wbm;
		tabstruct tab, wtab;
		PIXTYPE buf,wbuf, buft,wbuft;
		OFF_T fcurpos,wfcurpos, wfcurpos2,fcurpos2, bufshift, jumpsize;
		OFF_T2 fcurpos,wfcurpos, wfcurpos2,fcurpos2, bufshift, jumpsize;
		OFF_T2 currentElement, wcurrentElement, currentElement2, wcurrentElement2;
		size_t bufsize, bufsize2,
		size,meshsize;
		int i,j,k,m,n, step, nlines,
		@@ -65,7 +67,11 @@ void makeback(picstruct field, picstruct wfield, int wscale_flag)
		/* If the weight-map is not an external one, no stats are needed for it */
		if (wfield && wfield->flags&(INTERP_FIELD\|BACKRMS_FIELD))
		wfield= NULL;

		tab = field->tab;
		if (wfield)
		wtab = wfield->tab;
		else
		wtab = NULL; /* to avoid gcc -Wall warnings */
		w = field->width;
		bw = field->backw;
		bh = field->backh;
		@@ -83,12 +89,17 @@ void makeback(picstruct field, picstruct wfield, int wscale_flag)

		wfcurpos = wfcurpos2 = 0; /* to avoid gcc -Wall warnings */
		QFTELL(field->file, fcurpos, field->filename);
		currentElement = (tab->currentElement == 0) ? 1 : tab->currentElement; // CFITSIO

		if (wfield)
		{
		QFTELL(wfield->file, wfcurpos, wfield->filename);
		wcurrentElement = (wtab->currentElement == 0) ? 1 : wtab->currentElement; // CFITSIO
		}

		/* Allocate a correct amount of memory to store pixels */

		bufsize = (OFF_T)w*bh;
		bufsize = (OFF_T2)w*bh;
		meshsize = (size_t)bufsize;
		nlines = 0;
		if (bufsize > (size_t)BACK_BUFSIZE)
		@@ -96,8 +107,8 @@ void makeback(picstruct field, picstruct wfield, int wscale_flag)
		nlines = BACK_BUFSIZE/w;
		step = (field->backh-1)/nlines+1;
		bufsize = (size_t)(nlines = field->backh/step)*w;
		bufshift = (step/2)*(OFF_T)w;
		jumpsize = (step-1)*(OFF_T)w;
		bufshift = (step/2)*(OFF_T2)w;
		jumpsize = (step-1)*(OFF_T2)w;
		}
		else
		bufshift = jumpsize = 0; /* to avoid gcc -Wall warnings */
		@@ -172,16 +183,19 @@ void makeback(picstruct field, picstruct wfield, int wscale_flag)
		{
		/---- Image size too big, we have to skip a few data !/
		QFTELL(field->file, fcurpos2, field->filename);
		if (wfield)
		currentElement2 = (tab->currentElement == 0) ? 1 : tab->currentElement; // CFITSIO
		if (wfield){
		QFTELL(wfield->file, wfcurpos2, wfield->filename);
		wcurrentElement2 = (wtab->currentElement == 0) ? 1 : wtab->currentElement; // CFITSIO
		}
		if (j == ny-1 && (n=field->height%field->backh))
		{
		meshsize = n*(size_t)w;
		nlines = BACK_BUFSIZE/w;
		step = (n-1)/nlines+1;
		bufsize = (nlines = n/step)*(size_t)w;
		bufshift = (step/2)*(OFF_T)w;
		jumpsize = (step-1)*(OFF_T)w;
		bufshift = (step/2)*(OFF_T2)w;
		jumpsize = (step-1)*(OFF_T2)w;
		free(buf);
		QMALLOC(buf, PIXTYPE, bufsize); /* pixel buffer */
		if (wfield)
		@@ -192,35 +206,46 @@ void makeback(picstruct field, picstruct wfield, int wscale_flag)
		}

		/---- Read and skip, read and skip, etc... /
		QFSEEK(field->file, bufshift*(OFF_T)field->bytepix, SEEK_CUR,
		QFSEEK(field->file, bufshift*(OFF_T2)field->bytepix, SEEK_CUR,
		field->filename);
		tab->currentElement += bufshift; // CFITSIO

		buft = buf;
		for (i=nlines; i--; buft += w)
		{
		read_body(field->tab, buft, w);
		if (i)
		QFSEEK(field->file, jumpsize*(OFF_T)field->bytepix, SEEK_CUR,
		if (i) {
		QFSEEK(field->file, jumpsize*(OFF_T2)field->bytepix, SEEK_CUR,
		field->filename);
		tab->currentElement += jumpsize; // CFITSIO
		}
		}

		if (wfield)
		{
		/------ Read and skip, read and skip, etc... now on the weight-map /
		QFSEEK(wfield->file, bufshift*(OFF_T)wfield->bytepix, SEEK_CUR,
		QFSEEK(wfield->file, bufshift*(OFF_T2)wfield->bytepix, SEEK_CUR,
		wfield->filename);
		wtab->currentElement += bufshift; // CFITSIO

		wbuft = wbuf;
		for (i=nlines; i--; wbuft += w)
		{
		read_body(wfield->tab, wbuft, w);
		weight_to_var(wfield, wbuft, w);
		if (i)
		QFSEEK(wfield->file, jumpsize*(OFF_T)wfield->bytepix, SEEK_CUR,
		if (i){
		QFSEEK(wfield->file, jumpsize*(OFF_T2)wfield->bytepix, SEEK_CUR,
		wfield->filename);
		wtab->currentElement += jumpsize; // CFITSIO

		}
		}
		}
		backstat(backmesh, wbackmesh, buf, wbuf, bufsize, nx, w, bw,
		wfield?wfield->weight_thresh:0.0);
		QFSEEK(field->file, fcurpos2, SEEK_SET, field->filename);
		tab->currentElement = currentElement2; // CFITSIO

		bm = backmesh;
		for (m=nx; m--; bm++)
		if (bm->mean <= -BIG)
		@@ -230,6 +255,7 @@ void makeback(picstruct field, picstruct wfield, int wscale_flag)
		if (wfield)
		{
		QFSEEK(wfield->file, wfcurpos2, SEEK_SET, wfield->filename);
		wtab->currentElement = wcurrentElement2; // CFITSIO
		wbm = wbackmesh;
		for (m=nx; m--; wbm++)
		if (wbm->mean <= -BIG)
		@@ -284,8 +310,11 @@ void makeback(picstruct field, picstruct wfield, int wscale_flag)

		/* Go back to the original position */
		QFSEEK(field->file, fcurpos, SEEK_SET, field->filename);
		if (wfield)
		tab->currentElement = currentElement; // CFITSIO
		if (wfield) {
		QFSEEK(wfield->file, wfcurpos, SEEK_SET, wfield->filename);
		wfield->tab->currentElement = wcurrentElement; // CFITSIO
		}

		/* Median-filter and check suitability of the background map */
		NFPRINTF(OUTPUT, "Filtering background map(s)");

src/catout.c

+1 −2

Original line number	Diff line number	Diff line
		@@ -1053,7 +1053,7 @@ void reendcat()
		{
		keystruct *key;
		tabstruct *tab;
		OFF_T pos;
		OFF_T2 pos;
		char *head;

		switch(prefs.cat_type)
		@@ -1121,4 +1121,3 @@ void zerocat(void)
		}

src/field.c

+2 −0

Original line number	Diff line number	Diff line
		@@ -70,6 +70,7 @@ picstruct newfield(char filename, int flags, int ext)
		field->cat = cat;
		nok = 0;
		tab = cat->tab;
		if(tab->isTileCompressed) nok++;
		if (tab->naxis >= 2
		&& strncmp(tab->xtension, "BINTABLE", 8)
		&& strncmp(tab->xtension, "ASCTABLE", 8))
		@@ -78,6 +79,7 @@ picstruct newfield(char filename, int flags, int ext)
		for (ntab=cat->ntab; ext2-- && ntab--;)
		{
		tab=tab->nexttab;
		if(tab->isTileCompressed) nok++;
		if (tab->naxis >= 2
		&& strncmp(tab->xtension, "BINTABLE", 8)
		&& strncmp(tab->xtension, "ASCTABLE", 8))

src/fits/fitsbody.c

+108 −4

Original line number	Diff line number	Diff line
		@@ -85,7 +85,8 @@ PIXTYPE alloc_body(tabstruct tab, void (func)(PIXTYPE ptr, int npix))
		tab->extname);

		/* Decide if the data will go in physical memory or on swap-space */
		npix = tab->tabsize/tab->bytepix;
		//npix = tab->tabsize/tab->bytepix;
		npix = tab->naxisn[0] * tab->naxisn[1];
		size = npix*sizeof(PIXTYPE);
		if (size < body_ramleft)
		{
		@@ -93,6 +94,7 @@ PIXTYPE alloc_body(tabstruct tab, void (func)(PIXTYPE ptr, int npix))
		if ((tab->bodybuf = malloc(size)))
		{
		QFSEEK(tab->cat->file, tab->bodypos, SEEK_SET, tab->cat->filename);
		tab->currentElement = 1; // CFITSIO
		read_body(tab, (PIXTYPE *)tab->bodybuf, npix);
		/---- Apply pixel processing /
		if (func)
		@@ -120,6 +122,7 @@ PIXTYPE alloc_body(tabstruct tab, void (func)(PIXTYPE ptr, int npix))
		if (!spoonful)
		spoonful = DATA_BUFSIZE;
		QFSEEK(tab->cat->file, tab->bodypos, SEEK_SET, tab->cat->filename);
		tab->currentElement = 1; // CFITSIO
		read_body(tab, buffer, spoonful/sizeof(PIXTYPE));
		/-- Apply pixel processing /
		if (func)
		@@ -295,6 +298,61 @@ void free_body(tabstruct *tab)
		}


		void readTileCompressed(tabstruct tab, size_t spoonful, double bufdata0) {

		int status, hdutype;

		// first of all, move to correct HDU
		status = 0; fits_movabs_hdu(tab->infptr, tab->hdunum, &hdutype, &status);
		if (status != 0) {

		printf("Error moving to HDU %d\n", tab->hdunum);
		fits_report_error(stderr, status);
		}

		// pixels count from 1
		if (tab->currentElement == 0) tab->currentElement = 1;

		// now read section of image
		int datatype;
		switch(tab->bitpix){
		case BYTE_IMG:
		datatype = TBYTE;
		break;
		case SHORT_IMG:
		datatype = TSHORT;
		break;
		case LONG_IMG:
		datatype = TLONG;
		break;
		case FLOAT_IMG:
		datatype = TFLOAT;
		break;
		case DOUBLE_IMG:
		datatype = TDOUBLE;
		break;
		}

		int anynul;
		double bscale = 1.0, bzero = 0.0, nulval = 0.;

		// turn off any scaling so that we copy raw pixel values
		status = 0; fits_set_bscale(tab->infptr, bscale, bzero, &status);

		// now read the image
		status = 0; fits_read_img(tab->infptr, datatype, tab->currentElement, spoonful, &nulval, bufdata0, &anynul, &status);

		// report reading error
		if (status != 0) {

		printf("CFITSIO ERROR reading start=%d end=%d absolute end=%d\n", tab->currentElement, (tab->currentElement + spoonful) , (tab->naxisn[0]*tab->naxisn[1]));
		fits_report_error(stderr, status);
		}

		// update file 'pointer'
		tab->currentElement += spoonful;
		}

		/***** read_body **********************************************************
		PROTO read_body(tabstruct tab, PIXTYPE ptr, long size)
		PURPOSE Read floating point values from the body of a FITS table.
		@@ -323,14 +381,20 @@ void read_body(tabstruct tab, PIXTYPE ptr, size_t size)
		int curval, dval, blankflag, ival, iblank;

		size_t i, bowl, spoonful, npix;
		PIXTYPE bs,bz;
		//PIXTYPE bs,bz;
		double bs,bz;

		/* a NULL cat structure indicates that no data can be read */
		if (!(cat = tab->cat))
		return;

		bs = (PIXTYPE)tab->bscale;
		bz = (PIXTYPE)tab->bzero;
		// this cast from double to float loses precision
		//bs = (PIXTYPE)tab->bscale;
		//bz = (PIXTYPE)tab->bzero;

		bs = tab->bscale;
		bz = tab->bzero;

		blankflag = tab->blankflag;

		switch(tab->compress_type)
		@@ -344,7 +408,13 @@ void read_body(tabstruct tab, PIXTYPE ptr, size_t size)
		if (spoonful>size)
		spoonful = size;
		bufdata = (char *)bufdata0;

		// CFITSIO
		if (tab->isTileCompressed)
		readTileCompressed(tab, spoonful, bufdata0);
		else
		QFREAD(bufdata, spoonful*tab->bytepix, cat->file, cat->filename);

		switch(tab->bitpix)
		{
		case BP_BYTE:
		@@ -381,6 +451,7 @@ void read_body(tabstruct tab, PIXTYPE ptr, size_t size)
		break;

		case BP_SHORT:
		if (!tab->isTileCompressed)
		if (bswapflag)
		swapbytes(bufdata, 2, spoonful);
		if (blankflag)
		@@ -416,6 +487,7 @@ void read_body(tabstruct tab, PIXTYPE ptr, size_t size)
		break;

		case BP_LONG:
		if (!tab->isTileCompressed)
		if (bswapflag)
		swapbytes(bufdata, 4, spoonful);
		if (blankflag)
		@@ -452,6 +524,7 @@ void read_body(tabstruct tab, PIXTYPE ptr, size_t size)

		#ifdef HAVE_LONG_LONG_INT
		case BP_LONGLONG:
		if (!tab->isTileCompressed)
		if (bswapflag)
		swapbytes(bufdata, 8, spoonful);
		if (blankflag)
		@@ -487,6 +560,7 @@ void read_body(tabstruct tab, PIXTYPE ptr, size_t size)
		break;
		#endif
		case BP_FLOAT:
		if (!tab->isTileCompressed)
		if (bswapflag)
		swapbytes(bufdata, 4, spoonful);
		#pragma ivdep
		@@ -497,6 +571,7 @@ void read_body(tabstruct tab, PIXTYPE ptr, size_t size)
		case BP_DOUBLE:
		if (bswapflag)
		{
		if (!tab->isTileCompressed)
		swapbytes(bufdata, 8, spoonful);
		#pragma ivdep
		for (i=spoonful; i--; bufdata += sizeof(double))
		@@ -655,7 +730,13 @@ void read_ibody(tabstruct tab, FLAGTYPE ptr, size_t size)
		if (spoonful>size)
		spoonful = size;
		bufdata = (char *)bufdata0;

		// CFITSIO
		if (tab->isTileCompressed)
		readTileCompressed(tab, spoonful, bufdata0);
		else
		QFREAD(bufdata, spoonful*tab->bytepix, cat->file, cat->filename);

		switch(tab->bitpix)
		{
		case BP_BYTE:
		@@ -665,6 +746,7 @@ void read_ibody(tabstruct tab, FLAGTYPE ptr, size_t size)
		break;

		case BP_SHORT:
		if (!tab->isTileCompressed)
		if (bswapflag)
		swapbytes(bufdata, 2, spoonful);
		#pragma ivdep
		@@ -673,6 +755,7 @@ void read_ibody(tabstruct tab, FLAGTYPE ptr, size_t size)
		break;

		case BP_LONG:
		if (!tab->isTileCompressed)
		if (bswapflag)
		swapbytes(bufdata, 4, spoonful);
		#pragma ivdep
		@@ -682,6 +765,7 @@ void read_ibody(tabstruct tab, FLAGTYPE ptr, size_t size)

		#ifdef HAVE_LONG_LONG_INT
		case BP_LONGLONG:
		if (!tab->isTileCompressed)
		if (bswapflag)
		swapbytes(bufdata, 8, spoonful);
		#pragma ivdep
		@@ -922,6 +1006,10 @@ void write_body(tabstruct tab, PIXTYPE ptr, size_t size)
		#pragma ivdep
		for (i=spoonful; i--;)
		(bufdata++) = ((ptr++)-bz)/bs;


		// TODO not yet writing CFitsIO from SExtractor. CFITSIO - only perform byte-swap if we are NOT writing a tile-compressed format using cfitsio
		// if (tab->infptr == NULL) // TODO
		if (bswapflag)
		swapbytes(cbufdata0, 4, spoonful);
		}
		@@ -943,6 +1031,22 @@ void write_body(tabstruct tab, PIXTYPE ptr, size_t size)
		"read_body()");
		break;
		}

		// CFITSIO - if cfitsio output file has been set up, then proceed to write using cfitsio
		if (0 && tab->infptr != NULL) { // TODO

		int status = 0; fits_write_img(tab->infptr, TFLOAT, tab->currentElement, spoonful, cbufdata0, &status);

		if (status != 0) {

		printf("CFITSIO ERROR writing start=%d end=%d absolute end=%d\n", tab->currentElement, (tab->currentElement + spoonful) , (tab->naxisn[0]*tab->naxisn[1]));
		fits_report_error(stderr, status);
		}

		tab->currentElement = tab->currentElement + spoonful;
		}
		// otherwise, continue with usual AstrOmatic fits writing routine
		else
		QFWRITE(cbufdata0, spoonful*tab->bytepix, cat->file, cat->filename);
		}
		break;