source: trunk/src/RasterEdgeDetector.cpp@ 2736

Last change on this file since 2736 was 2626, checked in by Malte Marquarding, 12 years ago

Bug fix: rtrim doesn't exist in casacore!

File size: 4.0 KB
Line 
1//
2// C++ Implimentation: RasterEdgeDetector
3//
4// Description:
5//
6//
7// Author: Takeshi Nakazato <takeshi.nakazato@nao.ac.jp>, (C) 2012
8//
9// Copyright: See COPYING file that comes with this distribution
10//
11//
12#include <casa/Arrays/Vector.h>
13#include <casa/Arrays/ArrayMath.h>
14#include <casa/Arrays/ArrayIO.h>
15#include <casa/Containers/Record.h>
16//#include <casa/Utilities/GenSort.h>
17
18#include "RasterEdgeDetector.h"
19
20using namespace std ;
21using namespace casa ;
22
23namespace asap {
24
25RasterEdgeDetector::RasterEdgeDetector()
26 : EdgeDetector()
27{}
28
29RasterEdgeDetector::~RasterEdgeDetector()
30{}
31
32Vector<uInt> RasterEdgeDetector::detect()
33{
34 os_.origin(LogOrigin( "RasterEdgeDetector", "detect", WHERE )) ;
35
36 initDetect() ;
37
38 detectGap() ;
39 selection() ;
40
41 os_ << LogIO::DEBUGGING
42 << "Detected " << off_.nelements() << " integrations as OFF" << LogIO::POST ;
43
44 return off_ ;
45}
46
47void RasterEdgeDetector::parseOption( const Record &option )
48{
49 os_.origin(LogOrigin( "RasterEdgeDetector", "parseOption", WHERE )) ;
50
51 String name = "fraction" ;
52 if ( option.isDefined( name ) ) {
53 if ( option.dataType( name ) == TpString ) {
54 String fstr = option.asString( name ) ;
55 if ( fstr == "auto" ) {
56 fraction_ = -1.0 ; // indicates "auto"
57 }
58 else {
59 // should be "xx%" format
60 fstr.substr(0,fstr.size()-1) ;
61 fraction_ = String::toFloat( fstr ) * 0.01 ;
62 }
63 }
64 else {
65 fraction_ = option.asFloat( name ) ;
66 }
67 }
68 else {
69 fraction_ = 0.1 ; // default is 10%
70 }
71 name = "npts" ;
72 if ( option.isDefined( name ) ) {
73 npts_ = option.asInt( name ) ;
74 }
75 else {
76 npts_ = -1 ; // default is -1 (use fraction_)
77 }
78
79 os_ << "OPTION SUMMARY: " << endl
80 << " fraction=" << fraction_ << endl
81 << " npts=" << npts_ << LogIO::POST ;
82}
83
84void RasterEdgeDetector::detectGap()
85{
86 os_.origin(LogOrigin( "RasterEdgeDetector", "detectGap", WHERE )) ;
87
88 uInt n = time_.nelements() ;
89 Vector<Double> tdiff( n-1 ) ;
90 for ( uInt i = 0 ; i < n-1 ; i++ ) {
91 tdiff[i] = time_[i+1] - time_[i] ;
92 }
93 Double tmed = median( tdiff, False, True, False ) ;
94 Double threshold = tmed * 5.0 ;
95 uInt idx = 0 ;
96 tempuInt_[idx++] = 0 ;
97 for ( uInt i = 0 ; i < n-1 ; i++ ) {
98 if ( tdiff[i] > threshold ) {
99 tempuInt_[idx++] = i+1 ;
100 }
101 }
102 if ( tempuInt_[idx-1] != n ) {
103 tempuInt_[idx++] = n ;
104 }
105 gaplist_ = vectorFromTempStorage( idx ) ;
106
107 os_ << LogIO::DEBUGGING
108 << "Detected " << gaplist_.nelements() << " time gaps." << LogIO::POST ;
109}
110
111void RasterEdgeDetector::selection()
112{
113// os_.origin(LogOrigin( "RasterEdgeDetector", "selection", WHERE )) ;
114
115 uInt n = gaplist_.nelements() - 1 ;
116 uInt idx = 0 ;
117 for ( uInt i = 0 ; i < n ; i++ ) {
118 uInt start = gaplist_[i] ;
119 uInt end = gaplist_[i+1] ;
120 selectionPerRow( idx, start, end ) ;
121 }
122 off_ = vectorFromTempStorage( idx ) ;
123 //sort( off_, Sort::Ascending, Sort::QuickSort | Sort::NoDuplicates ) ;
124}
125
126void RasterEdgeDetector::selectionPerRow( uInt &idx,
127 const uInt &start,
128 const uint &end )
129{
130// os_.origin(LogOrigin( "RasterEdgeDetector", "selectionPerRow", WHERE )) ;
131
132 uInt len = end - start ;
133 uInt noff = numOff( len ) ;
134
135 if ( len > noff * 2 ) {
136 uInt n = start + noff ;
137 for ( uInt i = start ; i < n ; i++ ) {
138 tempuInt_[idx++] = i ;
139 }
140 n = end - noff ;
141 for ( uInt i = n ; i < end ; i++ ) {
142 tempuInt_[idx++] = i ;
143 }
144 }
145 else {
146 // add all elements to off position list
147 for ( uInt i = start ; i < end ; i++ ) {
148 tempuInt_[idx++] = i ;
149 }
150 }
151}
152
153uInt RasterEdgeDetector::numOff( const uInt &n )
154{
155 uInt ret = 0 ;
156 if ( npts_ > 0 ) {
157 ret = min( (uInt)npts_, n ) ;
158 }
159 else if ( fraction_ < 0 ) {
160 ret = optimumNumber( n ) ;
161 }
162 else {
163 ret = int( n * fraction_ ) ;
164 }
165 return ((ret > 0) ? ret : 1) ;
166}
167
168uInt RasterEdgeDetector::optimumNumber( const uInt &n )
169{
170 return uInt( sqrt( n + 1 ) ) - 1 ;
171}
172
173} // namespace asap
Note: See TracBrowser for help on using the repository browser.