[2613] | 1 | //
|
---|
| 2 | // C++ Implimentation: GenericEdgeDetector
|
---|
| 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 <math.h>
|
---|
| 13 |
|
---|
| 14 | #include <casa/Arrays/Vector.h>
|
---|
| 15 | #include <casa/Arrays/Matrix.h>
|
---|
| 16 | #include <casa/Arrays/ArrayMath.h>
|
---|
| 17 | #include <casa/Arrays/ArrayIO.h>
|
---|
| 18 | #include <casa/Utilities/GenSort.h>
|
---|
| 19 |
|
---|
| 20 | #include "GenericEdgeDetector.h"
|
---|
| 21 |
|
---|
| 22 | using namespace std ;
|
---|
| 23 | using namespace casa ;
|
---|
| 24 |
|
---|
| 25 | namespace asap {
|
---|
| 26 |
|
---|
| 27 | GenericEdgeDetector::GenericEdgeDetector()
|
---|
| 28 | : EdgeDetector()
|
---|
| 29 | {}
|
---|
| 30 |
|
---|
| 31 | GenericEdgeDetector::~GenericEdgeDetector()
|
---|
| 32 | {}
|
---|
| 33 |
|
---|
| 34 | Vector<uInt> GenericEdgeDetector::detect()
|
---|
| 35 | {
|
---|
[2615] | 36 | os_.origin(LogOrigin( "GenericEdgeDetector", "detect", WHERE )) ;
|
---|
| 37 |
|
---|
[2613] | 38 | initDetect() ;
|
---|
| 39 |
|
---|
| 40 | topixel() ;
|
---|
| 41 | countup() ;
|
---|
| 42 | thresholding() ;
|
---|
| 43 | labeling() ;
|
---|
| 44 | trimming() ;
|
---|
| 45 | selection() ;
|
---|
| 46 | tuning() ;
|
---|
| 47 |
|
---|
[2615] | 48 | os_ << LogIO::DEBUGGING
|
---|
| 49 | << "Detected " << off_.nelements() << " integrations as OFF" << LogIO::POST ;
|
---|
| 50 |
|
---|
[2613] | 51 | return off_ ;
|
---|
| 52 | }
|
---|
| 53 |
|
---|
| 54 | void GenericEdgeDetector::parseOption( const Record &option )
|
---|
| 55 | {
|
---|
| 56 | os_.origin(LogOrigin( "GenericEdgeDetector", "parseOption", WHERE )) ;
|
---|
| 57 |
|
---|
| 58 | String name = "fraction" ;
|
---|
| 59 | if ( option.isDefined( name ) ) {
|
---|
| 60 | if ( option.dataType( name ) == TpString ) {
|
---|
| 61 | // should be "xx%" format
|
---|
| 62 | String fstr = option.asString( name ) ;
|
---|
[2626] | 63 | fstr = fstr.substr(0,fstr.size()-1) ;
|
---|
[2613] | 64 | fraction_ = String::toFloat( fstr ) * 0.01 ;
|
---|
| 65 | }
|
---|
| 66 | else {
|
---|
| 67 | fraction_ = option.asFloat( name ) ;
|
---|
| 68 | }
|
---|
| 69 | }
|
---|
| 70 | else {
|
---|
| 71 | fraction_ = 0.1 ; // default is 10%
|
---|
| 72 | }
|
---|
| 73 |
|
---|
| 74 | name = "width" ;
|
---|
| 75 | if ( option.isDefined( name ) ) {
|
---|
| 76 | width_ = option.asFloat( name ) ;
|
---|
| 77 | }
|
---|
| 78 | else {
|
---|
| 79 | width_ = 0.5 ; // default is half of median separation
|
---|
| 80 | }
|
---|
| 81 |
|
---|
| 82 | name = "elongated" ;
|
---|
| 83 | if ( option.isDefined( name ) ) {
|
---|
[2614] | 84 | elongated_ = option.asBool( name ) ;
|
---|
[2613] | 85 | }
|
---|
| 86 | else {
|
---|
| 87 | elongated_ = False ; // default is two-dimensional processing
|
---|
| 88 | }
|
---|
| 89 |
|
---|
| 90 | os_ << "OPTION SUMMARY: " << endl
|
---|
| 91 | << " fraction=" << fraction_ << endl
|
---|
| 92 | << " width=" << width_ << endl
|
---|
| 93 | << " elongated=" << (elongated_ ? "True" : "False") << LogIO::POST ;
|
---|
| 94 | }
|
---|
| 95 |
|
---|
| 96 | void GenericEdgeDetector::topixel()
|
---|
| 97 | {
|
---|
[2615] | 98 | // os_.origin(LogOrigin( "GenericEdgeDetector", "topixel", WHERE )) ;
|
---|
[2613] | 99 |
|
---|
| 100 | setup() ;
|
---|
| 101 | Double *pdir_p = new Double[dir_.nelements()] ;
|
---|
| 102 | pdir_.takeStorage( dir_.shape(), pdir_p, TAKE_OVER ) ;
|
---|
| 103 | uInt len = time_.nelements() ;
|
---|
| 104 | Bool b ;
|
---|
| 105 | Double *px_p = pdir_p ;
|
---|
| 106 | Double *py_p = pdir_p + 1 ;
|
---|
| 107 | const Double *dir_p = dir_.getStorage( b ) ;
|
---|
| 108 | const Double *x_p = dir_p ;
|
---|
| 109 | const Double *y_p = dir_p + 1 ;
|
---|
| 110 | for ( uInt i = 0 ; i < len ; i++ ) {
|
---|
[2632] | 111 | *px_p = pcenx_ + ( *x_p - cenx_ ) * decCorr_ / dx_ ;
|
---|
[2613] | 112 | *py_p = pceny_ + ( *y_p - ceny_ ) / dy_ ;
|
---|
| 113 | px_p += 2 ;
|
---|
| 114 | py_p += 2 ;
|
---|
| 115 | x_p += 2 ;
|
---|
| 116 | y_p += 2 ;
|
---|
| 117 | }
|
---|
| 118 | dir_.freeStorage( dir_p, b ) ;
|
---|
| 119 | }
|
---|
| 120 |
|
---|
| 121 | void GenericEdgeDetector::setup()
|
---|
| 122 | {
|
---|
| 123 | os_.origin(LogOrigin( "GenericEdgeDetector", "setup", WHERE )) ;
|
---|
| 124 |
|
---|
[2632] | 125 | Double xmax, xmin, ymax, ymin ;
|
---|
| 126 | minMax( xmin, xmax, dir_.row( 0 ) ) ;
|
---|
| 127 | minMax( ymin, ymax, dir_.row( 1 ) ) ;
|
---|
| 128 | Double wx = ( xmax - xmin ) * 1.1 ;
|
---|
| 129 | Double wy = ( ymax - ymin ) * 1.1 ;
|
---|
| 130 |
|
---|
| 131 | cenx_ = 0.5 * ( xmin + xmax ) ;
|
---|
| 132 | ceny_ = 0.5 * ( ymin + ymax ) ;
|
---|
| 133 | decCorr_ = 1.0 / cos( ceny_ ) ;
|
---|
| 134 |
|
---|
[2613] | 135 | uInt len = time_.nelements() ;
|
---|
| 136 | Matrix<Double> dd = dir_.copy() ;
|
---|
| 137 | for ( uInt i = len-1 ; i > 0 ; i-- ) {
|
---|
[2632] | 138 | dd(0,i) = ( dd(0,i) - dd(0,i-1) ) * decCorr_ ;
|
---|
[2613] | 139 | dd(1,i) = dd(1,i) - dd(1,i-1) ;
|
---|
| 140 | }
|
---|
| 141 | Vector<Double> dr( len-1 ) ;
|
---|
| 142 | Bool b ;
|
---|
| 143 | const Double *dir_p = dd.getStorage( b ) ;
|
---|
| 144 | const Double *x_p = dir_p + 2 ;
|
---|
| 145 | const Double *y_p = dir_p + 3 ;
|
---|
| 146 | for ( uInt i = 0 ; i < len-1 ; i++ ) {
|
---|
| 147 | dr[i] = sqrt( (*x_p) * (*x_p) + (*y_p) * (*y_p) ) ;
|
---|
| 148 | x_p += 2 ;
|
---|
| 149 | y_p += 2 ;
|
---|
| 150 | }
|
---|
| 151 | dir_.freeStorage( dir_p, b ) ;
|
---|
| 152 | Double med = median( dr, False, True, True ) ;
|
---|
[2632] | 153 | dy_ = med * width_ ;
|
---|
| 154 | dx_ = dy_ * decCorr_ ;
|
---|
[2613] | 155 |
|
---|
| 156 | nx_ = uInt( ceil( wx / dx_ ) ) ;
|
---|
| 157 | ny_ = uInt( ceil( wy / dy_ ) ) ;
|
---|
| 158 |
|
---|
| 159 | pcenx_ = 0.5 * Double( nx_ - 1 ) ;
|
---|
| 160 | pceny_ = 0.5 * Double( ny_ - 1 ) ;
|
---|
| 161 |
|
---|
[2615] | 162 | os_ << LogIO::DEBUGGING
|
---|
| 163 | << "rangex=(" << xmin << "," << xmax << ")" << endl
|
---|
| 164 | << "rangey=(" << ymin << "," << ymax << ")" << endl
|
---|
| 165 | << "median separation between pointings: " << med << endl
|
---|
| 166 | << "dx=" << dx_ << ", dy=" << dy_ << endl
|
---|
| 167 | << "wx=" << wx << ", wy=" << wy << endl
|
---|
| 168 | << "nx=" << nx_ << ", ny=" << ny_ << LogIO::POST ;
|
---|
[2613] | 169 | }
|
---|
| 170 |
|
---|
| 171 | void GenericEdgeDetector::countup()
|
---|
| 172 | {
|
---|
| 173 | os_.origin(LogOrigin( "GenericEdgeDetector", "countup", WHERE )) ;
|
---|
| 174 |
|
---|
| 175 | uInt *a_p = new uInt[nx_*ny_] ;
|
---|
| 176 | apix_.takeStorage( IPosition(2,nx_,ny_), a_p, TAKE_OVER ) ;
|
---|
| 177 | apix_ = 0 ;
|
---|
| 178 |
|
---|
| 179 | uInt len = time_.nelements() ;
|
---|
| 180 | uInt ix ;
|
---|
| 181 | uInt iy ;
|
---|
| 182 | // pdir_ is always contiguous
|
---|
| 183 | const Double *pdir_p = pdir_.data() ;
|
---|
| 184 | const Double *px_p = pdir_p ;
|
---|
| 185 | const Double *py_p = pdir_p + 1 ;
|
---|
| 186 | long offset ;
|
---|
| 187 | for ( uInt i = 0 ; i < len ; i++ ) {
|
---|
| 188 | ix = uInt(round( *px_p )) ;
|
---|
| 189 | iy = uInt(round( *py_p )) ;
|
---|
| 190 | offset = ix + iy * nx_ ;
|
---|
| 191 | *(a_p+offset) += 1 ;
|
---|
| 192 | px_p += 2 ;
|
---|
| 193 | py_p += 2 ;
|
---|
| 194 | }
|
---|
| 195 |
|
---|
[2615] | 196 | os_ << LogIO::DEBUGGING
|
---|
| 197 | << "a.max()=" << max(apix_) << ",a.min()=" << min(apix_) << LogIO::POST ;
|
---|
[2613] | 198 | }
|
---|
| 199 |
|
---|
| 200 | void GenericEdgeDetector::thresholding()
|
---|
| 201 | {
|
---|
| 202 | uInt len = apix_.nelements() ;
|
---|
| 203 | uInt *a_p = apix_.data() ;
|
---|
| 204 | for ( uInt i = 0 ; i < len ; i++ ) {
|
---|
| 205 | *a_p = ((*a_p > 0) ? 1 : 0) ;
|
---|
| 206 | a_p++ ;
|
---|
| 207 | }
|
---|
| 208 | }
|
---|
| 209 |
|
---|
| 210 | void GenericEdgeDetector::labeling()
|
---|
| 211 | {
|
---|
| 212 | os_.origin(LogOrigin( "GenericEdgeDetector", "labeling", WHERE )) ;
|
---|
| 213 |
|
---|
| 214 | uInt n = 1 ;
|
---|
| 215 | uInt niter = 0 ;
|
---|
[2615] | 216 | const uInt maxiter = 100 ;
|
---|
[2613] | 217 | while ( n > 0 && niter < maxiter ) {
|
---|
[2617] | 218 | n = _labeling() ;
|
---|
[2615] | 219 | os_ << LogIO::DEBUGGING
|
---|
| 220 | << "cycle " << niter << ": labeled " << n << " pixels" << LogIO::POST ;
|
---|
[2613] | 221 | niter++ ;
|
---|
| 222 | }
|
---|
| 223 | if ( niter == maxiter ) {
|
---|
| 224 | // WARN
|
---|
| 225 | os_ << LogIO::WARN << "labeling not converged before maxiter=" << maxiter << LogIO::POST ;
|
---|
| 226 | }
|
---|
| 227 | }
|
---|
| 228 |
|
---|
[2617] | 229 | uInt GenericEdgeDetector::_labeling()
|
---|
[2613] | 230 | {
|
---|
| 231 | uInt n = 0 ;
|
---|
| 232 | for ( uInt ix = 0 ; ix < nx_ ; ix++ ) {
|
---|
[2617] | 233 | Vector<uInt> v( apix_.row( ix ) ) ;
|
---|
[2613] | 234 | uInt nx = __labeling( v ) ;
|
---|
| 235 | n += nx ;
|
---|
| 236 | }
|
---|
| 237 | for ( uInt iy = 0 ; iy < ny_ ; iy++ ) {
|
---|
[2617] | 238 | Vector<uInt> v( apix_.column( iy ) ) ;
|
---|
[2613] | 239 | uInt ny = __labeling( v ) ;
|
---|
| 240 | n += ny ;
|
---|
| 241 | }
|
---|
| 242 | return n ;
|
---|
| 243 | }
|
---|
| 244 |
|
---|
| 245 | uInt GenericEdgeDetector::__labeling( Vector<uInt> &a )
|
---|
| 246 | {
|
---|
| 247 | uInt n = 0 ;
|
---|
| 248 | if ( allEQ( a, n ) ) {
|
---|
| 249 | return n ;
|
---|
| 250 | }
|
---|
| 251 |
|
---|
| 252 | uInt start ;
|
---|
| 253 | uInt end ;
|
---|
| 254 | _search( start, end, a ) ;
|
---|
| 255 | for ( uInt i = start+1 ; i < end ; i++ ) {
|
---|
| 256 | if ( a[i] == 0 ) {
|
---|
| 257 | a[i] = 1 ;
|
---|
| 258 | n++ ;
|
---|
| 259 | }
|
---|
| 260 | }
|
---|
| 261 | return n ;
|
---|
| 262 | }
|
---|
| 263 |
|
---|
| 264 | void GenericEdgeDetector::_search( uInt &start,
|
---|
| 265 | uInt &end,
|
---|
| 266 | const Vector<uInt> &a )
|
---|
| 267 | {
|
---|
| 268 | uInt n = a.nelements() ;
|
---|
| 269 | start = 0 ;
|
---|
| 270 | while( a[start] == 0 ) {
|
---|
| 271 | start++ ;
|
---|
| 272 | }
|
---|
| 273 | end = n - 1 ;
|
---|
| 274 | while( a[end] == 0 ) {
|
---|
| 275 | end-- ;
|
---|
| 276 | }
|
---|
| 277 | }
|
---|
| 278 |
|
---|
| 279 | void GenericEdgeDetector::trimming()
|
---|
| 280 | {
|
---|
| 281 | os_.origin(LogOrigin( "GenericEdgeDetector", "trimming", WHERE )) ;
|
---|
| 282 |
|
---|
| 283 | const uInt n1 = sum( apix_ ) ;
|
---|
| 284 | const uInt nTrim = uInt(ceil( n1 * fraction_ )) ;
|
---|
[2615] | 285 | os_ << LogIO::DEBUGGING
|
---|
| 286 | << "number of nonzero pixel: " << n1 << endl
|
---|
| 287 | << "number of pixels to be trimmed: " << nTrim << LogIO::POST ;
|
---|
[2613] | 288 | uInt n = 0 ;
|
---|
| 289 | uInt niter = 0 ;
|
---|
| 290 | const uInt maxiter = 100 ;
|
---|
| 291 | if ( !elongated_ ) {
|
---|
| 292 | while ( n < nTrim && niter < maxiter ) {
|
---|
[2617] | 293 | uInt m = _trimming() ;
|
---|
[2615] | 294 | os_ << LogIO::DEBUGGING
|
---|
| 295 | << "cycle " << niter << ": trimmed " << m << " pixels" << LogIO::POST ;
|
---|
[2613] | 296 | n += m ;
|
---|
| 297 | niter++ ;
|
---|
| 298 | }
|
---|
| 299 | }
|
---|
| 300 | else if ( nx_ > ny_ ) {
|
---|
| 301 | os_ << "1D triming along x-axis" << LogIO::POST ;
|
---|
| 302 | while ( n < nTrim && niter < maxiter ) {
|
---|
[2617] | 303 | uInt m = _trimming1DX() ;
|
---|
[2615] | 304 | os_ << LogIO::DEBUGGING
|
---|
| 305 | << "cycle " << niter << ": trimmed " << m << " pixels" << LogIO::POST ;
|
---|
[2613] | 306 | n += m ;
|
---|
| 307 | niter++ ;
|
---|
| 308 | }
|
---|
| 309 | }
|
---|
| 310 | else { // nx_ < ny_
|
---|
| 311 | os_ << "1D triming along y-axis" << LogIO::POST ;
|
---|
| 312 | while ( n < nTrim && niter < maxiter ) {
|
---|
[2617] | 313 | uInt m = _trimming1DY() ;
|
---|
[2615] | 314 | os_ << LogIO::DEBUGGING
|
---|
| 315 | << "cycle " << niter << ": trimmed " << m << " pixels" << LogIO::POST ;
|
---|
[2613] | 316 | n += m ;
|
---|
| 317 | niter++ ;
|
---|
| 318 | }
|
---|
| 319 | }
|
---|
[2615] | 320 | os_ << LogIO::DEBUGGING
|
---|
| 321 | << "number of pixels actually trimmed: " << n << LogIO::POST ;
|
---|
[2613] | 322 |
|
---|
| 323 | if ( niter == maxiter ) {
|
---|
| 324 | // WARN
|
---|
| 325 | os_ << LogIO::WARN << "trimming not converged before maxiter=" << maxiter << LogIO::POST ;
|
---|
| 326 | }
|
---|
| 327 | }
|
---|
| 328 |
|
---|
[2617] | 329 | uInt GenericEdgeDetector::_trimming()
|
---|
[2613] | 330 | {
|
---|
| 331 | uInt n = 0 ;
|
---|
[2617] | 332 | Block<uInt> flatIdxList( apix_.nelements() ) ;
|
---|
[2613] | 333 | uInt start ;
|
---|
| 334 | uInt end ;
|
---|
| 335 | uInt flatIdx ;
|
---|
[2617] | 336 | for ( uInt ix = 0 ; ix < nx_ ; ix++ ) {
|
---|
| 337 | Vector<uInt> v( apix_.row( ix ) ) ;
|
---|
[2613] | 338 | if ( allEQ( v, (uInt)0 ) ) {
|
---|
| 339 | continue ;
|
---|
| 340 | }
|
---|
| 341 | _search( start, end, v ) ;
|
---|
[2617] | 342 | uInt offset = start * nx_ ;
|
---|
[2613] | 343 | flatIdx = offset + ix ;
|
---|
| 344 | flatIdxList[n++] = flatIdx ;
|
---|
| 345 | if ( start != end ) {
|
---|
[2617] | 346 | offset = end * nx_ ;
|
---|
[2613] | 347 | flatIdx = offset + ix ;
|
---|
| 348 | flatIdxList[n++] = flatIdx ;
|
---|
| 349 | }
|
---|
| 350 | }
|
---|
[2617] | 351 | for ( uInt iy = 0 ; iy < ny_ ; iy++ ) {
|
---|
| 352 | Vector<uInt> v( apix_.column( iy ) ) ;
|
---|
[2613] | 353 | if ( allEQ( v, (uInt)0 ) ) {
|
---|
| 354 | continue ;
|
---|
| 355 | }
|
---|
[2617] | 356 | uInt offset = iy * nx_ ;
|
---|
[2613] | 357 | _search( start, end, v ) ;
|
---|
| 358 | flatIdx = offset + start ;
|
---|
| 359 | flatIdxList[n++] = flatIdx ;
|
---|
| 360 | if ( start != end ) {
|
---|
| 361 | flatIdx = offset + end ;
|
---|
| 362 | flatIdxList[n++] = flatIdx ;
|
---|
| 363 | }
|
---|
| 364 | }
|
---|
| 365 | n = genSort( flatIdxList.storage(),
|
---|
| 366 | n,
|
---|
| 367 | Sort::Ascending,
|
---|
| 368 | Sort::QuickSort | Sort::NoDuplicates ) ;
|
---|
[2617] | 369 | Vector<uInt> v( IPosition(1,apix_.nelements()), apix_.data(), SHARE ) ;
|
---|
[2613] | 370 | const uInt *idx_p = flatIdxList.storage() ;
|
---|
| 371 | for ( uInt i = 0 ; i < n ; i++ ) {
|
---|
| 372 | v[*idx_p] = 0 ;
|
---|
| 373 | idx_p++ ;
|
---|
| 374 | }
|
---|
| 375 |
|
---|
| 376 | return n ;
|
---|
| 377 | }
|
---|
| 378 |
|
---|
[2617] | 379 | uInt GenericEdgeDetector::_trimming1DX()
|
---|
[2613] | 380 | {
|
---|
| 381 | uInt n = 0 ;
|
---|
[2617] | 382 | const uInt nx = apix_.nrow() ;
|
---|
[2613] | 383 | Vector<uInt> v1, v2 ;
|
---|
| 384 | uInt ix, jx ;
|
---|
[2617] | 385 | for ( ix = 0 ; ix < nx_ ; ix++ ) {
|
---|
| 386 | v1.reference( apix_.row( ix ) ) ;
|
---|
[2613] | 387 | if ( anyNE( v1, n ) ) break ;
|
---|
| 388 | }
|
---|
[2614] | 389 | for ( jx = nx-1 ; jx > ix ; jx-- ) {
|
---|
[2617] | 390 | v2.reference( apix_.row( jx ) ) ;
|
---|
[2613] | 391 | if ( anyNE( v2, n ) ) break ;
|
---|
| 392 | }
|
---|
| 393 | n += _trimming1D( v1 ) ;
|
---|
| 394 | if ( ix != jx )
|
---|
| 395 | n+= _trimming1D( v2 ) ;
|
---|
| 396 |
|
---|
| 397 | return n ;
|
---|
| 398 | }
|
---|
| 399 |
|
---|
[2617] | 400 | uInt GenericEdgeDetector::_trimming1DY()
|
---|
[2613] | 401 | {
|
---|
| 402 | uInt n = 0 ;
|
---|
[2617] | 403 | const uInt ny = apix_.ncolumn() ;
|
---|
[2613] | 404 | Vector<uInt> v1, v2 ;
|
---|
| 405 | uInt iy, jy ;
|
---|
[2617] | 406 | for ( iy = 0 ; iy < ny_ ; iy++ ) {
|
---|
| 407 | v1.reference( apix_.column( iy ) ) ;
|
---|
[2613] | 408 | if ( anyNE( v1, n ) ) break ;
|
---|
| 409 | }
|
---|
[2614] | 410 | for ( jy = ny-1 ; jy > iy ; jy-- ) {
|
---|
[2617] | 411 | v2.reference( apix_.column( jy ) ) ;
|
---|
[2613] | 412 | if ( anyNE( v2, n ) ) break ;
|
---|
| 413 | }
|
---|
| 414 | n += _trimming1D( v1 ) ;
|
---|
| 415 | if ( iy != jy )
|
---|
| 416 | n+= _trimming1D( v2 ) ;
|
---|
| 417 |
|
---|
| 418 | return n ;
|
---|
| 419 | }
|
---|
| 420 |
|
---|
| 421 | uInt GenericEdgeDetector::_trimming1D( Vector<uInt> &a )
|
---|
| 422 | {
|
---|
| 423 | uInt len = a.nelements() ;
|
---|
| 424 | uInt n = 0 ;
|
---|
| 425 | for ( uInt i = 0 ; i < len ; i++ ) {
|
---|
| 426 | if ( a[i] == 1 ) {
|
---|
| 427 | a[i] = 0 ;
|
---|
| 428 | n++ ;
|
---|
| 429 | }
|
---|
| 430 | }
|
---|
| 431 |
|
---|
| 432 | return n ;
|
---|
| 433 | }
|
---|
| 434 |
|
---|
| 435 | void GenericEdgeDetector::selection()
|
---|
| 436 | {
|
---|
[2615] | 437 | // os_.origin(LogOrigin( "GenericEdgeDetector", "selection", WHERE )) ;
|
---|
[2613] | 438 |
|
---|
| 439 | uInt nrow = pdir_.shape()[1] ;
|
---|
| 440 | const Double *px_p = pdir_.data() ;
|
---|
| 441 | const Double *py_p = px_p + 1 ;
|
---|
| 442 | Vector<uInt> v( IPosition(1,apix_.nelements()), apix_.data(), SHARE ) ;
|
---|
| 443 | uInt n = 0 ;
|
---|
| 444 | for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
|
---|
| 445 | uInt idx = int(round(*px_p)) + int(round(*py_p)) * nx_ ;
|
---|
| 446 | if ( v[idx] == 0 ) {
|
---|
| 447 | tempuInt_[n++] = irow ;
|
---|
| 448 | }
|
---|
| 449 | px_p += 2 ;
|
---|
| 450 | py_p += 2 ;
|
---|
| 451 | }
|
---|
| 452 | off_ = vectorFromTempStorage( n ) ;
|
---|
| 453 | }
|
---|
| 454 |
|
---|
| 455 | void GenericEdgeDetector::tuning()
|
---|
| 456 | {
|
---|
| 457 | os_.origin(LogOrigin( "GenericEdgeDetector", "tuning", WHERE )) ;
|
---|
| 458 |
|
---|
| 459 | const uInt len = off_.nelements() ;
|
---|
| 460 | if ( len == 0 )
|
---|
| 461 | return ;
|
---|
| 462 |
|
---|
[2617] | 463 | Block<uInt> diff( len-1 ) ;
|
---|
[2615] | 464 | for ( uInt i = 0 ; i < len-1 ; i++ ) {
|
---|
| 465 | diff[i] = off_[i+1] - off_[i] ;
|
---|
[2613] | 466 | }
|
---|
| 467 | const uInt threshold = 3 ;
|
---|
| 468 | uInt n = 0 ;
|
---|
| 469 | for ( uInt i = 0 ; i < len ; i++ ) {
|
---|
| 470 | tempuInt_[n++] = off_[i] ;
|
---|
| 471 | }
|
---|
| 472 | for ( uInt i = 1 ; i < len ; i++ ) {
|
---|
[2615] | 473 | uInt ii = i - 1 ;
|
---|
| 474 | if ( diff[ii] != 1 && diff[ii] < threshold ) {
|
---|
| 475 | uInt t = off_[ii]+1 ;
|
---|
[2613] | 476 | uInt u = off_[i] ;
|
---|
| 477 | for ( uInt j = t ; j < u ; j++ ) {
|
---|
[2615] | 478 | os_ << LogIO::DEBUGGING
|
---|
| 479 | << "move " << j << " from ON to OFF" << LogIO::POST ;
|
---|
[2613] | 480 | tempuInt_[n++] = j ;
|
---|
| 481 | }
|
---|
| 482 | }
|
---|
| 483 | }
|
---|
| 484 | if ( n > len ) {
|
---|
| 485 | off_.resize() ;
|
---|
| 486 | off_ = vectorFromTempStorage( n ) ;
|
---|
| 487 | }
|
---|
| 488 | }
|
---|
| 489 |
|
---|
| 490 | } // namespace asap
|
---|