source: trunk/src/Detection/columns.cc @ 543

Last change on this file since 543 was 541, checked in by MatthewWhiting, 15 years ago

Changing all calls of uint to unsigned int, as there are sometimes compilers that don't know about that typedef. Also added an include call for stdlib.h to fitsHeader.cc so that it knows about calloc.

File size: 25.5 KB
Line 
1// -----------------------------------------------------------------------
2// columns.cc: Define a set of columns for Duchamp output.
3// -----------------------------------------------------------------------
4// Copyright (C) 2006, Matthew Whiting, ATNF
5//
6// This program is free software; you can redistribute it and/or modify it
7// under the terms of the GNU General Public License as published by the
8// Free Software Foundation; either version 2 of the License, or (at your
9// option) any later version.
10//
11// Duchamp is distributed in the hope that it will be useful, but WITHOUT
12// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13// FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
14// for more details.
15//
16// You should have received a copy of the GNU General Public License
17// along with Duchamp; if not, write to the Free Software Foundation,
18// Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307, USA
19//
20// Correspondence concerning Duchamp may be directed to:
21//    Internet email: Matthew.Whiting [at] atnf.csiro.au
22//    Postal address: Dr. Matthew Whiting
23//                    Australia Telescope National Facility, CSIRO
24//                    PO Box 76
25//                    Epping NSW 1710
26//                    AUSTRALIA
27// -----------------------------------------------------------------------
28#include <iostream>
29#include <sstream>
30
31#include <math.h>
32#include <duchamp/duchamp.hh>
33#include <duchamp/fitsHeader.hh>
34#include <vector>
35#include <string>
36#include <duchamp/Detection/detection.hh>
37#include <duchamp/Detection/columns.hh>
38
39namespace duchamp
40{
41
42  namespace Column
43  {
44
45
46    Col::Col()
47    {
48      /// @brief Set the default values for the parameters.
49      this->width=1;
50      this->precision=0;
51      this->name=" ";
52      this->units=" ";
53      this->type=UNKNOWN;
54    }
55
56    Col::~Col(){}
57
58    Col::Col(const Col& c)
59    {
60      operator=(c);
61    }
62
63    Col& Col::operator=(const Col& c)
64    {
65      if(this==&c) return *this;
66      this->width = c.width;
67      this->precision = c.precision;
68      this->name = c.name;
69      this->units = c.units;
70      this->type = c.type;
71      return *this;
72    }
73
74    Col::Col(int num, int prec)
75    {
76      /// A specialised constructor that defines one of the default
77      ///  columns, as defined in the Column namespace
78      /// \param num The number of the column to be constructed.
79      ///            Corresponds to the order of the columns in the const
80      ///            arrays in the Column namespace.
81      /// \param prec The precision to use, if not the default. A value
82      /// <0 or no value given results in the default being used.
83
84      if((num>=0)&&(num<numColumns)){
85        this->width     = defaultWidth[num];
86        if(prec<0) this->precision = defaultPrec[num];
87        else this->precision = prec;
88        this->name      = defaultName[num];
89        this->units     = defaultUnits[num];
90        this->type = COLNAME(num);
91      }
92      else{
93        std::stringstream errmsg;
94        errmsg << "Incorrect value for Col(num) --> num="<<num
95               << ", should be between 0 and " << numColumns-1 << ".\n";
96        duchampError("Column constructor", errmsg.str());
97        this->width = 1;
98        this->precision = 0;
99        this->name = " ";
100        this->units = " ";
101        this->type=UNKNOWN;
102      }
103    }
104
105    Col::Col(std::string name, std::string units, int width, int prec)
106    {
107      /// A generic constructor designed to make a Column other than
108      /// the predefined ones listed in the Column namespace
109      /// \param name The title of the column
110      /// \param units What units the column takes
111      /// \param width The starting width
112      /// \param prec The starting precision
113
114      this->width = width;
115      this->precision = prec;
116      this->name = name;
117      this->units = units;
118      this->type = UNKNOWN;
119    }
120
121    //------------------------------------------------------------
122
123    void Col::changePrec(int p)
124    {
125      /// A direct way to alter the precision without having to use
126      /// Col::upPrec(). If the precision increases, the width will
127      /// increase by the same amount. If it decreases, no change is
128      /// made to the width.
129      /// \param p The new precision value.
130
131      if(p > this->precision) this->width += (p-this->precision);
132      this->precision = p;
133
134    }
135   
136    //------------------------------------------------------------
137    void Col::printTitle(std::ostream &stream)
138    {
139      stream << std::setw(this->width)
140             << std::setfill(' ')
141             << this->name;
142    }
143
144    void Col::printUnits(std::ostream &stream)
145    {
146      stream << std::setw(this->width)
147             << std::setfill(' ')
148             << this->units;
149    }
150 
151    void Col::printDash (std::ostream &stream)
152    {
153      stream << std::setw(this->width)
154             << std::setfill('-')
155             << ""
156             << std::setfill(' ');
157    }
158
159    void Col::printBlank(std::ostream &stream)
160    {
161      stream << std::setw(this->width)
162             << std::setfill(' ')
163             << "";
164    }
165 
166    template <class T> void Col::printEntry(std::ostream &stream, T value)
167    {
168      ///  \param stream Where the printing is done.
169      ///  \param value  The value to be printed.
170
171      stream << std::setprecision(this->precision)
172             << std::setw(this->width)
173             << std::setfill(' ')
174             << value;
175    }
176    template void Col::printEntry<int>(std::ostream &stream, int value);
177    template void Col::printEntry<long>(std::ostream &stream, long value);
178    template void Col::printEntry<unsigned>(std::ostream &stream, unsigned value);
179    template void Col::printEntry<float>(std::ostream &stream, float value);
180    template void Col::printEntry<double>(std::ostream &stream, double value);
181    template void Col::printEntry<std::string>(std::ostream &stream, std::string value);
182
183 
184    bool Col::doCol(std::string tableType, bool flagFint)
185    {
186      ///  @details Uses the info in the isFile etc arrays to
187      ///  determine whether a given column, referenced by the
188      ///  enumeration COLNAME, is used for a given table type.
189      /// \param tableType The type of table: one of file, screen, log, votable.
190      /// \param flagFint Whether to use FINT (true) or FTOT
191      /// (false). This defaults to true, so need not be given. It only
192      /// applies to the screen and votable cases -- both are written for the
193      /// results file case.
194      /// \return True if column is used for given table type. False
195      /// otherwise. False if tableType not one of four listed.
196     
197      if(tableType == "file") return isFile[this->type];
198      else if(tableType == "screen"){
199        if(flagFint && this->type==FTOT) return false;
200        if(!flagFint && this->type==FINT) return false;
201        return isScreen[this->type];
202      }
203      else if(tableType == "log") return isLog[this->type];
204      else if(tableType == "votable"){
205        if(flagFint && this->type==FTOT) return false;
206        if(!flagFint && this->type==FINT) return false;
207        return isVOTable[this->type];
208      }
209      else return false;
210    }
211
212 
213
214    std::vector<Col> getFullColSet(std::vector<Detection> &objectList,
215                                   FitsHeader &head)
216    {
217      ///  @details A function that returns a std::vector of Col
218      ///  objects containing information on the columns necessary for
219      ///  output to the results file:
220      ///    Obj#,NAME,X,Y,Z,RA,DEC,VEL,w_RA,w_DEC,w_VEL,F_tot,F_int,F_peak,
221      ///                X1,X2,Y1,Y2,Z1,Z2,Npix,Flag,
222      ///                XAV,YAV,ZAV,XCENT,YCENT,ZCENT,XPEAK,YPEAK,ZPEAK
223      ///
224      ///   Each object in the provided objectList is checked to see if it
225      ///   requires any column to be widened, or for that column to have
226      ///   its precision increased.
227      ///
228      ///   Both Ftot and Fint are provided -- it is up to the calling
229      ///   function to determine which to use.
230      ///
231      /// \param objectList A std::vector list of Detection objects that the
232      /// columns need to fit.
233      /// \param head The FitsHeader object defining the World Coordinate
234      /// System.
235      /// \return A std::vector list of Col definitions.
236
237      std::vector<Col> newset;
238
239      // desired precisions for fluxes, velocities and SNR value
240      int precVel,precFlux,precSNR;
241      if(objectList.size()>0){
242        precVel=objectList[0].getVelPrec();
243        precFlux=objectList[0].getFpeakPrec(); // same as FintPrec at this point.
244        precSNR=objectList[0].getSNRPrec();
245      }
246      else {
247        precVel = prVEL;
248        precFlux = prFLUX;
249        precSNR = prSNR;
250      }
251
252      // set up the default columns
253      newset.push_back( Col(NUM) );
254      newset.push_back( Col(NAME) );
255      newset.push_back( Col(X) );
256      newset.push_back( Col(Y) );
257      newset.push_back( Col(Z) );
258      newset.push_back( Col(RA) );
259      newset.push_back( Col(DEC) );
260      newset.push_back( Col(RAJD) );
261      newset.push_back( Col(DECJD) );
262      newset.push_back( Col(VEL, precVel) );
263      newset.push_back( Col(WRA) );
264      newset.push_back( Col(WDEC) );
265      newset.push_back( Col(W50, precVel) );
266      newset.push_back( Col(W20, precVel) );
267      newset.push_back( Col(WVEL, precVel) );
268      newset.push_back( Col(FINT, precFlux) );
269      newset.push_back( Col(FTOT, precFlux) );
270      newset.push_back( Col(FPEAK, precFlux) );
271      newset.push_back( Col(SNRPEAK, precSNR) );
272      newset.push_back( Col(X1) );
273      newset.push_back( Col(X2) );
274      newset.push_back( Col(Y1) );
275      newset.push_back( Col(Y2) );
276      newset.push_back( Col(Z1) );
277      newset.push_back( Col(Z2) );
278      newset.push_back( Col(NPIX) );
279      newset.push_back( Col(FLAG) );
280      newset.push_back( Col(XAV) );
281      newset.push_back( Col(YAV) );
282      newset.push_back( Col(ZAV) );
283      newset.push_back( Col(XCENT) );
284      newset.push_back( Col(YCENT) );
285      newset.push_back( Col(ZCENT) );
286      newset.push_back( Col(XPEAK) );
287      newset.push_back( Col(YPEAK) );
288      newset.push_back( Col(ZPEAK) );
289
290
291      // Now test each object against each new column
292      for( unsigned int obj = 0; obj < objectList.size(); obj++){
293        std::string tempstr;
294        int tempwidth;
295        float val,minval;
296
297        // Obj#
298        tempwidth = int( log10(objectList[obj].getID()) + 1) + 1;
299        for(int i=newset[NUM].getWidth();i<tempwidth;i++) newset[NUM].widen();
300
301        // Name
302        tempwidth = objectList[obj].getName().size() + 1;
303        for(int i=newset[NAME].getWidth();i<tempwidth;i++) newset[NAME].widen();
304
305        // X position
306        val = objectList[obj].getXcentre() + objectList[obj].getXOffset();
307        if((val<1.)&&(val>0.)){
308          minval = pow(10, -1. * (newset[X].getPrecision()+1));
309          if(val < minval) newset[X].upPrec();
310        }
311        tempwidth = int( log10(val) + 1) + newset[X].getPrecision() + 2;
312        for(int i=newset[X].getWidth();i<tempwidth;i++) newset[X].widen();
313        // Y position
314        val = objectList[obj].getYcentre() + objectList[obj].getYOffset();
315        if((val<1.)&&(val>0.)){
316          minval = pow(10, -1. * (newset[Y].getPrecision()+1));
317          if(val < minval) newset[Y].upPrec();
318        }
319        tempwidth = int( log10(val) + 1) + newset[Y].getPrecision() + 2;
320        for(int i=newset[Y].getWidth();i<tempwidth;i++) newset[Y].widen();
321        // Z position
322        val = objectList[obj].getZcentre() + objectList[obj].getZOffset();
323        if((val<1.)&&(val>0.)){
324          minval = pow(10, -1. * (newset[Z].getPrecision()+1));
325          if((val>0.)&&(val < minval)) newset[Z].upPrec();
326        }
327        tempwidth = int( log10(val) + 1) + newset[Z].getPrecision() + 2;
328        for(int i=newset[Z].getWidth();i<tempwidth;i++) newset[Z].widen();
329
330        if(head.isWCS()){ 
331          // RA -- assign correct title. Check width but should be ok by definition
332          tempstr = head.WCS().lngtyp;
333          newset[RA].setName(tempstr);
334          tempwidth = objectList[obj].getRAs().size() + 1;
335          for(int i=newset[RA].getWidth();i<tempwidth;i++) newset[RA].widen();
336     
337          // Dec -- assign correct title. Check width, should be ok by definition
338          tempstr = head.WCS().lattyp;
339          newset[DEC].setName(tempstr);
340          tempwidth = objectList[obj].getDecs().size() + 1;
341          for(int i=newset[DEC].getWidth();i<tempwidth;i++) newset[DEC].widen();
342
343          // RA decimal degrees -- assign correct title. Check width but should be OK
344          tempstr = head.WCS().lngtyp;
345          newset[RAJD].setName(tempstr);
346          val = objectList[obj].getRA();
347          tempwidth = int( log10(fabs(val)) + 1) + newset[RAJD].getPrecision() + 2;
348          for(int i=newset[RAJD].getWidth();i<tempwidth;i++) newset[RAJD].widen();
349     
350          // Dec decimal degrees -- assign correct title. Check width but should be OK
351          tempstr = head.WCS().lattyp;
352          newset[DECJD].setName(tempstr);
353          val = objectList[obj].getDec();
354          tempwidth = int( log10(fabs(val)) + 1) + newset[DECJD].getPrecision() + 2;
355          for(int i=newset[DECJD].getWidth();i<tempwidth;i++) newset[DECJD].widen();
356
357          // Vel -- check width, title and units.
358          //if(head.isSpecOK()){
359          if(head.canUseThirdAxis()){
360            if(head.WCS().spec < 0)  // if it's not a spectral axis
361              newset[VEL].setName( head.WCS().ctype[2] );
362//          else if(head.WCS().restfrq == 0) // using frequency, not velocity
363            else if(head.getSpectralDescription()==duchamp::duchampSpectralDescription[FREQUENCY]) // using frequency, not velocity
364              newset[VEL].setName("FREQ");
365            if(head.getSpectralUnits().size()>0)
366              newset[VEL].setUnits("[" + head.getSpectralUnits() + "]");
367            tempwidth = newset[VEL].getUnits().size() + 1;
368            for(int i=newset[VEL].getWidth();i<tempwidth;i++) newset[VEL].widen();
369       
370            val = objectList[obj].getVel();
371            if((fabs(val) < 1.)&&(val>0.)){
372              minval = pow(10, -1. * (newset[VEL].getPrecision()+1));
373              if(val < minval) newset[VEL].upPrec();
374            }
375            tempwidth = int(log10(fabs(val)) + 1) + newset[VEL].getPrecision() + 2;
376            if(val<0) tempwidth++;
377            for(int i=newset[VEL].getWidth();i<tempwidth;i++) newset[VEL].widen();
378          }
379
380          // w_RA -- check width & title. leave units for the moment.
381          tempwidth = newset[RA].getUnits().size() + 1;
382          for(int i=newset[RA].getWidth();i<tempwidth;i++) newset[RA].widen();
383          val = objectList[obj].getRAWidth();
384          if((fabs(val) < 1.)&&(val>0.)){
385            minval = pow(10, -1. * (newset[WRA].getPrecision()+1));
386            if(val < minval) newset[WRA].upPrec();
387          }
388          tempwidth = int( log10(fabs(val)) + 1) + newset[WRA].getPrecision() + 2;
389          if(val<0) tempwidth++;
390          for(int i=newset[WRA].getWidth();i<tempwidth;i++) newset[WRA].widen();
391
392          // w_DEC -- check width & title. leave units for the moment.
393          tempwidth = newset[DEC].getUnits().size() + 1;
394          for(int i=newset[DEC].getWidth();i<tempwidth;i++) newset[DEC].widen();
395          val = objectList[obj].getDecWidth();
396          if((fabs(val) < 1.)&&(val>0.)){
397            minval = pow(10, -1. * (newset[WDEC].getPrecision()+1));
398            if(val < minval) newset[WDEC].upPrec();
399          }
400          tempwidth = int( log10(fabs(val)) + 1) + newset[WDEC].getPrecision() + 2;
401          if(val<0) tempwidth++;
402          for(int i=newset[WDEC].getWidth();i<tempwidth;i++) newset[WDEC].widen();
403
404          // w_50 -- check width, title and units.
405          if(head.canUseThirdAxis()){
406            if(head.getSpectralUnits().size()>0)
407              newset[W50].setUnits("[" + head.getSpectralUnits() + "]");
408            tempwidth = newset[W50].getUnits().size() + 1;
409            for(int i=newset[W50].getWidth();i<tempwidth;i++) newset[W50].widen();
410            val = objectList[obj].getW50();
411            if((fabs(val) < 1.)&&(val>0.)){
412              minval = pow(10, -1. * (newset[W50].getPrecision()+1));
413              if(val < minval) newset[W50].upPrec();
414            }
415            tempwidth = int( log10(fabs(val)) + 1) + newset[W50].getPrecision() + 2;
416            if(val<0) tempwidth++;
417            for(int i=newset[W50].getWidth();i<tempwidth;i++) newset[W50].widen();
418          }
419
420          // w_20 -- check width, title and units.
421          if(head.canUseThirdAxis()){
422            if(head.getSpectralUnits().size()>0)
423              newset[W20].setUnits("[" + head.getSpectralUnits() + "]");
424            tempwidth = newset[W20].getUnits().size() + 1;
425            for(int i=newset[W20].getWidth();i<tempwidth;i++)newset[W20].widen();
426            val = objectList[obj].getW20();
427            if((fabs(val) < 1.)&&(val>0.)){
428              minval = pow(10, -1. * (newset[W20].getPrecision()+1));
429              if(val < minval) newset[W20].upPrec();
430            }
431            tempwidth = int( log10(fabs(val)) + 1) + newset[W20].getPrecision() + 2;
432            if(val<0) tempwidth++;
433            for(int i=newset[W20].getWidth();i<tempwidth;i++) newset[W20].widen();
434          }
435
436          // w_Vel -- check width, title and units.
437          if(head.canUseThirdAxis()){
438            if(head.WCS().spec < 0) // if it's not a spectral axis
439              newset[WVEL].setName( std::string("w_") + head.WCS().ctype[2] );
440//          else if(head.WCS().restfrq == 0) // using frequency, not velocity
441            else if(head.getSpectralDescription()==duchamp::duchampSpectralDescription[FREQUENCY]) // using frequency, not velocity
442              newset[WVEL].setName("w_FREQ");
443            if(head.getSpectralUnits().size()>0)
444              newset[WVEL].setUnits("[" + head.getSpectralUnits() + "]");
445            tempwidth = newset[WVEL].getUnits().size() + 1;
446            for(int i=newset[WVEL].getWidth();i<tempwidth;i++)newset[WVEL].widen();
447            tempwidth = newset[WVEL].getName().size() + 1;
448            for(int i=newset[WVEL].getWidth();i<tempwidth;i++) newset[WVEL].widen();
449            val = objectList[obj].getVelWidth();
450            if((fabs(val) < 1.)&&(val>0.)){
451              minval = pow(10, -1. * (newset[WVEL].getPrecision()+1));
452              if(val < minval) newset[WVEL].upPrec();
453            }
454            tempwidth = int( log10(fabs(val)) + 1) + newset[WVEL].getPrecision() + 2;
455            if(val<0) tempwidth++;
456            for(int i=newset[WVEL].getWidth();i<tempwidth;i++) newset[WVEL].widen();
457          }
458
459          // F_int -- check width & units
460          if(head.getIntFluxUnits().size()>0)
461            newset[FINT].setUnits("[" + head.getIntFluxUnits() + "]");
462          tempwidth = newset[FINT].getUnits().size() + 1;
463          for(int i=newset[FINT].getWidth();i<tempwidth;i++) newset[FINT].widen();
464          val = objectList[obj].getIntegFlux();
465          if((fabs(val) < 1.)// &&(val>0.)
466             ){
467            int minprec = int(fabs(log10(fabs(val))))+2;
468            for(int i=newset[FINT].getPrecision();i<minprec;i++) newset[FINT].upPrec();
469          }
470          tempwidth = int( log10(fabs(val)) + 1) + newset[FINT].getPrecision() + 2;
471          if(val<0) tempwidth++;
472          for(int i=newset[FINT].getWidth();i<tempwidth;i++) newset[FINT].widen();
473     
474        }
475
476        // F_tot
477        newset[FTOT].setUnits("[" + head.getFluxUnits() + "]");
478        tempwidth = newset[FTOT].getUnits().size() + 1;
479        for(int i=newset[FTOT].getWidth();i<tempwidth;i++) newset[FTOT].widen();
480        val = objectList[obj].getTotalFlux();
481        //     std::cerr << val << "\n";
482        if((fabs(val) < 1.)// &&(val>0.)
483           ){
484          int minprec = int(fabs(log10(fabs(val))))+2;
485          for(int i=newset[FTOT].getPrecision();i<minprec;i++) newset[FTOT].upPrec();
486        }
487        tempwidth = int( log10(fabs(val)) + 1) + newset[FTOT].getPrecision() + 2;
488        if(val<0) tempwidth++;
489        for(int i=newset[FTOT].getWidth();i<tempwidth;i++) newset[FTOT].widen();
490
491        // F_peak
492        newset[FPEAK].setUnits("[" + head.getFluxUnits() + "]");
493        tempwidth = newset[FPEAK].getUnits().size() + 1;
494        for(int i=newset[FPEAK].getWidth();i<tempwidth;i++) newset[FPEAK].widen();
495        val = objectList[obj].getPeakFlux();
496        if((fabs(val) < 1.)// &&(val>0.)
497           ){
498          int minprec = int(fabs(log10(fabs(val))))+2;
499          for(int i=newset[FPEAK].getPrecision();i<minprec;i++) newset[FPEAK].upPrec();
500        }
501        tempwidth = int( log10(fabs(val)) + 1) + newset[FPEAK].getPrecision() + 2;
502        if(val<0) tempwidth++;
503        for(int i=newset[FPEAK].getWidth();i<tempwidth;i++) newset[FPEAK].widen();
504
505        // S/N_peak
506        val = objectList[obj].getPeakSNR();
507        if((fabs(val) < 1.)&&(val>0.)){
508          minval = pow(10, -1. * (newset[SNRPEAK].getPrecision()+1));
509          if(val < minval) newset[SNRPEAK].upPrec();
510        }
511        tempwidth = int( log10(fabs(val)) + 1) + newset[SNRPEAK].getPrecision() +2;
512        if(val<0) tempwidth++;
513        for(int i=newset[SNRPEAK].getWidth();i<tempwidth;i++) newset[SNRPEAK].widen();
514
515        // X1 position
516        val = objectList[obj].getXmin() + objectList[obj].getXOffset();
517        tempwidth = int( log10(val) + 1) + newset[X1].getPrecision() + 1;
518        for(int i=newset[X1].getWidth();i<tempwidth;i++) newset[X1].widen();
519        // X2 position
520        val = objectList[obj].getXmax() + objectList[obj].getXOffset();
521        tempwidth = int( log10(val) + 1) + newset[X2].getPrecision() + 1;
522        for(int i=newset[X2].getWidth();i<tempwidth;i++) newset[X2].widen();
523        // Y1 position
524        val = objectList[obj].getYmin() + objectList[obj].getYOffset();
525        tempwidth = int( log10(val) + 1) + newset[Y1].getPrecision() + 1;
526        for(int i=newset[Y1].getWidth();i<tempwidth;i++) newset[Y1].widen();
527        // Y2 position
528        val = objectList[obj].getYmax() + objectList[obj].getYOffset();
529        tempwidth = int( log10(val) + 1) + newset[Y2].getPrecision() + 1;
530        for(int i=newset[Y2].getWidth();i<tempwidth;i++) newset[Y2].widen();
531        // Z1 position
532        val = objectList[obj].getZmin() + objectList[obj].getZOffset();
533        tempwidth = int( log10(val) + 1) + newset[Z1].getPrecision() + 1;
534        for(int i=newset[Z1].getWidth();i<tempwidth;i++) newset[Z1].widen();
535        // Z2 position
536        val = objectList[obj].getZmax() + objectList[obj].getZOffset();
537        tempwidth = int( log10(val) + 1) + newset[Z2].getPrecision() + 1;
538        for(int i=newset[Z2].getWidth();i<tempwidth;i++) newset[Z2].widen();
539
540        // Npix
541        tempwidth = int( log10(objectList[obj].getSize()) + 1) +
542          newset[NPIX].getPrecision() + 1;
543        for(int i=newset[NPIX].getWidth();i<tempwidth;i++) newset[NPIX].widen();
544   
545        // average X position
546        val = objectList[obj].getXAverage() + objectList[obj].getXOffset();
547        if((val<1.)&&(val>0.)){
548          minval = pow(10, -1. * (newset[XAV].getPrecision()+1));
549          if(val < minval) newset[XAV].upPrec();
550        }
551        tempwidth = int( log10(val) + 1) + newset[XAV].getPrecision() + 2;
552        for(int i=newset[XAV].getWidth();i<tempwidth;i++) newset[XAV].widen();
553        // average Y position
554        val = objectList[obj].getYAverage() + objectList[obj].getYOffset();
555        if((val<1.)&&(val>0.)){
556          minval = pow(10, -1. * (newset[YAV].getPrecision()+1));
557          if(val < minval) newset[YAV].upPrec();
558        }
559        tempwidth = int( log10(val) + 1) + newset[YAV].getPrecision() + 2;
560        for(int i=newset[YAV].getWidth();i<tempwidth;i++) newset[YAV].widen();
561        // average Z position
562        val = objectList[obj].getZAverage() + objectList[obj].getZOffset();
563        if((val<1.)&&(val>0.)){
564          minval = pow(10, -1. * (newset[ZAV].getPrecision()+1));
565          if((val>0.)&&(val < minval)) newset[ZAV].upPrec();
566        }
567        tempwidth = int( log10(val) + 1) + newset[ZAV].getPrecision() + 2;
568        for(int i=newset[ZAV].getWidth();i<tempwidth;i++) newset[ZAV].widen();
569   
570        // X position of centroid
571        val = objectList[obj].getXCentroid() + objectList[obj].getXOffset();
572        if((val<1.)&&(val>0.)){
573          minval = pow(10, -1. * (newset[XCENT].getPrecision()+1));
574          if(val < minval) newset[XCENT].upPrec();
575        }
576        tempwidth = int( log10(val) + 1) + newset[XCENT].getPrecision() + 2;
577        for(int i=newset[XCENT].getWidth();i<tempwidth;i++) newset[XCENT].widen();
578        // Y position of centroid
579        val = objectList[obj].getYCentroid() + objectList[obj].getYOffset();
580        if((val<1.)&&(val>0.)){
581          minval = pow(10, -1. * (newset[YCENT].getPrecision()+1));
582          if(val < minval) newset[YCENT].upPrec();
583        }
584        tempwidth = int( log10(val) + 1) + newset[YCENT].getPrecision() + 2;
585        for(int i=newset[YCENT].getWidth();i<tempwidth;i++) newset[YCENT].widen();
586        // Z position of centroid
587        val = objectList[obj].getZCentroid() + objectList[obj].getZOffset();
588        if((val<1.)&&(val>0.)){
589          minval = pow(10, -1. * (newset[ZCENT].getPrecision()+1));
590          if((val>0.)&&(val < minval)) newset[ZCENT].upPrec();
591        }
592        tempwidth = int( log10(val) + 1) + newset[ZCENT].getPrecision() + 2;
593        for(int i=newset[ZCENT].getWidth();i<tempwidth;i++) newset[ZCENT].widen();
594   
595        // X position of peak flux
596        val = objectList[obj].getXPeak() + objectList[obj].getXOffset();
597        if((val<1.)&&(val>0.)){
598          minval = pow(10, -1. * (newset[XPEAK].getPrecision()+1));
599          if(val < minval) newset[XPEAK].upPrec();
600        }
601        tempwidth = int( log10(val) + 1) + newset[XPEAK].getPrecision() + 2;
602        for(int i=newset[XPEAK].getWidth();i<tempwidth;i++) newset[XPEAK].widen();
603        // Y position of peak flux
604        val = objectList[obj].getYPeak() + objectList[obj].getYOffset();
605        if((val<1.)&&(val>0.)){
606          minval = pow(10, -1. * (newset[YPEAK].getPrecision()+1));
607          if(val < minval) newset[YPEAK].upPrec();
608        }
609        tempwidth = int( log10(val) + 1) + newset[YPEAK].getPrecision() + 2;
610        for(int i=newset[YPEAK].getWidth();i<tempwidth;i++) newset[YPEAK].widen();
611        // Z position of peak flux
612        val = objectList[obj].getZPeak() + objectList[obj].getZOffset();
613        if((val<1.)&&(val>0.)){
614          minval = pow(10, -1. * (newset[ZPEAK].getPrecision()+1));
615          if((val>0.)&&(val < minval)) newset[ZPEAK].upPrec();
616        }
617        tempwidth = int( log10(val) + 1) + newset[ZPEAK].getPrecision() + 2;
618        for(int i=newset[ZPEAK].getWidth();i<tempwidth;i++) newset[ZPEAK].widen();
619      }
620
621      return newset;
622
623    }
624
625    std::vector<Col> getLogColSet(std::vector<Detection> &objectList,
626                                  FitsHeader &head)
627    {
628      ///  @details A function that returns a std::vector of Col
629      ///    objects containing information on the columns necessary
630      ///    for logfile output:
631      ///    Obj#,X,Y,Z,F_tot,F_peak,X1,X2,Y1,Y2,Z1,Z2,Npix
632      ///
633      ///   Each object in the provided objectList is checked to see if
634      ///    it requires any column to be widened, or for that column to
635      ///    have its precision increased.
636      ///
637      /// \param objectList A std::vector list of Detection objects that the columns need to fit.
638      /// \param head The FitsHeader object defining the World Coordinate System.
639      /// \return A std::vector list of Col definitions.
640
641      std::vector<Col> newset,tempset;
642 
643      // set up the default columns:
644      //  get from FullColSet, and select only the ones we want.
645      tempset = getFullColSet(objectList,head);
646
647      newset.push_back( tempset[NUM] );
648      newset.push_back( tempset[X] );
649      newset.push_back( tempset[Y] );
650      newset.push_back( tempset[Z] );
651      newset.push_back( tempset[FTOT] );
652      newset.push_back( tempset[FPEAK] );
653      newset.push_back( tempset[SNRPEAK] );
654      newset.push_back( tempset[X1] );
655      newset.push_back( tempset[X2] );
656      newset.push_back( tempset[Y1] );
657      newset.push_back( tempset[Y2] );
658      newset.push_back( tempset[Z1] );
659      newset.push_back( tempset[Z2] );
660      newset.push_back( tempset[NPIX] );
661
662      return newset;
663
664    }
665
666  }
667
668}
Note: See TracBrowser for help on using the repository browser.