Changeset 1610 for branches/alma/src
- Timestamp:
- 08/03/09 11:34:53 (15 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/alma/src/Scantable.cpp
r1603 r1610 1308 1308 int oldsize = abcissa.size() ; 1309 1309 double olddnu = abcissa[1] - abcissa[0] ; 1310 int refChan = 0 ; 1311 double frac = 0.0 ; 1312 double wedge = 0.0 ; 1313 double pile = 0.0 ; 1310 //int refChan = 0 ; 1311 //double frac = 0.0 ; 1312 //double wedge = 0.0 ; 1313 //double pile = 0.0 ; 1314 int ichan = 0 ; 1314 1315 double wsum = 0.0 ; 1315 /*** 1316 * ichan = 0 1317 ***/ 1318 //ofs << "olddnu = " << olddnu << ", dnu = " << dnu << endl ; 1319 pile += dnu ; 1320 wedge = olddnu * ( refChan + 1 ) ; 1321 while ( wedge < pile ) { 1322 newspec[0] += olddnu * oldspec[refChan] ; 1323 newflag[0] = newflag[0] || oldflag[refChan] ; 1324 //ofs << "channel " << refChan << " is included in new channel 0" << endl ; 1325 refChan++ ; 1326 wedge += olddnu ; 1327 wsum += olddnu ; 1328 //ofs << "newspec[0] = " << newspec[0] << " wsum = " << wsum << endl ; 1329 } 1330 frac = ( wedge - pile ) / olddnu ; 1331 wsum += ( 1.0 - frac ) * olddnu ; 1332 newspec[0] += ( 1.0 - frac ) * olddnu * oldspec[refChan] ; 1333 newflag[0] = newflag[0] || oldflag[refChan] ; 1334 //ofs << "channel " << refChan << " is partly included in new channel 0" << " with fraction of " << ( 1.0 - frac ) << endl ; 1335 //ofs << "newspec[0] = " << newspec[0] << " wsum = " << wsum << endl ; 1336 newspec[0] /= wsum ; 1337 //ofs << "newspec[0] = " << newspec[0] << endl ; 1338 //ofs << "wedge = " << wedge << ", pile = " << pile << endl ; 1339 1340 /*** 1341 * ichan = 1 - nChan-2 1342 ***/ 1343 for ( int ichan = 1 ; ichan < nChan - 1 ; ichan++ ) { 1344 pile += dnu ; 1345 newspec[ichan] += frac * olddnu * oldspec[refChan] ; 1346 newflag[ichan] = newflag[ichan] || oldflag[refChan] ; 1347 //ofs << "channel " << refChan << " is partly included in new channel " << ichan << " with fraction of " << frac << endl ; 1348 refChan++ ; 1349 wedge += olddnu ; 1350 wsum = frac * olddnu ; 1351 //ofs << "newspec[" << ichan << "] = " << newspec[ichan] << " wsum = " << wsum << endl ; 1352 while ( wedge < pile ) { 1353 newspec[ichan] += olddnu * oldspec[refChan] ; 1354 newflag[ichan] = newflag[ichan] || oldflag[refChan] ; 1355 //ofs << "channel " << refChan << " is included in new channel " << ichan << endl ; 1356 refChan++ ; 1357 wedge += olddnu ; 1358 wsum += olddnu ; 1359 //ofs << "newspec[" << ichan << "] = " << newspec[ichan] << " wsum = " << wsum << endl ; 1360 } 1361 frac = ( wedge - pile ) / olddnu ; 1362 wsum += ( 1.0 - frac ) * olddnu ; 1363 newspec[ichan] += ( 1.0 - frac ) * olddnu * oldspec[refChan] ; 1364 newflag[ichan] = newflag[ichan] || oldflag[refChan] ; 1365 //ofs << "channel " << refChan << " is partly included in new channel " << ichan << " with fraction of " << ( 1.0 - frac ) << endl ; 1366 //ofs << "wedge = " << wedge << ", pile = " << pile << endl ; 1367 //ofs << "newspec[" << ichan << "] = " << newspec[ichan] << " wsum = " << wsum << endl ; 1368 newspec[ichan] /= wsum ; 1369 //ofs << "newspec[" << ichan << "] = " << newspec[ichan] << endl ; 1370 } 1371 1372 /*** 1373 * ichan = nChan-1 1374 ***/ 1375 // NOTE: Assumed that all spectra have the same bandwidth 1376 pile += dnu ; 1377 newspec[nChan-1] += frac * olddnu * oldspec[refChan] ; 1378 newflag[nChan-1] = newflag[nChan-1] || oldflag[refChan] ; 1379 //ofs << "channel " << refChan << " is partly included in new channel " << nChan-1 << " with fraction of " << frac << endl ; 1380 refChan++ ; 1381 wedge += olddnu ; 1382 wsum = frac * olddnu ; 1383 //ofs << "newspec[" << nChan - 1 << "] = " << newspec[nChan-1] << " wsum = " << wsum << endl ; 1384 for ( int jchan = refChan ; jchan < oldsize ; jchan++ ) { 1385 newspec[nChan-1] += olddnu * oldspec[jchan] ; 1386 newflag[nChan-1] = newflag[nChan-1] || oldflag[jchan] ; 1387 wsum += olddnu ; 1388 //ofs << "channel " << jchan << " is included in new channel " << nChan-1 << " with fraction of " << frac << endl ; 1389 //ofs << "newspec[" << nChan - 1 << "] = " << newspec[nChan-1] << " wsum = " << wsum << endl ; 1390 } 1391 //ofs << "wedge = " << wedge << ", pile = " << pile << endl ; 1392 //ofs << "newspec[" << nChan - 1 << "] = " << newspec[nChan-1] << " wsum = " << wsum << endl ; 1393 newspec[nChan-1] /= wsum ; 1394 //ofs << "newspec[" << nChan - 1 << "] = " << newspec[nChan-1] << endl ; 1316 Vector<Float> z( nChan ) ; 1317 z[0] = abcissa[0] - 0.5 * olddnu + 0.5 * dnu ; 1318 for ( int ii = 1 ; ii < nChan ; ii++ ) 1319 z[ii] = z[ii-1] + dnu ; 1320 Vector<Float> zi( nChan+1 ) ; 1321 Vector<Float> yi( oldsize + 1 ) ; 1322 zi[0] = z[0] - 0.5 * dnu ; 1323 zi[1] = z[0] + 0.5 * dnu ; 1324 for ( int ii = 2 ; ii < nChan ; ii++ ) 1325 zi[ii] = zi[ii-1] + dnu ; 1326 zi[nChan] = z[nChan-1] + 0.5 * dnu ; 1327 yi[0] = abcissa[0] - 0.5 * olddnu ; 1328 yi[1] = abcissa[1] + 0.5 * olddnu ; 1329 for ( int ii = 2 ; ii < oldsize ; ii++ ) 1330 yi[ii] = abcissa[ii-1] + olddnu ; 1331 yi[oldsize] = abcissa[oldsize-1] + 0.5 * olddnu ; 1332 if ( dnu > 0.0 ) { 1333 for ( int ii = 0 ; ii < nChan ; ii++ ) { 1334 double zl = zi[ii] ; 1335 double zr = zi[ii+1] ; 1336 for ( int j = ichan ; j < oldsize ; j++ ) { 1337 double yl = yi[j] ; 1338 double yr = yi[j+1] ; 1339 if ( yl <= zl ) { 1340 if ( yr <= zl ) { 1341 continue ; 1342 } 1343 else if ( yr <= zr ) { 1344 newspec[ii] += oldspec[j] * ( yr - zl ) ; 1345 newflag[ii] = newflag[ii] || oldflag[j] ; 1346 wsum += ( yr - zl ) ; 1347 } 1348 else { 1349 newspec[ii] += oldspec[j] * dnu ; 1350 newflag[ii] = newflag[ii] || oldflag[j] ; 1351 wsum += dnu ; 1352 ichan = j ; 1353 break ; 1354 } 1355 } 1356 else if ( yl < zr ) { 1357 if ( yr <= zr ) { 1358 newspec[ii] += oldspec[j] * ( yr - yl ) ; 1359 newflag[ii] = newflag[ii] || oldflag[j] ; 1360 wsum += ( yr - yl ) ; 1361 } 1362 else { 1363 newspec[ii] += oldspec[j] * ( zr - yl ) ; 1364 newflag[ii] = newflag[ii] || oldflag[j] ; 1365 wsum += ( zr - yl ) ; 1366 ichan = j ; 1367 break ; 1368 } 1369 } 1370 else { 1371 ichan = j - 1 ; 1372 break ; 1373 } 1374 } 1375 newspec[ii] /= wsum ; 1376 wsum = 0.0 ; 1377 } 1378 } 1379 else if ( dnu < 0.0 ) { 1380 for ( int ii = 0 ; ii < nChan ; ii++ ) { 1381 double zl = zi[ii] ; 1382 double zr = zi[ii+1] ; 1383 for ( int j = ichan ; j < oldsize ; j++ ) { 1384 double yl = yi[j] ; 1385 double yr = yi[j+1] ; 1386 if ( yl >= zl ) { 1387 if ( yr >= zl ) { 1388 continue ; 1389 } 1390 else if ( yr >= zr ) { 1391 newspec[ii] += oldspec[j] * abs( yr - zl ) ; 1392 newflag[ii] = newflag[ii] || oldflag[j] ; 1393 wsum += abs( yr - zl ) ; 1394 } 1395 else { 1396 newspec[ii] += oldspec[j] * abs( dnu ) ; 1397 newflag[ii] = newflag[ii] || oldflag[j] ; 1398 wsum += abs( dnu ) ; 1399 ichan = j ; 1400 break ; 1401 } 1402 } 1403 else if ( yl > zr ) { 1404 if ( yr >= zr ) { 1405 newspec[ii] += oldspec[j] * abs( yr - yl ) ; 1406 newflag[ii] = newflag[ii] || oldflag[j] ; 1407 wsum += abs( yr - yl ) ; 1408 } 1409 else { 1410 newspec[ii] += oldspec[j] * abs( zr - yl ) ; 1411 newflag[ii] = newflag[ii] || oldflag[j] ; 1412 wsum += abs( zr - yl ) ; 1413 ichan = j ; 1414 break ; 1415 } 1416 } 1417 else { 1418 ichan = j - 1 ; 1419 break ; 1420 } 1421 } 1422 newspec[ii] /= wsum ; 1423 wsum = 0.0 ; 1424 } 1425 } 1426 // * ichan = 0 1427 // ***/ 1428 // //ofs << "olddnu = " << olddnu << ", dnu = " << dnu << endl ; 1429 // pile += dnu ; 1430 // wedge = olddnu * ( refChan + 1 ) ; 1431 // while ( wedge < pile ) { 1432 // newspec[0] += olddnu * oldspec[refChan] ; 1433 // newflag[0] = newflag[0] || oldflag[refChan] ; 1434 // //ofs << "channel " << refChan << " is included in new channel 0" << endl ; 1435 // refChan++ ; 1436 // wedge += olddnu ; 1437 // wsum += olddnu ; 1438 // //ofs << "newspec[0] = " << newspec[0] << " wsum = " << wsum << endl ; 1439 // } 1440 // frac = ( wedge - pile ) / olddnu ; 1441 // wsum += ( 1.0 - frac ) * olddnu ; 1442 // newspec[0] += ( 1.0 - frac ) * olddnu * oldspec[refChan] ; 1443 // newflag[0] = newflag[0] || oldflag[refChan] ; 1444 // //ofs << "channel " << refChan << " is partly included in new channel 0" << " with fraction of " << ( 1.0 - frac ) << endl ; 1445 // //ofs << "newspec[0] = " << newspec[0] << " wsum = " << wsum << endl ; 1446 // newspec[0] /= wsum ; 1447 // //ofs << "newspec[0] = " << newspec[0] << endl ; 1448 // //ofs << "wedge = " << wedge << ", pile = " << pile << endl ; 1449 1450 // /*** 1451 // * ichan = 1 - nChan-2 1452 // ***/ 1453 // for ( int ichan = 1 ; ichan < nChan - 1 ; ichan++ ) { 1454 // pile += dnu ; 1455 // newspec[ichan] += frac * olddnu * oldspec[refChan] ; 1456 // newflag[ichan] = newflag[ichan] || oldflag[refChan] ; 1457 // //ofs << "channel " << refChan << " is partly included in new channel " << ichan << " with fraction of " << frac << endl ; 1458 // refChan++ ; 1459 // wedge += olddnu ; 1460 // wsum = frac * olddnu ; 1461 // //ofs << "newspec[" << ichan << "] = " << newspec[ichan] << " wsum = " << wsum << endl ; 1462 // while ( wedge < pile ) { 1463 // newspec[ichan] += olddnu * oldspec[refChan] ; 1464 // newflag[ichan] = newflag[ichan] || oldflag[refChan] ; 1465 // //ofs << "channel " << refChan << " is included in new channel " << ichan << endl ; 1466 // refChan++ ; 1467 // wedge += olddnu ; 1468 // wsum += olddnu ; 1469 // //ofs << "newspec[" << ichan << "] = " << newspec[ichan] << " wsum = " << wsum << endl ; 1470 // } 1471 // frac = ( wedge - pile ) / olddnu ; 1472 // wsum += ( 1.0 - frac ) * olddnu ; 1473 // newspec[ichan] += ( 1.0 - frac ) * olddnu * oldspec[refChan] ; 1474 // newflag[ichan] = newflag[ichan] || oldflag[refChan] ; 1475 // //ofs << "channel " << refChan << " is partly included in new channel " << ichan << " with fraction of " << ( 1.0 - frac ) << endl ; 1476 // //ofs << "wedge = " << wedge << ", pile = " << pile << endl ; 1477 // //ofs << "newspec[" << ichan << "] = " << newspec[ichan] << " wsum = " << wsum << endl ; 1478 // newspec[ichan] /= wsum ; 1479 // //ofs << "newspec[" << ichan << "] = " << newspec[ichan] << endl ; 1480 // } 1481 1482 // /*** 1483 // * ichan = nChan-1 1484 // ***/ 1485 // // NOTE: Assumed that all spectra have the same bandwidth 1486 // pile += dnu ; 1487 // newspec[nChan-1] += frac * olddnu * oldspec[refChan] ; 1488 // newflag[nChan-1] = newflag[nChan-1] || oldflag[refChan] ; 1489 // //ofs << "channel " << refChan << " is partly included in new channel " << nChan-1 << " with fraction of " << frac << endl ; 1490 // refChan++ ; 1491 // wedge += olddnu ; 1492 // wsum = frac * olddnu ; 1493 // //ofs << "newspec[" << nChan - 1 << "] = " << newspec[nChan-1] << " wsum = " << wsum << endl ; 1494 // for ( int jchan = refChan ; jchan < oldsize ; jchan++ ) { 1495 // newspec[nChan-1] += olddnu * oldspec[jchan] ; 1496 // newflag[nChan-1] = newflag[nChan-1] || oldflag[jchan] ; 1497 // wsum += olddnu ; 1498 // //ofs << "channel " << jchan << " is included in new channel " << nChan-1 << " with fraction of " << frac << endl ; 1499 // //ofs << "newspec[" << nChan - 1 << "] = " << newspec[nChan-1] << " wsum = " << wsum << endl ; 1500 // } 1501 // //ofs << "wedge = " << wedge << ", pile = " << pile << endl ; 1502 // //ofs << "newspec[" << nChan - 1 << "] = " << newspec[nChan-1] << " wsum = " << wsum << endl ; 1503 // newspec[nChan-1] /= wsum ; 1504 // //ofs << "newspec[" << nChan - 1 << "] = " << newspec[nChan-1] << endl ; 1395 1505 1396 specCol_.put( irow, newspec ) ; 1397 flagsCol_.put( irow, newflag ) ; 1398 1399 // ofs.close() ; 1506 // specCol_.put( irow, newspec ) ; 1507 // flagsCol_.put( irow, newflag ) ; 1508 1509 // // ofs.close() ; 1510 1400 1511 1401 1512 return ;
Note:
See TracChangeset
for help on using the changeset viewer.