Added hyb2onccal GTests (#4037) (8b0364e7) · Commits · aflab / astrogeology / ISIS3

isis/src/hayabusa2/apps/hyb2onc2isis/Hayabusa2OncInstrument.trn

+1 −1

Original line number	Diff line number	Diff line
		@@ -98,7 +98,7 @@ Group = SpacecraftClockEndCount
		Auto
		InputKey = P_SCCEC
		InputPosition = FitsLabels
		OutputName = SpacecraftClockStartCount
		OutputName = SpacecraftClockEndCount
		OutputPosition = (Object, IsisCube, Group, Instrument)
		Translation = (, )
		End_Group

isis/src/hayabusa2/apps/hyb2onccal/Hyb2OncCalUtils.h

+10 −3

Original line number	Diff line number	Diff line
		@@ -221,10 +221,17 @@ namespace Isis {
		}

		double smear = 0;
		for (int j = 0; j < imageIn.size(); j++) {
		int skipCount = 0;
		for (int i = 0; i < imageIn.size(); i++) {
		if (!IsSpecial(imageIn[i])) {
		// Left out g_darkCurrent subtraction for now.
		smear += ((imageIn[j] * pow(2.0, 12 - g_bitDepth) - g_bias) / imageIn.size());
		smear += (imageIn[i] * pow(2.0, 12 - g_bitDepth) - g_bias);
		}
		else {
		skipCount++;
		}
		}
		smear /= (imageIn.size() - skipCount);
		smear *= g_timeRatio;

		// Iterate over the line space

isis/src/hayabusa2/apps/hyb2onccal/hyb2onccal.cpp

+18 −9

Original line number	Diff line number	Diff line
		@@ -38,7 +38,17 @@ using namespace std;

		namespace Isis {

		void hyb2onccal(UserInterface &ui) {
		void hyb2onccal(UserInterface &ui, Pvl *log) {
		Cube icube;
		CubeAttributeInput inAtt = ui.GetInputAttribute("FROM");
		if (inAtt.bands().size() != 0) {
		icube.setVirtualBands(inAtt.bands());
		}
		icube.open(ui.GetFileName("FROM"));
		hyb2onccal(&icube, ui, log);
		}

		void hyb2onccal(Cube icube, UserInterface &ui, Pvl log) {
		g_calStep = ui.GetString("UNITS");

		const QString hyb2cal_program = "hyb2onccal";
		@@ -47,7 +57,7 @@ namespace Isis {
		QString hyb2cal_runtime = Application::DateTime();

		ProcessBySample p;
		Cube *icube = p.SetInputCube("FROM");
		p.SetInputCube(icube);

		// Basic assurances...
		if (icube->bandCount() != 1) {
		@@ -62,7 +72,7 @@ namespace Isis {
		g_filter = bandbin["FilterName"][0];
		}
		catch(IException &e) {
		QString msg = "Unable to read FilterName keyword in the BandBin group "
		QString msg = "Unable to read [FilterName] keyword in the BandBin group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io, msg, _FILEINFO_);
		}
		@@ -73,7 +83,7 @@ namespace Isis {
		instrument = inst["InstrumentId"][0];
		}
		catch(IException &e) {
		QString msg = "Unable to read InstrumentId keyword in the Instrument group "
		QString msg = "Unable to read [InstrumentId] keyword in the Instrument group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io, msg, _FILEINFO_);
		}
		@@ -188,9 +198,6 @@ namespace Isis {
		if (smearCorrection == "ONBOARD") {
		g_onBoardSmearCorrection = true;
		}
		else {
		qDebug() << icube->fileName();
		}

		QString compmode = inst["Compression"];
		// TODO: verify that the compression factor/scale is actually 16 for compressed Hayabusa2 images.
		@@ -268,7 +275,8 @@ namespace Isis {

		} //Finished setting flatfield file

		Cube *ocube = p.SetOutputCube("TO");
		Cube *ocube = p.SetOutputCube(ui.GetFileName("TO"), ui.GetOutputAttribute("TO"),
		icube->sampleCount(), icube->lineCount(), icube->bandCount());
		QString calfile = loadCalibrationVariables(ui.GetAsString("CONFIG"));
		g_timeRatio = g_Tvct/(g_texp + g_Tvct);

		@@ -281,6 +289,7 @@ namespace Isis {
		g_calibrationScale = 1.0 / (g_texp * g_sensitivity);
		units = "W / (m**2 micrometer sr)";
		}

		// Output units of I/F
		else if (g_calStep == "IOF") {
		// Convert to radiance
		@@ -369,7 +378,7 @@ namespace Isis {

		// Write Calibration group to output file
		ocube->putGroup(calibrationLog);
		Application::Log(calibrationLog);
		log->addGroup(calibrationLog);
		p.EndProcess();

		}

isis/src/hayabusa2/apps/hyb2onccal/hyb2onccal.h

+7 −8

Original line number	Diff line number	Diff line
		#ifndef HYB2ONCCAL_H
		#define HYB2ONCCAL_H
		#ifndef hyb2onccal_h
		#define hyb2onccal_h

		#include "Cube.h"
		#include "Pvl.h"
		#include "UserInterface.h"

		namespace Isis {

		extern void hyb2onccal(UserInterface &ui);

		extern void hyb2onccal(Cube icube, UserInterface &ui, Pvl log);
		extern void hyb2onccal(UserInterface &ui, Pvl *log);
		}

		#endif //hyb2onccal_h



		#endif // HYB2ONCCAL_H

isis/src/hayabusa2/apps/hyb2onccal/main.cpp

+9 −585

Original line number	Diff line number	Diff line
		@@ -135,595 +135,19 @@ static double g_v_standard(1.0);
		// static double g_v_standard(3.42E-3);//!< Base conversion for all filters (Tbl. 9)

		void IsisMain() {

		UserInterface &ui = Application::GetUserInterface();
		// g_iofCorrection = ui.GetString("UNITS");

		const QString hyb2cal_program = "hyb2onccal";
		const QString hyb2cal_version = "1.1";
		const QString hyb2cal_revision = "$Revision$";
		QString hyb2cal_runtime = Application::DateTime();

		ProcessBySample p;

		Cube *icube = p.SetInputCube("FROM");

		// Basic assurances...
		if (icube->bandCount() != 1) {
		throw IException(IException::User,
		"ONC images may only contain one band", _FILEINFO_);
		}

		PvlGroup &inst = icube->group("Instrument");
		PvlGroup &bandbin = icube->group("BandBin");

		try{
		g_filter = bandbin["FilterName"][0];
		}

		catch(IException &e) {
		QString msg = "Unable to read FilterName keyword in the BandBin group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io,msg, _FILEINFO_);

		}

		QString instrument("");

		try {
		instrument = inst["InstrumentId"][0];
		}
		catch(IException &e) {
		QString msg = "Unable to read InstrumentId keyword in the Instrument group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io,msg, _FILEINFO_);

		}

		if ( instrument=="ONC-W1" ) {
		g_instrument = InstrumentType::ONCW1;
		}
		else if ( instrument=="ONC-W2" ) {
		g_instrument = InstrumentType::ONCW2;
		}
		else if ( instrument == "ONC-T" ) {
		g_instrument = InstrumentType::ONCT;
		}
		else {
		QString msg = "Unidentified instrument key in the "
		"InstrumentId key of the Instrument Pvl group.";
		throw IException(IException::Io,msg, _FILEINFO_);
		}



		//Set up binning and image subarea mapping
		binning = inst["Binning"];
		int startLine = inst["SelectedImageAreaY1"];
		int startSample = inst["SelectedImageAreaX1"];
		int lastLine = inst["SelectedImageAreaY2"];
		int lastSample = inst["SelectedImageAreaX2"];

		AlphaCube myAlpha(1024,1024,icube->sampleCount(), icube->lineCount(),
		startSample,startLine,lastSample,lastLine);

		g_bitDepth = inst["BitDepth"];

		alpha = &myAlpha;

		try {
		g_texp = inst["ExposureDuration"][0].toDouble();
		}
		catch(IException &e) {
		QString msg = "Unable to read [ExposureDuration] keyword in the Instrument group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io,msg, _FILEINFO_);
		}


		try {
		g_AEtemperature = inst["ONCAETemperature"][0].toDouble();
		}
		catch(IException &e) {
		QString msg = "Unable to read [ONCAETemperature] keyword in the Instrument group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io,msg, _FILEINFO_);

		}

		try {
		g_CCD_T_temperature = inst["ONCTCCDTemperature"][0].toDouble();
		}
		catch(IException &e) {
		QString msg = "Unable to read [ONCTCCDTemperature] keyword in the Instrument group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io,msg, _FILEINFO_);
		}

		try {
		g_ECT_T_temperature = inst["ONCTElectricCircuitTemperature"][0].toDouble();
		}
		catch(IException &e) {
		QString msg = "Unable to read [ONCTElectricCircuitTemperature] keyword in the Instrument group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io,msg, _FILEINFO_);
		}


		try {

		g_startTime=inst["SpacecraftClockStartCount"][0];
		}

		catch (IException &e) {
		QString msg = "Unable to read [SpacecraftClockStartCount] keyword in the Instrument group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io,msg, _FILEINFO_);
		}


		QString smearCorrection;

		try {
		smearCorrection = inst["SmearCorrection"][0];
		}

		catch (IException &e) {
		QString msg = "Unable to read [SmearCorrection] keyword in the Instrument group "
		"from input file [" + icube->fileName() + "]";
		throw IException(e, IException::Io,msg, _FILEINFO_);
		}

		if (smearCorrection=="ONBOARD") {
		g_onBoardSmearCorrection=true;
		}

		QString compmode = inst["Compression"];
		// TODO: verify that the compression factor/scale is actually 16 for compressed Hayabusa2 images.
		g_compfactor = ( "lossy" == compmode.toLower() ) ? 16.0 : 1.0;

		QString target = inst["TargetName"];
		g_target = target;


		// NOTE we do not have a valid flat-field for the W1 or W2 images.
		FileName flatfile = "NONE";
		if (g_instrument == InstrumentType::ONCT) {
		QScopedPointer<Cube, TemporaryCubeDeleter> flatcube;
		flatfile = DetermineFlatFieldFile(g_filter);

		QString reducedFlat(flatfile.expanded());

		// Image is not cropped
		if (startLine == 1 && startSample == 1) {

		// Determine if we need to subsample the flat field if pixel binning occurred
		// TODO: test a binned image (add test case).
		if (binning > 1) {
		QString scale(toString(binning));
		FileName newflat = FileName::createTempFile("$TEMPORARY/" +
		flatfile.baseName() + "_reduced.cub");
		reducedFlat = newflat.expanded();
		QString parameters = "FROM=" + flatfile.expanded() +
		" TO=" + newflat.expanded() +
		" MODE=SCALE" +
		" LSCALE=" + scale +
		" SSCALE=" + scale;

		Pvl appLog;
		try {
		ProgramLauncher::RunIsisProgram("reduce", parameters);
		hyb2onccal(ui, &appLog);
		}
		catch (IException& ie) {
		remove(reducedFlat.toLatin1().data());
		throw ie;
		catch (...) {
		for (auto grpIt = appLog.beginGroup(); grpIt!= appLog.endGroup(); grpIt++) {
		Application::Log(*grpIt);
		}
		QScopedPointer<Cube, TemporaryCubeDeleter> reduced(new Cube(reducedFlat, "r"));
		flatcube.swap( reduced );
		throw;
		}

		// Set up processing for flat field as a second input file
		CubeAttributeInput att;
		p.SetInputCube(reducedFlat, att);
		for (auto grpIt = appLog.beginGroup(); grpIt!= appLog.endGroup(); grpIt++) {
		Application::Log(*grpIt);
		}
		else {
		// Image is cropped so we have to deal with it
		FileName transFlat =
		FileName::createTempFile("$TEMPORARY/" + flatfile.baseName() + "_translated.cub");

		Cube *flatOriginal = new Cube(flatfile.expanded() );

		int transform[5] = {binning,startSample,startLine,lastSample,lastLine};

		// Translates and scales the flatfield image. Scaling
		// might be necessary in the event that the raw image was also binned.

		translate(flatOriginal,transform,transFlat.expanded());

		QScopedPointer<Cube, TemporaryCubeDeleter> translated(new Cube(transFlat.expanded(), "r"));
		flatcube.swap(translated);

		CubeAttributeInput att;
		std::cout<<" before second setinputcube()...\n"<<std::endl;
		p.SetInputCube(transFlat.expanded(),att);
		std::cout<<" after second setinputcube()...\n"<<std::endl;
		}
		}

		Cube *ocube = p.SetOutputCube("TO");
		//QString fname = ocube->fileName();

		QString calfile = loadCalibrationVariables(ui.GetAsString("CONFIG"));

		g_timeRatio = g_Tvct/(g_texp + g_Tvct);

		g_iof = 1.0; // Units of DN

		QString g_units = "DN";
		// if ( "radiance" == g_iofCorrection.toLower() ) {
		// // Units of RADIANCE
		// g_iof = g_iof * g_iofScale;
		// g_units = "W / (m**2 micrometer sr)";
		// }
		//
		// if ( !sunDistanceAU(startTime, target, g_solarDist) ) {
		// throw IException(IException::Programmer, "Cannot calculate distance to sun!",
		// _FILEINFO_);
		// }
		//
		// if ( "iof" == g_iofCorrection.toLower() ) {
		// // Units of I/F
		// // TODO: Note, we do not have a correct g_iofScale (== 1 right now). This was provided for
		// // Hayabusa AMICA's v-band and all other bands were normalized according to this value. We do
		// // not have this value for Hayabusa2 ONC. Need to correct this value.
		// g_iof = pi_c() * (g_solarDist * g_solarDist) *
		// g_iofScale / g_solarFlux / g_texp;
		// g_units = "I over F";
		// }

		// Calibrate!
		try {
		p.Progress()->SetText("Calibrating Hayabusa2 Cube");
		p.StartProcess(Calibrate);
		}
		catch (IException &ie) {
		throw IException(ie, IException::Programmer,
		"Radiometric calibration failed!", _FILEINFO_);
		}

		// Log calibration activity performed so far
		PvlGroup calibrationLog("RadiometricCalibration");
		calibrationLog.addKeyword(PvlKeyword("SoftwareName", hyb2cal_program));
		calibrationLog.addKeyword(PvlKeyword("SoftwareVersion", hyb2cal_version));
		calibrationLog.addKeyword(PvlKeyword("ProcessDate", hyb2cal_runtime));
		calibrationLog.addKeyword(PvlKeyword("CalibrationFile", calfile));
		calibrationLog.addKeyword(PvlKeyword("FlatFieldFile", flatfile.originalPath()
		+ "/" + flatfile.name()));
		calibrationLog.addKeyword(PvlKeyword("CompressionFactor", toString(g_compfactor, 2)));

		// Parameters
		PvlKeyword bn("Bias_Bn");
		bn.addValue(toString(g_b0, 8));
		bn.addValue(toString(g_b1, 8));
		bn.addValue(toString(g_b2, 8));
		calibrationLog.addKeyword(bn);

		PvlKeyword bnae("Bias_AECorrection");
		bnae.addValue(toString(g_bae0, 8));
		bnae.addValue(toString(g_bae1, 8));
		calibrationLog.addKeyword(bnae);

		PvlKeyword linearityCoefs("Linearity");
		linearityCoefs.addValue(toString(g_L0,8));
		linearityCoefs.addValue(toString(g_L1,8));
		linearityCoefs.addValue(toString(g_L2,8));
		calibrationLog.addKeyword(linearityCoefs);



		static double g_AEtemperature(0.0);

		static double g_CCD_T_temperature(0.0);
		static double g_ECT_T_temperature(0.0);

		calibrationLog.addKeyword(PvlKeyword("Bias_AETemp", toString(g_AEtemperature, 16)));


		switch (g_instrument) {

		case InstrumentType::ONCT:
		calibrationLog.addKeyword(PvlKeyword("Bias_CCDTemp", toString(g_CCD_T_temperature, 16)));
		calibrationLog.addKeyword(PvlKeyword("Bias_ECTTemp", toString(g_ECT_T_temperature, 16)));
		break;

		case InstrumentType::ONCW1:
		calibrationLog.addKeyword(PvlKeyword("Bias_CCDTemp", toString(g_CCD_W1_temperature, 16)));
		calibrationLog.addKeyword(PvlKeyword("Bias_ECTTemp", toString(g_ECT_W1_temperature, 16)));
		break;

		case InstrumentType::ONCW2:
		calibrationLog.addKeyword(PvlKeyword("Bias_CCDTemp", toString(g_CCD_W2_temperature, 16)));
		calibrationLog.addKeyword(PvlKeyword("Bias_ECTTemp", toString(g_ECT_W2_temperature, 16)));
		break;
		}



		calibrationLog.addKeyword(PvlKeyword("Bias", toString(g_bias, 16), "DN"));
		calibrationLog.addKeyword(PvlKeyword("Smear_Tvct", toString(g_Tvct, 16)));
		calibrationLog.addKeyword(PvlKeyword("Smear_texp", toString(g_texp, 16)));



		calibrationLog.addKeyword(PvlKeyword("CalibrationUnits", g_iofCorrection));
		calibrationLog.addKeyword(PvlKeyword("RadianceStandard", toString(g_v_standard, 16)));
		calibrationLog.addKeyword(PvlKeyword("RadianceScaleFactor", toString(g_iofScale, 16)));
		calibrationLog.addKeyword(PvlKeyword("SolarDistance", toString(g_solarDist, 16), "AU"));
		calibrationLog.addKeyword(PvlKeyword("SolarFlux", toString(g_solarFlux, 16)));
		calibrationLog.addKeyword(PvlKeyword("IOFFactor", toString(g_iof, 16)));
		calibrationLog.addKeyword(PvlKeyword("Units", g_units));

		// Write Calibration group to output file
		ocube->putGroup(calibrationLog);
		Application::Log(calibrationLog);
		//configFile.clear();
		p.EndProcess();

		}


		/**
		* @brief Determine name of flat field file to apply
		* @author 2016-03-30 Kris Becker
		* @param filter Name of ONC filter
		* @return FileName Path and name of flat file file
		* @internal
		* @history 2017-07-27 Ian Humphrey & Kaj Williams - Adapted from amicacal.
		*/
		FileName DetermineFlatFieldFile(const QString &filter) {

		QString fileName = "$hayabusa2/calibration/flatfield/";
		// FileName consists of binned/notbinned, camera, and filter
		fileName += "flat_" + filter.toLower() + ".cub";
		FileName final(fileName);
		return final;
		}



		/**
		* @brief Loads the calibration variables into the program.
		* @param config QString Name of the calibration file to load.
		*
		* Loads g_b0-g_b2,g_bae0-g_bae1,g_d0-g_d1
		*
		*
		*/
		QString loadCalibrationVariables(const QString &config) {



		// UserInterface& ui = Application::GetUserInterface();

		FileName calibFile(config);
		if ( config.contains("?") ) calibFile = calibFile.highestVersion();

		// Pvl configFile;
		g_configFile.read(calibFile.expanded());

		// Load the groups
		PvlGroup &Bias = g_configFile.findGroup("Bias");
		PvlGroup &DarkCurrent = g_configFile.findGroup("DarkCurrent");
		PvlGroup &Smear = g_configFile.findGroup("SmearRemoval");
		PvlGroup &solar = g_configFile.findGroup("SOLARFLUX");
		PvlGroup &linearity = g_configFile.findGroup("Linearity");
		// PvlGroup &iof = g_configFile.findGroup("RAD");

		// Load Smear Removal Variables
		g_Tvct = Smear["Tvct"];



		// Load DarkCurrent variables and calculate the dark current
		g_d0 = DarkCurrent["D"][0].toDouble();
		g_d1 = DarkCurrent["D"][1].toDouble();
		double CCDTemp(0.0);

		switch (g_instrument) {
		case InstrumentType::ONCT:
		CCDTemp = g_CCD_T_temperature;
		break;
		case InstrumentType::ONCW1:
		CCDTemp = g_CCD_W1_temperature;
		break;
		case InstrumentType::ONCW2:
		CCDTemp = g_CCD_W2_temperature;
		break;
		default:
		CCDTemp = g_CCD_T_temperature;
		}

		g_darkCurrent = g_texpexp(0.10CCDTemp+0.52);

		// Load Bias variables
		g_b0 = Bias["B"][0].toDouble();
		g_b1 = Bias["B"][1].toDouble();
		g_b2 = Bias["B"][2].toDouble();
		g_bae0 = Bias["B_AE"][0].toDouble();
		g_bae1 = Bias["B_AE"][1].toDouble();



		// Compute BIAS correction factor (it's a constant so do it once!)
		g_bias = g_b0+g_b1g_CCD_T_temperature+g_b2g_ECT_T_temperature;
		g_bias = (g_bae0 + g_bae1g_AEtemperature); //correction factor
		qDebug() << "Bias: " << g_bias;

		// Load the Solar Flux for the specific filter
		g_solarFlux = solar[g_filter.toLower()];

		//Load the linearity variables
		g_L0 = linearity["L"][0].toDouble();
		g_L1 = linearity["L"][1].toDouble();
		g_L2 = linearity["L"][2].toDouble();


		// radiance = g_v_standard * g_iofScale
		// iof = radiance * pi *dist_au^2
		// g_iofScale = iof[g_filter];

		return ( calibFile.original() );
		}


		/**
		* @brief Apply radiometric correction to each line of an AMICA image.
		* @author 2016-03-30 Kris Becker
		* @param in Raw image and flat field
		* @param out Radometrically corrected image
		* @internal
		* @history 2017-07-2017 Ian Humphrey & Kaj Williams - Adapted from amicacal.
		*/
		void Calibrate(vector<Buffer >& in, vector<Buffer >& out) {

		Buffer& imageIn = *in[0];
		Buffer& flatField = *in[1];
		Buffer& imageOut = *out[0];

		int pixelsToNull = 0;

		// Note that is isn't currently tested, as we do not have a test with a hayabusa2 image that
		// has been on-board cropped.
		int currentSample = imageIn.Sample();
		int alphaSample = alpha->AlphaSample(currentSample);

		if ( (alphaSample <= pixelsToNull) \|\| (alphaSample >= (1024 - pixelsToNull ) ) ) {

		for (int i = 0; i < imageIn.size(); i++ ) {
		imageOut[i] = Isis::Null;
		}
		return;
		}


		// TODO: Once smear model and readout time are known, we can add smear correction.
		// Compute smear component here as its a constant for the entire sample
		// double t1 = g_timeRatio/imageIn.size();
		// double b = binning;
		// double c1(1.0); //default if no binning
		//
		// if (binning > 1) {
		// c1 = 1.0/(1.0 + t1((b -1.0)/(2.0b) ) );
		// }
		//
		// double smear = 0;
		// for (int j = 0; j < imageIn.size(); j++ ) {
		// if ( !IsSpecial(imageIn[j]) ) {
		// smear += t1 * ( (imageIn[j] * g_compfactor) - g_bias);
		// }
		// }





		// Iterate over the line space
		for (int i = 0; i < imageIn.size(); i++) {
		//qDebug() << "imageIn: " << imageIn[i];
		imageOut[i] = imageIn[i]*pow(2.0,12-g_bitDepth);
		//qDebug() << "imageIOut: " << imageOut[i];


		// Check for special pixel in input image and pass through
		if ( IsSpecial(imageOut[i]) ) {
		imageOut[i] = imageIn[i];
		continue;
		}

		// Apply compression factor here to raise LOSSY dns to proper response
		imageOut[i] *= g_compfactor;



		// 1) BIAS Removal - Only needed if not on-board corrected
		// TODO: find image with > 1 subimages



		if ( !g_onBoardSmearCorrection ) {


		if ( (imageOut[i] - g_bias) <= 0.0) {
		imageOut[i] = Null;
		continue;
		}
		else {
		imageOut[i] = imageOut[i] - g_bias;
		}
		}
		#if 0
		double linearCorrection;
		linearCorrection = g_L0+g_L1pow(imageOut[i],2.0)+g_L2pow(imageOut[i],3.0);
		imageOut[i]*=linearCorrection;

		#endif


		// DARK Current
		imageOut[i] = imageOut[i] - g_darkCurrent;




		// READOUT Smear Removal - Not needed if on-board corrected. Binning is
		// accounted for in computation of c1 before loop.
		// if (nsubImages <= 1) {
		// imageOut[i] = c1*(imageOut[i] - smear);
		// }
		if (g_onBoardSmearCorrection) {
		double smear = 0;
		for (int j=0;j < imageIn.size();j++) {
		smear += (imageOut[j]/imageIn.size() );
		}
		smear*=g_timeRatio;
		imageOut[i] = imageOut[i] - smear;

		}

		//Linearity Correction
		//In the SIS this adjustment is made just after the bias, but
		//in the Calibration paper it happens just before the flat field correction.

		#if 0
		double linearCorrection;
		linearCorrection = g_L0+g_L1pow(imageOut[i],2.0)+g_L2pow(imageOut[i],3.0);
		qDebug() << "linearCorrection=" << linearCorrection;
		imageOut[i]*=linearCorrection;
		#endif


		// FLATFIELD correction
		// Check for any special pixels in the flat field (unlikely)
		// If we have only one input cube, that means that we do not have a flat-field (W1/W2).
		if (in.size() == 2) {
		// Note that our current flat-fields to not have special pixel values.
		if ( IsSpecial(flatField[i]) ) {
		imageOut[i] = Isis::Null;
		continue;
		}
		else {
		if (flatField[i] != 0) {
		imageOut[i] /= flatField[i];
		}
		}
		}

		// TODO: once the radiance values are known for each band, we can correctly compute I/F.
		// For now, g_iof is 1, so output will be in DNs.
		// 7) I/F or Radiance Conversion (or g_iof might = 1, in which case the output will be in DNs)
		imageOut[i] *= g_iof;

		}


		return;
		}