Removed PSF correction from AMICACAL until some bugs related to that are... (c2abace6) · Commits · aflab / astrogeology / ISIS3

isis/src/hayabusa/apps/amicacal/AmicaCalUtils.h

+58 −89

Original line number	Diff line number	Diff line
		@@ -31,9 +31,10 @@
		*
		*
		*/

		using namespace std;
		using namespace cv;
		using namespace std;





		@@ -42,12 +43,12 @@ namespace Isis {


		/**
		* @brief Load required NAIF kernels required for timing needs
		*
		* @description Load required NAIF kernels required for timing needs.
		* This method maintains the loading of kernels for HAYABUSA timing and
		* planetary body ephemerides to support time and relative positions of planet
		* bodies.
		*/

		static void loadNaifTiming() {
		static bool naifLoaded = false;
		if (!naifLoaded) {
		@@ -72,8 +73,6 @@ static void loadNaifTiming() {
		QString pckName3(pck3.expanded());
		QString pckName4(pck4.expanded());



		furnsh_c(leapsecondsName.toAscii().data());
		furnsh_c(sclkName.toAscii().data());

		@@ -91,13 +90,12 @@ static void loadNaifTiming() {


		/**
		* @brief Computes the distance from the Sun to the observed body
		*
		* @description Computes the distance from the Sun to the observed body.
		* This method requires the appropriate NAIK kernels to be loaded that
		* provides instrument time support, leap seconds and planet body ephemeris.
		*
		* @return double Distance in AU between Sun and observed body
		* @return @b double Distance in AU between Sun and observed body.
		*/

		static bool sunDistanceAU(const QString &scStartTime,
		const QString &target,
		double &sunDist) {
		@@ -110,15 +108,14 @@ static bool sunDistanceAU(const QString &scStartTime,
		SpiceInt tcode;
		SpiceBoolean found;
		bodn2c_c(target.toAscii().data(), &tcode, &found);
		//cout << "tcode = " << tcode << endl;
		//cout << "found = " << found << endl;

		if (!found) return (false);

		// Convert starttime to et
		double obsStartTime;
		scs2e_c(-130, scStartTime.toAscii().data(), &obsStartTime);

		//cout << "obsStartTime = " << obsStartTime << endl;


		// Get the vector from target to sun and determine its length
		double sunv[3];
		@@ -129,7 +126,7 @@ static bool sunDistanceAU(const QString &scStartTime,

		double sunkm = vnorm_c(sunv);

		//cout << "lt = " << lt << endl;

		// Return in AU units
		sunDist = sunkm / 1.49597870691E8;

		@@ -145,7 +142,7 @@ static bool sunDistanceAU(const QString &scStartTime,
		*
		* @param icube A pointer to the input cube
		*
		* @return Mat A pointer to the OpenMat object
		* @return @b Mat A pointer to the OpenMat object
		*/

		Mat * isis2mat(Cube *icube) {
		@@ -170,8 +167,9 @@ return matrix;

		}


		/**
		* @brief Translates an OpenMat object to an Isis::Cube with one band
		* @brief Translates an OpenMat object to an ISIS::Cube with one band
		*
		* @author 2016-04-19 Tyler Wilson
		*
		@@ -206,8 +204,6 @@ void mat2isis(Mat *matrix, QString cubeName) {

		}



		}


		@@ -218,11 +214,10 @@ void mat2isis(Mat *matrix, QString cubeName) {
		*
		* @param matrix A pointer to the OpenMat object
		*
		* @param cubeName The name of the Isis::Cube that is being created.
		* @param cubeName The name of the ISIS::Cube that is being created.
		*
		*/


		void translate(Cube flatField, int transform, QString fname) {

		Mat * originalMat = isis2mat(flatField);
		@@ -233,12 +228,9 @@ void translate(Cube flatField, int transform, QString fname) {
		int lastsample = transform[3];
		int lastline = transform[4];



		int width = (lastsample-startsample);
		int height = (lastline-startline);


		Size sz(flatField->lineCount()/scale,flatField->sampleCount()/scale);

		Mat * resizedMatrix = new Mat();
		@@ -258,9 +250,9 @@ void translate(Cube flatField, int transform, QString fname) {
		mat2isis(resizedMatrix,fname);
		}


		}


		/**
		* @brief This function characterizes the point spread function
		* for near light sources. See equation [1] in:
		@@ -276,23 +268,20 @@ void translate(Cube flatField, int transform, QString fname) {
		* @param r The distance (in pixels) from the optical center of the
		* point source
		*
		* @return
		* @return @b double An estimate of the diffusion of light from near light sources
		* at pixel with coordinates (x,y) relative to the central pixel (with coordinates (0,0) ).
		*
		*/



		static double f_focused(double alpha,int binning,double x,double y) {

		double X = x*binning;
		double Y = y*binning;
		double r = sqrt(XX+YY);
		return exp(-alpha*r);

		}



		/**
		* @brief This function characterizes an attenuation function modeled as a sum
		* of Gaussian functions. See equation [2] in:
		@@ -307,44 +296,29 @@ static double f_focused(double alpha,int binning,double x,double y) {
		*
		* @param sigma An empirically derived vector of standard deviations of length N
		* also stored in $hayabusa/calibration/amica/amicaCalibration????.trn"
		* @param N The number of filter-bands (N=6)
		* @param binning The number of lines/samples which are binned together.
		* @param x The x-coordinate in pixels of the current pixel we are evaluating.
		* @param y The y-coordinate in pixels of the current pixel we are evaluating.
		*
		* @return
		* @return @b double An estimate of diffuse light at the point (x,y) relative to the
		* central pixel (at coordinates (0,0) ).
		*
		*/

		/*
		* static double f_diffuse(double A,double sigma,int N,double r){

		double sum = 0;

		for (int i = 0; i < N; i ++) {
		sum += (A[i]/(sigma[i]sqrt(2.0pi_c() ) ) )exp(-(rr)/(2sigma[i]sigma[i]) );
		}


		return sum;


		}
		*/


		static double f_unfocused(double * A,double * sigma, int N,int binning,double x,double y) {

		//double X = abs(binning*x-512);
		//double Y = abs(binning*y-512);

		double X = binning*x;
		double Y = binning*y;

		double r = sqrt(XX+YY);


		double sum = 0;

		for (int i = 0; i < N; i ++) {
		sum += (A[i]/(sigma[i]sqrt(2.0pi_c() ) ) )exp(-(rr)/(2sigma[i]sigma[i]) );


		}


		@@ -353,51 +327,51 @@ static double f_unfocused(double * A,double * sigma, int N,int binning,double x,
		}


		/*

		static double psf(double A, double sigma, int N,double alpha,double Ixy,int binning,double x, double y) {



		//double psfVal = f_unfocused(A,sigma,N,x,y);
		//double psfVal = Ixy*f_unfocused(A,sigma,N,binning,x,y);

		//double psfVal = Ixy*(f_unfocused(A,sigma,N,binning,x,y) + f_focused(alpha,x,y) );

		double psfVal = Ixy*f_unfocused(A,sigma,N,binning,x,y);



		return psfVal;



		}


		/**
		* @brief This function returns a matrix of [size x size] of light distribution values.
		* The center pixel value is at the center of the matrix, and the values around the central
		* pixel represent the fraction of light intensity from the central pixel that seeps into the
		* neighboring pixels.
		* @param size The dimension of the matrix (in pixels).
		* @param A An empirically derived vector of constants of length N which is stored in
		* $hayabusa/calibration/amica/amicaCalibration????.trn"
		* Each filter has a different A vector.
		* @param sigma sigma An empirically derived vector of standard deviations of length N
		* also stored in $hayabusa/calibration/amica/amicaCalibration????.trn"
		* @param alpha alpha An empirically derived constant which is stored in
		* $hayabusa/calibration/amica/amicaCalibration????.trn" Each filter has a different
		* alpha value.
		* @param N
		* @param binning
		* @return @b double * A pointer to a [size x size] matrix of light distribution values.
		*/





		double * setPSFFilter(int size, double A,double sigma, double alpha,int N,int binning) {


		double * psfVals = new double[size*size];


		int i = 0;


		for(double y = -(size / 2) ; y <= (size / 2) ; y++) {
		for(double x = -(size / 2) ; x <= (size / 2) ; x++) {

		if (x == 0 && y ==0) {


		psfVals[i] = 0;
		i++;

		}
		else {

		psfVals[i]=f_unfocused(A,sigma,N,binning,x,y) +f_focused(alpha,binning,x,y);
		//psfVals[i]=f_unfocused(A,sigma,N,binning,x,y);
		//psfVals[i] = f_focused(alpha,binning,x,y);

		i++;

		}
		//i++;
		}
		}

		@@ -405,14 +379,9 @@ double * setPSFFilter(int size, double A,double sigma, double alpha,int N,int
		}




		}






		#endif

isis/src/hayabusa/apps/amicacal/amicacal.cpp

+213 −262

File changed.

Preview size limit exceeded, changes collapsed.

isis/src/hayabusa/apps/amicacal/amicacal.xml

+78 −25

Original line number	Diff line number	Diff line
		@@ -17,14 +17,15 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		<li><a href="#HotPixel">Step 3: Remove the hot pixels.</a ></li>
		<li><a href="#Smear">Step 4: Readout smear correction.</a ></li>
		<li><a href="#Flatfield">Step 5: Uniformity (flat field) correction.</a ></li>
		<li><a href="#ConvertIOF">Step 6: Conversion to I over F</a ></li>
		<li><a href="#ConvertIOF">Step 6: Conversion to I/F</a ></li>
		<li><a href="#PSFCorrection">Step 7: PSF correction</a ></li>
		<li><a href="#Notes">Notes</a></li>
		<li><a href="#References">References</a ></li>
		</ul>

		<p>
		<h1><a id="BiasRemoval">Step 1: Bias removal.</a ></h1>
		The signal level of AMICA is determined by a combination of an of an electronic offset
		<h1><a id="BiasRemoval:">Step 1: Bias removal.</a ></h1>
		The signal level of AMICA is determined by a combination of an electronic offset
		which defines the "zero" level of the analog-to-digital converter (the bias) as well
		as any additional signal due to thermal electrons (the dark current). In general, the bias
		offset depends upon the electronics and CCD temperature. The bias for AMICA was empirically
		@@ -32,8 +33,8 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica

		<div style="font-size:90%; font-family:monospace; margin:5px; margin-left:30px;">

		<center><strong>BIAS(t) = B<sub>0</sub> + B<sub>1</sub>t + B<sub>2</sub>t<sup>2</sup>
		</strong>, where <strong>t</strong> is in units of days.</center>
		<center><strong>BIAS(t) = B<sub>0</sub> + B<sub>1</sub>t + B<sub>2</sub>t<sup>2</sup>,
		</strong> where <strong>t</strong> is in units of days.</center>
		</div>

		The constants <strong>B<sub>0</sub></strong>, <strong>B<sub>1</sub></strong> and
		@@ -66,7 +67,7 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		<p>
		<h1><a id="Linearity">Step 2: Linearity correction:</a ></h1>
		The DN level stored in raw images is approximately proportional to the number
		of photons detected in each CCD pixel over the exposure time. Linearity is a meausure
		of photons detected in each CCD pixel over the exposure time. Linearity is a measure
		of the degree to which the CCD response is proportional to the incident flux.
		Linearity was tested during the pre-flight and mission phase of the mission. During the
		pre-flight phase, DN values taken at different exposure times showed a linear intensity with
		@@ -117,12 +118,12 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		</p>
		<p>
		<h1><a id="HotPixel">Step 3: Hot pixel removal:</a ></h1>
		During the mission the CCD was exposed to intense radition from cosmic rays. This resulted
		During the mission the CCD was exposed to intense radiation from cosmic rays. This resulted
		in some pixels with an aberrant dark current value that was higher than their neighbors
		(the hot pixels). The number of hot pixels increased during the mission. This made
		accurate measurement of faint objects such as stars difficult. The following hot pixels
		were identified and removed from the output images during calibration by setting
		their values to Isis::Null:
		their values to ISIS::Null:

		<br/><br/>

		@@ -171,21 +172,47 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		is produced. The read-out smear is predominant in images with exposure times
		on the order of 100 microseconds. Originally the mission plan was to correct
		for this smear on-board the spacecraft. However, an anomaly in one of the reaction
		wheels on October 2, 2005 made this impossible. So some of the images returned
		by the AMICA instrument are smear-corrected, while others are not. This necessitated
		wheels on October 2, 2005 made this impossible. As a result of this, some of the images
		returned by the AMICA instrument are smear-corrected, while others are not. This necessitated
		creating a smear model to remove smears from images taken during the descending and ascending
		period. The read-out smear brightness <strong>I<sub>smear</sub></strong>
		for unbinned images is modeled from the observed images as follows:

		<div style="font-size:90%; font-family:monospace; margin:5px; margin-left:30px;">
		<br/><br/>
		<center><strong>I<sub>smear</sub>(H) = SUM(H=0...N<sub>v</sub> -1) [K*[(I<sub>raw</sub>(H,V) -
		I<sub>sky</sub>(H,V)]/N<sub>v</sub>]</strong></center>
		<center><strong>I<sub>smear</sub>(H) = SUM(H=0...N<sub>v</sub> -1) { [K*[(I<sub>raw</sub>(H,V) -
		I<sub>sky</sub>(H,V)]/N<sub>v</sub>] } </strong></center>
		<br/><br/>
		<center>Where: <strong>K = t<sub>VCT</sub>/(t<sub>VCT</sub>+t<sub>exp</sub>)
		</strong></center>
		<br/><br/>



		For binned images a different formula is applied:


		<center><strong>I<sub>smear</sub>(H) = SUM(H=0...N<sub>v</sub> - 1)
		{ [ K1*I<sub>raw</sub>(H,V) ] }</strong>
		</center>

		<br/><br/>
		<center>Where: <strong>K1 = (1/N<sub>v</sub>)*[t<sub>VCT</sub>/
		(t<sub>VCT</sub>+t<sub>exp</sub>)]
		</strong></center>

		<br/><br/>

		A correction factor is then applied after the smear component is subtracted
		from the input image (I<sub>raw</sub>):

		<center> <strong> I<sub>out</sub> = C*(I<sub>raw></sub> - I<sub>smear</sub>)</strong></center>
		<br/><br/>

		<center>Where: <strong>C = 1/( 1 + K1 * ( (B - 1)/2*B ) ) </strong></center>

		<br/><br/>

		<center><table border = "1" cellpadding="5">
		<tr>
		<td>Variable</td><td>Description</td>
		@@ -216,18 +243,27 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		<tr>
		<td>L</td><td>The sample number of the pixel.</td>
		</tr>
		<tr>
		<td>B</td><td>The binning number.</td>
		</tr>
		</table></center>



		</div>






		</p>
		<p>
		<h1><a id="Flatfield">Step 5: Flat-field correction:</a ></h1>
		Performs a correction for pixel-to-pixel variation in CCD response and vignetting
		(reduction of image brightness near the periphery compared with the center).
		Flat-field images for all bands were acquired using an integrating sphere
		at NEC Toshiba Space Systems Ltd at room temperature (around 30 degrees Celsius).
		at NEC Space Technologies, Ltd. at room temperature (around 30 degrees Celsius).
		A flat-field image is one which has constant uniform brightness everywhere. AMICA
		was pointed into the integrating sphere to acquire images of a field that is known be
		be spatially uniform to an accuracy of approximately 2%. The correction is accomplished
		@@ -236,7 +272,7 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		</p>
		<p>
		<h1><a id="ConvertIOF">Step 6: Convert units to I over F:</a ></h1>
		This converts the raw DN values to I/F (radiance) units. The formula is:
		This converts the raw DN values to I/F (radiance) units. The formula is as follows:
		<br/><br/>
		<div style="font-size:90%; font-family:monospace; margin:5px; margin-left:30px;">
		<center><strong>IoF = Raw[pi(R<sub>s</sub>)<sup>2</sup>]/F</strong></center>
		@@ -264,6 +300,7 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica

		</p>
		<p>

		<h1><a id="PSFCorrection">Step 7: PSF Correction:</a ></h1>
		The point-spread function (PSF) is a model describing the two-dimensional distribution of
		light emitted from an infinitely small point source. For point light sources, there
		@@ -321,7 +358,7 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		Both the focused and the unfocused components are calculated for each pixel in the image, and
		subtracted out of the final calibrated image.
		<strong>A<sub>i</sub></strong> and <strong>sigma<sub>i</sub></strong> are constants
		derived empircally from images of Itokawa and the Moon and enumerated below:
		derived empirically from images of Itokawa and the Moon and enumerated below:
		<br/><br/>


		@@ -384,14 +421,30 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		</table></center>
		</div>
		</p>
		<p>
		<h1><a id="Notes">Notes</a></h1>
		<ol>
		<li>Any in-flight dark current noise for this mission was buried in the read-out noise (~ 4 DN).
		When observing Itokawa, the image intensity had a range of 1000-3000 DN. Therefore,
		the dark current was considered negligible, and is not accounted for in this calibration
		application (see Ishiguro (2010) ).
		</li>
		<li>The derivation for the formula for smear-removal in the case of unbinned images
		does not appear in any of the references given below. Questions concerning
		the derivation steps for this formula should be directed towards the USGS.
		</li>
		</ol>
		</p>

		<p>
		<h1><a id="References">References:</a></h1>
		<ol>
		<li>Ishiguro, M., et al. The Hayabusa Spacecraft Asteroid Multi-band Imaging Camera (AMICA).
		<li>Ishiguro, Masateru, et al. <i>"The Hayabusa Spacecraft Asteroid Multi-band Imaging
		Camera (AMICA)"</i>.
		Icarus 207(2010) 714-731.</li>
		<li>Ishiguro, M. Scattered light correction of Hayabusa/AMICA data and quantitative
		spectral comparisons of Itokawa. Publ. Astron. Soc. Japan (2014) 66(3), 55 (1-9).</li>
		<li>Ishiguro, Masateru. <i>"Scattered light correction of Hayabusa/AMICA data and
		quantitative spectral comparisons of Itokawa"</i>.
		Publ. Astron. Soc. Japan (2014) 66(3), 55 (1-9).</li>
		</ol>

		</p>
		@@ -426,9 +479,9 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		Added calibration correction routines for linearity, PSF, smear removal for binned
		and unbinned images, and converting to I/F units. Also added hot pixel removal.
		Updated the application documentation. Added routines for converting
		a single-band Isis cube to an OpenCV::OpenMat matrix and vice versa for the purpose
		a single-band ISIS cube to an OpenCV::OpenMat matrix and vice versa for the purpose
		of translating and scaling flat-field images so their dimensions match
		some AMICA images which have been cropped.
		some AMICA images which have been cropped. Fixes #3880.

		</change>
		</history>
		@@ -442,7 +495,7 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		Input cube
		</brief>
		<description>
		This is the Isis Cube to be calibrated.
		This is the ISIS cube to be calibrated.
		</description>
		<filter>
		*.cub
		@@ -463,7 +516,7 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		</filter>
		</parameter>


		<!--
		<parameter name="PSFCORRECTED">
		<type>cube</type>
		<fileMode>output</fileMode>
		@@ -478,7 +531,7 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		*.cub
		</filter>
		</parameter>

		-->



		@@ -488,7 +541,7 @@ xsi:noNamespaceSchemaLocation="http://isis.astrogeology.usgs.gov/Schemas/Applica
		<type>boolean</type>
		<default><item>True</item></default>
		<brief>
		Set polarized pixels to Isis::Null.
		Set polarized pixels to ISIS::Null.
		</brief>
		<description>
		Four position angle glass polarizers (of dimension 200 x 200 pixels each)

isis/src/hayabusa/apps/amicacal/tsts/iof/Makefile

+1 −3

Original line number	Diff line number	Diff line
		@@ -5,10 +5,8 @@ include $(ISISROOT)/make/isismake.tsts
		commands:
		$(APPNAME) \
		FROM=$(INPUT)/st_2459265790_w.cub \
		TO=$(OUTPUT)/corrected_iof_st_2459265790_w.cub \
		PSFCORRECTED=$(OUTPUT)/corrected_psfcorrected_iof_st_2459265790_w.cub >/dev/null;\
		TO=$(OUTPUT)/corrected_iof_st_2459265790_w.cub > /dev/null;\
		$(APPNAME) \
		FROM=$(INPUT)/st_2459265790_w.cub \
		TO=$(OUTPUT)/corrected_noiof_st_2459265790_w.cub \
		PSFCORRECTED=$(OUTPUT)/corrected_psfcorrected_noiof_st_2459265790_w.cub \
		iof=no > /dev/null;

isis/src/hayabusa/apps/amicacal/tsts/nullpolarpix/Makefile

+1 −3

Original line number	Diff line number	Diff line
		@@ -5,10 +5,8 @@ include $(ISISROOT)/make/isismake.tsts
		commands:
		$(APPNAME) \
		FROM=$(INPUT)/st_2459265790_w.cub \
		TO=$(OUTPUT)/corrected_polarpixnulled_st_2459265790_w.cub \
		PSFCORRECTED=$(OUTPUT)/corrected_polarpixnulled_psfcorrected_st_2459265790_w.cub >/dev/null;\
		TO=$(OUTPUT)/corrected_polarpixnulled_st_2459265790_w.cub >/dev/null;\
		$(APPNAME) \
		FROM=$(INPUT)/st_2459265790_w.cub \
		TO=$(OUTPUT)/corrected_polarpixnotnulled_st_2459265790_w.cub \
		PSFCORRECTED=$(OUTPUT)/corrected_polarpixnotnulled_psfcorrected_st_2459265790_w.cub \
		nullpolarpix=no > /dev/null;