USGS

Isis 3.0 Application Source Code Reference

Home

HiCalUtil.h

Go to the documentation of this file.
00001 #ifndef HiCalUtil_h
00002 #define HiCalUtil_h
00003 /**
00004  * @file
00005  * $Revision: 5086 $
00006  * $Date: 2013-03-15 11:35:41 -0700 (Fri, 15 Mar 2013) $
00007  *
00008  *   Unless noted otherwise, the portions of Isis written by the USGS are
00009  *   public domain. See individual third-party library and package descriptions
00010  *   for intellectual property information, user agreements, and related
00011  *   information.
00012  *
00013  *   Although Isis has been used by the USGS, no warranty, expressed or
00014  *   implied, is made by the USGS as to the accuracy and functioning of such
00015  *   software and related material nor shall the fact of distribution
00016  *   constitute any such warranty, and no responsibility is assumed by the
00017  *   USGS in connection therewith.
00018  *
00019  *   For additional information, launch
00020  *   $ISISROOT/doc//documents/Disclaimers/Disclaimers.html
00021  *   in a browser or see the Privacy & Disclaimers page on the Isis website,
00022  *   http://isis.astrogeology.usgs.gov, and the USGS privacy and disclaimers on
00023  *   http://www.usgs.gov/privacy.html.
00024  */
00025 
00026 #include <cmath>
00027 #include <string>
00028 #include <vector>
00029 #include <iostream>
00030 #include <sstream>
00031 
00032 #include "HiCalTypes.h"
00033 #include "DbProfile.h"
00034 #include "IString.h"
00035 #include "Statistics.h"
00036 #include "PvlKeyword.h"
00037 #include "NumericalApproximation.h"
00038 #include "FileName.h"
00039 #include "Pvl.h"
00040 #include "IException.h"
00041 
00042 namespace Isis {
00043 
00044  template <typename T> inline T MIN(const T &A, const T &B) {
00045    if ( A < B ) { return (A); }
00046    else         { return (B); }
00047  }
00048 
00049  template <typename T> inline T MAX(const T &A, const T &B) {
00050    if ( A > B ) { return (A); }
00051    else         { return (B); }
00052  }
00053 
00054 /**
00055  * @brief Counts number of valid pixels in vector
00056  * @param v Vector to inspect
00057  * @return int Number valid pixels in vector
00058  */
00059 inline int ValidCount(const HiVector &v) {
00060   int n = 0;
00061   for (int i = 0 ; i < v.dim() ; i++) {
00062     if (!IsSpecial(v[i])) n++;
00063   }
00064   return (n);
00065 }
00066 
00067 /**
00068  * @brief Counts number of invalid pixels in vector
00069  * @param v Vector to inspect
00070  * @return int Number invalid (special) pixels in vector
00071  */
00072 inline int InValidCount(const HiVector &v) {
00073   int n = 0;
00074   for (int i = 0 ; i < v.dim() ; i++) {
00075     if (IsSpecial(v[i])) n++;
00076   }
00077   return (n);
00078 }
00079 
00080 /**
00081  * @brief Convert HiRISE Cpmm number to Ccd number
00082  *
00083  * @param cpmm Cpmm number
00084  */
00085 inline int CpmmToCcd(int cpmm) {
00086   const int cpmm2ccd[] = {0,1,2,3,12,4,10,11,5,13,6,7,8,9};
00087   if ( (cpmm < 0) || (cpmm >= (int)(sizeof(cpmm2ccd)/sizeof(int))) ) {
00088     QString mess = "CpmmToCdd: Bad CpmmNumber (" + toString(cpmm) + ")";
00089     throw IException(IException::User, mess, _FILEINFO_);
00090   }
00091   return (cpmm2ccd[cpmm]);
00092 }
00093 
00094 
00095 /**
00096  * @brief Convert HiRISE Ccd number to string filter name
00097  *
00098  * @param ccd Ccd number of device
00099  */
00100 inline QString CcdToFilter(int ccd) {
00101   if ( (ccd < 0) || (ccd > 13) ) {
00102     QString mess = "CcdToFilter: Bad Ccd Number (" + QString(ccd) + ")";
00103     throw IException(IException::User, mess, _FILEINFO_);
00104   }
00105 
00106   QString filter;
00107   if ( ccd <= 9 )     { filter = "RED"; }
00108   else if (ccd <= 11) { filter = "IR";  }
00109   else                { filter = "BG";  }
00110   return (filter);
00111 }
00112 
00113 /**
00114  * @brief Crop specified lines from a buffer
00115  *
00116  * This function extracts lines from a buffer and returns a new buffer with the
00117  * extracted lines.
00118  *
00119  * @param m      Buffer to extract lines from
00120  * @param sline  Starting line number (first line is 0)
00121  * @param eline  Last line to extract
00122  *
00123  * @return HiMatrix Buffer containing cropped lines
00124  */
00125 inline HiMatrix cropLines(const HiMatrix &m, int sline, int eline) {
00126   int nlines(eline-sline+1);
00127   HiMatrix mcrop(nlines, m.dim2());
00128   for (int l =  0 ; l < nlines ; l++) {
00129     for (int s = 0 ; s < m.dim2() ; s++) {
00130       mcrop[l][s] = m[l+sline][s];
00131     }
00132   }
00133   return (mcrop);
00134 }
00135 
00136 /**
00137  * @brief Crop specified samples from a buffer
00138  *
00139  * This function extracts samples from a buffer and returns a new buffer with
00140  * the extracted samples.
00141  *
00142  * @param m      Buffer to extract samples from
00143  * @param ssamp  Startling sample (first sample 0)
00144  * @param esamp  Ending sample to extract
00145  *
00146  * @return HiMatrix  Buffer with cropped samples
00147  */
00148 inline HiMatrix cropSamples(const HiMatrix &m, int ssamp, int esamp) {
00149   int nsamps(esamp-ssamp+1);
00150   HiMatrix mcrop(m.dim1(), nsamps);
00151   for (int l =  0 ; l < m.dim1() ; l++) {
00152     for (int s = 0 ; s < nsamps ; s++) {
00153       mcrop[l][s] = m[l][s+ssamp];
00154     }
00155   }
00156   return (mcrop);
00157 }
00158 
00159 /**
00160  * @brief Reduces by averaging specified lines from a buffer
00161  *
00162  * This function produces a vector from a 2-D buffer of averaged lines at each
00163  * vertical sample location
00164  *
00165  * @param m      Buffer to reduce
00166  * @param sline  Starting line number (first line is 0)
00167  * @param eline  Last line to average (-1 means use all lines)
00168  *
00169  * @return HiVector Buffer containing averaged lines
00170  */
00171 inline HiVector averageLines(const HiMatrix &m, int sline = 0, int eline = -1) {
00172   eline = (eline == -1) ? m.dim1() - 1 : eline;
00173   HiVector v(m.dim2());
00174   for (int s =  0 ; s < m.dim2() ; s++) {
00175     Statistics stats;
00176     for (int l = sline ; l <= eline ; l++) {
00177       stats.AddData(&m[l][s], 1);
00178     }
00179     v[s] = stats.Average();
00180   }
00181   return (v);
00182 }
00183 
00184 /**
00185  * @brief Reduces by averaging specified samples from a buffer
00186  *
00187  * This function produces a vector from a 2-D buffer of averaged samples at each
00188  * horizontal line location
00189  *
00190  * @param m      Buffer to reduce
00191  * @param ssamp  Starting sample number (first sample is 0)
00192  * @param esamp  Last sample to average (-1 means use all samples)
00193  *
00194  * @return HiVector Buffer containing averaged samples
00195  */
00196 inline HiVector averageSamples(const HiMatrix &m, int ssamp = 0,
00197                                int esamp = -1) {
00198   esamp = (esamp == -1) ? m.dim2() - 1 : esamp;
00199   HiVector v(m.dim1());
00200   for (int l =  0 ; l < m.dim1() ; l++) {
00201     Statistics stats;
00202     for (int s = ssamp ; s <= esamp ; s++) {
00203       stats.AddData(&m[l][s], 1);
00204     }
00205     v[l] = stats.Average();
00206   }
00207   return (v);
00208 }
00209 
00210 /**
00211  * @brief Find a keyword in a profile using default for non-existant keywords
00212  *
00213  * This template function will extract keyword values from a profile and provide
00214  * the values at the indicated index.  If the keyword does not exist or the
00215  * indicated value at index, the provided default will be used instead.
00216  *
00217  */
00218 template <typename T>
00219   T ConfKey(const DbProfile &conf, const QString &keyname, const T &defval,
00220             int index = 0) {
00221   if (!conf.exists(keyname)) { return (defval); }
00222   if (conf.count(keyname) < index) { return (defval); }
00223   QString iValue(conf.value(keyname, index));
00224   T value = iValue;  // This makes it work with a string?
00225   return (value);
00226 }
00227 
00228 /**
00229  * @brief Helper function to convert values to Integers
00230  *
00231  * @param T Type of value to convert
00232  * @param value Value to convert
00233  *
00234  * @return int Converted value
00235  */
00236 template <typename T> int ToInteger(const T &value) {
00237     return (QString(value).trimmed().toInt());
00238 }
00239 
00240 /**
00241  * @brief Helper function to convert values to doubles
00242  *
00243  * @param T Type of value to convert
00244  * @param value Value to convert
00245  *
00246  * @return double Converted value
00247  */
00248 template <typename T> double ToDouble(const T &value) {
00249     return (QString(value).trimmed().toDouble());
00250 }
00251 
00252 /**
00253  * @brief Helper function to convert values to strings
00254  *
00255  * @param T Type of value to convert
00256  * @param value Value to convert
00257  *
00258  * @return string Converted value
00259  */
00260 template <typename T> QString ToString(const T &value) {
00261     return (toString(value).trimmed());
00262 }
00263 
00264 /**
00265  * @brief Shortened string equality test
00266  *
00267  * @param v1 First value
00268  * @param v2 Second value
00269  *
00270  * @return bool True if they are equal w/o regard to case
00271  */
00272 inline bool IsEqual(const QString &v1, const QString &v2 = "TRUE") {
00273     return (v1.toUpper() == v2.toUpper());
00274 }
00275 
00276 /**
00277  * @brief Determines if the keyword value is the expected value
00278  *
00279  * This function checks the existance of a keyword in a profile, extracts the
00280  * first value and tests it for equivelance to the expected value.  Note this
00281  * test is case insensitive.
00282  *
00283  * @param prof Profile to find the expected keyword in
00284  * @param key  Name of keyword in profile to test
00285  * @param value Value to test for in keyword
00286  *
00287  * @return bool Returns true only if the keyword exists in the profile and it is
00288  *              the value expected.
00289  */
00290 inline bool IsTrueValue(const DbProfile &prof, const QString &key,
00291                    const QString &value = "TRUE") {
00292   if ( prof.exists(key) ) {
00293     return (IsEqual(prof(key), value));
00294   }
00295   return (false);
00296 }
00297 
00298 /**
00299  * @brief Checks profile flag to skip the current Module
00300  *
00301  * This function looks for the keyword \b Debug::SkipModule and checks its
00302  * value. True is returned if the value is TRUE (case insensentive).
00303  *
00304  * @param prof Module profile from config file
00305  *
00306  * @return bool  True if the value of the Debug::SkipModule keyword is TRUE,
00307  *         otherwise it returns false for all other values.
00308  */
00309 inline bool SkipModule(const DbProfile &prof) {
00310   return (IsTrueValue(prof, "Debug::SkipModule", "TRUE"));
00311 }
00312 
00313 
00314 inline HiMatrix appendLines(const HiMatrix &top, const HiMatrix &bottom) {
00315   //  Ensure same number samples
00316   if (top.dim2() != bottom.dim2()) {
00317        std::ostringstream mess;
00318         mess << "Top buffer samples (" << top.dim2()
00319              << ") do not match bottom buffer samples (" << bottom.dim2()
00320              << ")";
00321         throw IException(IException::User, mess.str(), _FILEINFO_);
00322   }
00323 
00324   int nlines(top.dim1()+bottom.dim1());
00325   HiMatrix mat(nlines, top.dim2());
00326   for (int lt = 0 ; lt < top.dim1() ; lt++) {
00327     for (int s = 0 ; s < top.dim2() ; s++) {
00328       mat[lt][s] = top[lt][s];
00329     }
00330   }
00331 
00332   int topl = top.dim1();
00333   for (int lb = 0 ; lb < bottom.dim1() ; lb++) {
00334     for (int s = 0 ; s < top.dim2() ; s++) {
00335       mat[topl+lb][s] = bottom[lb][s];
00336     }
00337   }
00338   return (mat);
00339 }
00340 
00341 inline HiMatrix appendSamples(const HiMatrix &left, const HiMatrix &right) {
00342   //  Ensure same number samples
00343   if (right.dim1() != right.dim1()) {
00344        std::ostringstream mess;
00345         mess << "Left buffer lines (" << left.dim1()
00346              << ") do not match right buffer lines (" << right.dim1()
00347              << ")";
00348         throw IException(IException::User, mess.str(), _FILEINFO_);
00349   }
00350 
00351   int nsamps(left.dim2()+right.dim2());
00352   HiMatrix mat(left.dim1(), nsamps);
00353   for (int ll = 0 ; ll < left.dim1() ; ll++) {
00354     for (int s = 0 ; s < left.dim2() ; s++) {
00355       mat[ll][s] = left[ll][s];
00356     }
00357   }
00358 
00359   int lefts = left.dim2();
00360   for (int lr = 0 ; lr < right.dim1() ; lr++) {
00361     for (int s = 0 ; s < right.dim2() ; s++) {
00362       mat[lr][lefts+s] = right[lr][s];
00363     }
00364   }
00365   return (mat);
00366 }
00367 
00368 /**
00369  * @brief Compute HiRISE line times
00370  *
00371  * This class will compute the time (in seconds) for a HiRISE observation line
00372  * based upon the binning mode and line number.  It is assumed that the first
00373  * line will be time 0, but that is up to the caller to keep that straight.
00374  *
00375  * @author ????-??-?? Unknown
00376  *
00377  * @internal
00378  */
00379 class HiLineTimeEqn {
00380   public:
00381     HiLineTimeEqn() : _bin(1), _ltime(1.0) { }
00382     HiLineTimeEqn(int bin, double ltime) : _bin(bin), _ltime(ltime) { }
00383     virtual ~HiLineTimeEqn() { }
00384 
00385     void setLineTime(double ltime) { _ltime = ltime; }
00386     void setBin(int bin) { _bin = bin; }
00387     double Time(const double line) const {
00388       return (line * (_bin * _ltime * 1.0E-6));
00389     }
00390     double operator()(const double line) const { return (Time(line)); }
00391 
00392   private:
00393     double _bin;
00394     double _ltime;
00395 };
00396 
00397 /**
00398  * @brief Implements (classic) HiRISE temperature equation
00399  *
00400  * This function computes the dark current temperature and returns the results
00401  * in electrons/sec/pixel.
00402  *
00403  * @param temperature Temperature (typically of the FPGA)
00404  * @param napcm2      Dark current for silicone diodes (nano-ampres/cm^2)
00405  * @param px          Pixel size in microns
00406  *
00407  * @return double     Dark current temperature (electrons/sec/pixel)
00408  */
00409 inline double HiTempEqn(const double temperature, const double napcm2 = 2.0,
00410                         const double px = 12.0) {
00411   double temp = temperature + 273.0;
00412   double eg   = 1.1557 - (7.021E-4 * temp * temp) / (1108.0 + temp);
00413   const double K = 1.38E-23;
00414   const double Q = 1.6E-19;
00415   return (napcm2*(px*px)*2.55E7*std::pow(temp,1.5)*std::exp(-eg*Q/2.0/K/temp));
00416 }
00417 
00418 /**
00419  * @brief Rebins a vector to a different size
00420  *
00421  * This function can rebin to both larger and smaller sizes.  It fits the data
00422  * to a cubic spline and then computes the value at the rebin pixel index.  One
00423  * advantage to this approach is that on input, special pixels are ignored - on
00424  * output there will never be special pixels unless there are not enough points
00425  * to conpute the cubic spline on which case this function throws an exception.
00426  *
00427  * @param v  Input vector to rebin
00428  * @param n  Size of the new output vector
00429  *
00430  * @return HiVector
00431  * @history 2008-11-05 Jeannie Walldren Replaced references to
00432  *          DataInterp class with NumericalApproximation.
00433  */
00434 inline HiVector rebin(const HiVector &v, int n) {
00435   if ( n == v.dim() ) { return (v); }
00436   NumericalApproximation nterp(NumericalApproximation::CubicNatural);
00437   double mag((double) v.dim()/ (double) n);
00438 
00439   for ( int i = 0 ; i < v.dim() ; i++ ) {
00440     double x((double) i);
00441     if ( !IsSpecial(v[i]) ) {
00442       nterp.AddData(x,v[i]);
00443     }
00444   }
00445 
00446   // Compute the spline and fill the output vector
00447   HiVector vout(n);
00448   for ( int j = 0 ; j < n ; j++ ) {
00449     double x = (double) j * mag;
00450     vout[j] = nterp.Evaluate(x, NumericalApproximation::NearestEndpoint);
00451   }
00452   return (vout);
00453 }
00454 
00455 /**
00456  * @brief Deletes HiRISE specific BLOBS from cube file
00457  *
00458  * Ths function removes only the HiRISE specific
00459  *
00460  * @param label Input label associated with file from which to remove the HiRISE
00461  *              blobs
00462  */
00463 inline void RemoveHiBlobs(Pvl &label) {
00464   for ( int blob = 0 ; blob < label.objects() ; blob++ ) {
00465     if ( label.object(blob).isNamed("Table") ) {
00466       QString name = label.object(blob)["Name"][0];
00467       if ( name.toLower() == "hirise calibration ancillary" ||
00468            name.toLower() == "hirise calibration image" ||
00469            name.toLower() == "hirise ancillary" ) {
00470         label.deleteObject(blob);
00471         blob--;
00472       }
00473     }
00474   }
00475   return;
00476 }
00477 
00478 /**
00479  * @brief Return the statistics for a vector of data
00480  *
00481  * The default statistic returned is the median of the dataset but can be
00482  * changed with a compile time change.  The vector passed in will be sorted so
00483  * that the median can be determined.  If the vector has an even number of
00484  * elements in it, the average of the two center values will be returned.  The
00485  * array is assumed to be clean data, no ISIS special pixels.
00486  *
00487  * @param data    Vector containing data compute the statistic.  It will be
00488  *                altered (sorted) upon return to the caller.
00489  *
00490  * @return double The median (default) of the data
00491  */
00492 inline double GainLineStat(std::vector<double> &data) {
00493 
00494   //  Check for special conditions
00495   if (data.size() == 0)  return (Null);
00496   if (data.size() == 1) return (data[0]);
00497 
00498 #if defined(USE_AVERAGE)
00499   Statistics stats;
00500   stats.AddData(&data[0], data.size());
00501   return (stat.Average());
00502 #else
00503   std::sort(data.begin(), data.end());
00504   int meanIndex = data.size()/2;
00505   if ((data.size() % 2) == 1)  return data[meanIndex];
00506   return ((data[meanIndex-1]+data[meanIndex])/2.0);
00507 #endif
00508 }
00509 
00510 };     // namespace Isis
00511 #endif