Changeset 699
- Timestamp:
- 10/20/05 16:35:24 (19 years ago)
- Location:
- branches/Release-2-fixes/src
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/Release-2-fixes/src/SDMath.cc
r653 r699 110 110 111 111 112 113 112 SDMemTable* SDMath::frequencyAlignment(const SDMemTable& in, 114 113 const String& refTime, … … 116 115 Bool perFreqID) const 117 116 { 118 // Get frame info from Table 119 117 // Get frame info from Table 120 118 std::vector<std::string> info = in.getCoordInfo(); 121 119 122 // Parse frequency system 123 120 // Parse frequency system 124 121 String systemStr(info[1]); 125 122 String baseSystemStr(info[3]); 126 123 if (baseSystemStr==systemStr) { 127 124 throw(AipsError("You have not set a frequency frame different from the initial - use function set_freqframe")); 128 125 } 129 // 126 130 127 MFrequency::Types freqSystem; 131 128 MFrequency::getType(freqSystem, systemStr); 132 129 133 // Do it134 135 130 return frequencyAlign(in, freqSystem, refTime, method, perFreqID); 136 131 } … … 138 133 139 134 140 CountedPtr<SDMemTable> SDMath::average(const std::vector<CountedPtr<SDMemTable> >& in,141 const Vector<Bool>& mask, Bool scanAv,142 const String& weightStr, Bool alignFreq) const143 // 135 CountedPtr<SDMemTable> 136 SDMath::average(const std::vector<CountedPtr<SDMemTable> >& in, 137 const Vector<Bool>& mask, Bool scanAv, 138 const String& weightStr, Bool alignFreq) const 144 139 // Weighted averaging of spectra from one or more Tables. 145 // 146 { 147 148 // Convert weight type 149 140 { 141 // Convert weight type 150 142 WeightType wtType = NONE; 151 143 convertWeightString(wtType, weightStr, True); 152 144 153 // Create output Table by cloning from the first table 154 145 // Create output Table by cloning from the first table 155 146 SDMemTable* pTabOut = new SDMemTable(*in[0],True); 156 147 if (in.size() > 1) { … … 159 150 } 160 151 } 161 // Setup 162 152 // Setup 163 153 IPosition shp = in[0]->rowAsMaskedArray(0).shape(); // Must not change 164 154 Array<Float> arr(shp); … … 166 156 const Bool useMask = (mask.nelements() == shp(asap::ChanAxis)); 167 157 168 // Columns from Tables 169 158 // Columns from Tables 170 159 ROArrayColumn<Float> tSysCol; 171 160 ROScalarColumn<Double> mjdCol; … … 175 164 ROScalarColumn<Int> scanIDCol; 176 165 177 // Create accumulation MaskedArray. We accumulate for each channel,if,pol,beam 178 // Note that the mask of the accumulation array will ALWAYS remain ALL True. 179 // The MA is only used so that when data which is masked Bad is added to it, 180 // that data does not contribute. 166 // Create accumulation MaskedArray. We accumulate for each 167 // channel,if,pol,beam Note that the mask of the accumulation array 168 // will ALWAYS remain ALL True. The MA is only used so that when 169 // data which is masked Bad is added to it, that data does not 170 // contribute. 181 171 182 172 Array<Float> zero(shp); … … 186 176 MaskedArray<Float> sum(zero,good); 187 177 188 // Counter arrays 189 178 // Counter arrays 190 179 Array<Float> nPts(shp); // Number of points 191 180 nPts = 0.0; … … 193 182 nInc = 1.0; 194 183 195 // Create accumulation Array for variance. We accumulate for 196 // each if,pol,beam, but average over channel. So we need 197 // a shape with one less axis dropping channels. 198 184 // Create accumulation Array for variance. We accumulate for each 185 // if,pol,beam, but average over channel. So we need a shape with 186 // one less axis dropping channels. 199 187 const uInt nAxesSub = shp.nelements() - 1; 200 188 IPosition shp2(nAxesSub); … … 209 197 IPosition pos2(nAxesSub,0); // For indexing 210 198 211 // Time-related accumulators 212 199 // Time-related accumulators 213 200 Double time; 214 201 Double timeSum = 0.0; … … 216 203 Double interval = 0.0; 217 204 218 // To get the right shape for the Tsys accumulator we need to 219 // access a column from the first table. The shape of this 220 // array must not change. Note however that since the TSysSqSum 221 // array is used in a normalization process, and that I ignore the 222 // channel axisreplication of values for now, it loses a dimension205 // To get the right shape for the Tsys accumulator we need to access 206 // a column from the first table. The shape of this array must not 207 // change. Note however that since the TSysSqSum array is used in a 208 // normalization process, and that I ignore the channel axis 209 // replication of values for now, it loses a dimension 223 210 224 211 Array<Float> tSysSum, tSysSqSum; … … 227 214 tSysCol.attach(tabIn,"TSYS"); 228 215 tSysSum.resize(tSysCol.shape(0)); 229 //230 216 tSysSqSum.resize(shp2); 231 217 } 232 tSysSum = 0.0;218 tSysSum = 0.0; 233 219 tSysSqSum = 0.0; 234 220 Array<Float> tSys; 235 221 236 // Scan and row tracking 237 222 // Scan and row tracking 238 223 Int oldScanID = 0; 239 224 Int outScanID = 0; … … 243 228 Int tableStart = 0; 244 229 245 // Source and FreqID 246 230 // Source and FreqID 247 231 String sourceName, oldSourceName, sourceNameStart; 248 232 Vector<uInt> freqID, freqIDStart, oldFreqID; 249 233 250 // Loop over tables 251 234 // Loop over tables 252 235 Float fac = 1.0; 253 236 const uInt nTables = in.size(); 254 237 for (uInt iTab=0; iTab<nTables; iTab++) { 255 238 256 // Should check that the frequency tables don't change if doing FreqAlignment 257 258 // Attach columns to Table 259 239 // Should check that the frequency tables don't change if doing 240 // FreqAlignment 241 242 // Attach columns to Table 260 243 const Table& tabIn = in[iTab]->table(); 261 244 tSysCol.attach(tabIn, "TSYS"); … … 266 249 scanIDCol.attach(tabIn, "SCANID"); 267 250 268 // Loop over rows in Table 269 251 // Loop over rows in Table 270 252 const uInt nRows = in[iTab]->nRow(); 271 253 for (uInt iRow=0; iRow<nRows; iRow++) { 272 273 // Check conformance 274 254 // Check conformance 275 255 IPosition shp2 = in[iTab]->rowAsMaskedArray(iRow).shape(); 276 256 if (!shp.isEqual(shp2)) { 257 delete pTabOut; 277 258 throw (AipsError("Shapes for all rows must be the same")); 278 259 } 279 260 280 // If we are not doing scan averages, make checks for source and 281 // frequency setup and warn if averaging across them 282 261 // If we are not doing scan averages, make checks for source 262 // and frequency setup and warn if averaging across them 283 263 scanIDCol.getScalar(iRow, scanID); 284 264 285 // Get quantities from columns 286 265 // Get quantities from columns 287 266 srcNameCol.getScalar(iRow, sourceName); 288 267 mjdCol.get(iRow, time); … … 291 270 fqIDCol.get(iRow, freqID); 292 271 293 // Initialize first source and freqID 294 272 // Initialize first source and freqID 295 273 if (iRow==0 && iTab==0) { 296 274 sourceNameStart = sourceName; … … 298 276 } 299 277 300 // If we are doing scan averages, see if we are at the end of an 301 // accumulation period (scan). We must check soutce names too, 302 // since we might have two tables with one scan each but different303 // source names; we shouldn't average different sources together 304 278 // If we are doing scan averages, see if we are at the end of 279 // an accumulation period (scan). We must check soutce names 280 // too, since we might have two tables with one scan each but 281 // different source names; we shouldn't average different 282 // sources together 305 283 if (scanAv && ( (scanID != oldScanID) || 306 284 (iRow==0 && iTab>0 && sourceName!=oldSourceName))) { 307 285 308 // Normalize data in 'sum' accumulation array according to weighting scheme 309 310 normalize(sum, sumSq, tSysSqSum, nPts, intSum, wtType, asap::ChanAxis, nAxesSub);311 312 // Get ScanContainer for the first row of this averaged Scan 313 286 // Normalize data in 'sum' accumulation array according to 287 // weighting scheme 288 normalize(sum, sumSq, tSysSqSum, nPts, intSum, wtType, 289 asap::ChanAxis, nAxesSub); 290 291 // Get ScanContainer for the first row of this averaged Scan 314 292 SDContainer scOut = in[iTab]->getSDContainer(rowStart); 315 293 316 // Fill scan container. The source and freqID come from the317 // first row of the first table that went into this average ( 318 //should be the same for all rows in the scan average)294 // Fill scan container. The source and freqID come from the 295 // first row of the first table that went into this average 296 // ( should be the same for all rows in the scan average) 319 297 320 298 Float nR(nAccum); 321 299 fillSDC(scOut, sum.getMask(), sum.getArray(), tSysSum/nR, outScanID, 322 timeSum/nR, intSum, sourceNameStart, freqIDStart); 323 324 // Write container out to Table 325 300 timeSum/nR, intSum, sourceNameStart, freqIDStart); 301 302 // Write container out to Table 326 303 pTabOut->putSDContainer(scOut); 327 328 // Reset accumulators 329 304 305 // Reset accumulators 330 306 sum = 0.0; 331 307 sumSq = 0.0; 332 308 nAccum = 0; 333 // 309 334 310 tSysSum =0.0; 335 311 tSysSqSum =0.0; … … 338 314 nPts = 0.0; 339 315 340 // Increment 341 316 // Increment 342 317 rowStart = iRow; // First row for next accumulation 343 318 tableStart = iTab; // First table for next accumulation 344 sourceNameStart = sourceName; // First source name for next accumulation 319 sourceNameStart = sourceName; // First source name for next 320 // accumulation 345 321 freqIDStart = freqID; // First FreqID for next accumulation 346 // 322 347 323 oldScanID = scanID; 348 outScanID += 1; // Scan ID for next accumulation period 324 outScanID += 1; // Scan ID for next 325 // accumulation period 349 326 } 350 327 351 // Accumulate352 353 accumulate(timeSum, intSum, nAccum, sum, sumSq, nPts, tSysSum, tSysSqSum, 354 tSys,nInc, mask, time, interval, in, iTab, iRow, asap::ChanAxis,328 // Accumulate 329 accumulate(timeSum, intSum, nAccum, sum, sumSq, nPts, 330 tSysSum, tSysSqSum, tSys, 331 nInc, mask, time, interval, in, iTab, iRow, asap::ChanAxis, 355 332 nAxesSub, useMask, wtType); 356 // 357 oldSourceName = sourceName; 358 oldFreqID = freqID; 333 oldSourceName = sourceName; 334 oldFreqID = freqID; 359 335 } 360 336 } 361 337 362 // OK at this point we have accumulation data which is either 363 // - accumulated from all tables into one row 364 // or 365 // - accumulated from the last scan average 366 // 367 // Normalize data in 'sum' accumulation array according to weighting scheme 368 369 normalize(sum, sumSq, tSysSqSum, nPts, intSum, wtType, asap::ChanAxis, nAxesSub); 370 371 // Create and fill container. The container we clone will be from 372 // the last Table and the first row that went into the current 373 // accumulation. It probably doesn't matter that much really... 374 338 // OK at this point we have accumulation data which is either 339 // - accumulated from all tables into one row 340 // or 341 // - accumulated from the last scan average 342 // 343 // Normalize data in 'sum' accumulation array according to weighting 344 // scheme 345 346 normalize(sum, sumSq, tSysSqSum, nPts, intSum, wtType, 347 asap::ChanAxis, nAxesSub); 348 349 // Create and fill container. The container we clone will be from 350 // the last Table and the first row that went into the current 351 // accumulation. It probably doesn't matter that much really... 375 352 Float nR(nAccum); 376 353 SDContainer scOut = in[tableStart]->getSDContainer(rowStart); 377 354 fillSDC(scOut, sum.getMask(), sum.getArray(), tSysSum/nR, outScanID, 378 355 timeSum/nR, intSum, sourceNameStart, freqIDStart); 379 356 pTabOut->putSDContainer(scOut); 380 357 pTabOut->resetCursor(); 381 // 358 382 359 return CountedPtr<SDMemTable>(pTabOut); 383 360 } … … 385 362 386 363 387 CountedPtr<SDMemTable> SDMath::binaryOperate(const CountedPtr<SDMemTable>& 388 left, 389 const CountedPtr<SDMemTable>& 390 right, 391 const String& op, Bool preserve, 392 Bool doTSys) const 393 { 394 395 // Check operator 396 364 CountedPtr<SDMemTable> 365 SDMath::binaryOperate(const CountedPtr<SDMemTable>& left, 366 const CountedPtr<SDMemTable>& right, 367 const String& op, Bool preserve, Bool doTSys) const 368 { 369 370 // Check operator 397 371 String op2(op); 398 372 op2.upcase(); … … 413 387 } 414 388 415 // Check rows 416 389 // Check rows 417 390 const uInt nRowLeft = left->nRow(); 418 391 const uInt nRowRight = right->nRow(); 419 Bool ok = (nRowRight==1 &&nRowLeft>0) ||420 (nRowRight>=1 &&nRowLeft==nRowRight);392 Bool ok = (nRowRight==1 && nRowLeft>0) || 393 (nRowRight>=1 && nRowLeft==nRowRight); 421 394 if (!ok) { 422 395 throw (AipsError("The right Scan Table can have one row or the same number of rows as the left Scan Table")); 423 396 } 424 397 425 // Input Tables 426 398 // Input Tables 427 399 const Table& tLeft = left->table(); 428 400 const Table& tRight = right->table(); 429 401 430 // TSys columns 431 402 // TSys columns 432 403 ROArrayColumn<Float> tSysLeftCol, tSysRightCol; 433 404 if (doTSys) { 434 tSysLeftCol.attach(tLeft, "TSYS"); 435 tSysRightCol.attach(tRight, "TSYS"); 436 } 437 438 // First row for right 439 405 tSysLeftCol.attach(tLeft, "TSYS"); 406 tSysRightCol.attach(tRight, "TSYS"); 407 } 408 409 // First row for right 440 410 Array<Float> tSysLeftArr, tSysRightArr; 441 411 if (doTSys) tSysRightCol.get(0, tSysRightArr); 442 MaskedArray<Float>* pMRight = new MaskedArray<Float>(right->rowAsMaskedArray(0)); 412 MaskedArray<Float>* pMRight = 413 new MaskedArray<Float>(right->rowAsMaskedArray(0)); 414 443 415 IPosition shpRight = pMRight->shape(); 444 416 445 // Output Table cloned from left 446 417 // Output Table cloned from left 447 418 SDMemTable* pTabOut = new SDMemTable(*left, True); 448 419 pTabOut->appendToHistoryTable(right->getHistoryTable()); 449 // Loop over rows 450 420 421 // Loop over rows 451 422 for (uInt i=0; i<nRowLeft; i++) { 452 453 // Get data 454 455 MaskedArray<Float> mLeft(left->rowAsMaskedArray(i)); 456 IPosition shpLeft = mLeft.shape(); 457 if (doTSys) tSysLeftCol.get(i, tSysLeftArr); 458 // 459 if (nRowRight>1) { 460 delete pMRight; 461 pMRight = new MaskedArray<Float>(right->rowAsMaskedArray(i)); 462 shpRight = pMRight->shape(); 463 if (doTSys) tSysRightCol.get(i, tSysRightArr); 464 } 465 // 466 if (!shpRight.isEqual(shpLeft)) { 467 throw(AipsError("left and right scan tables are not conformant")); 468 } 469 if (doTSys) { 470 if (!tSysRightArr.shape().isEqual(tSysRightArr.shape())) { 471 throw(AipsError("left and right Tsys data are not conformant")); 472 } 473 if (!shpRight.isEqual(tSysRightArr.shape())) { 474 throw(AipsError("left and right scan tables are not conformant")); 475 } 476 } 477 478 // Make container 479 423 424 // Get data 425 MaskedArray<Float> mLeft(left->rowAsMaskedArray(i)); 426 IPosition shpLeft = mLeft.shape(); 427 if (doTSys) tSysLeftCol.get(i, tSysLeftArr); 428 429 if (nRowRight>1) { 430 delete pMRight; 431 pMRight = new MaskedArray<Float>(right->rowAsMaskedArray(i)); 432 shpRight = pMRight->shape(); 433 if (doTSys) tSysRightCol.get(i, tSysRightArr); 434 } 435 436 if (!shpRight.isEqual(shpLeft)) { 437 delete pTabOut; 438 delete pMRight; 439 throw(AipsError("left and right scan tables are not conformant")); 440 } 441 if (doTSys) { 442 if (!tSysRightArr.shape().isEqual(tSysRightArr.shape())) { 443 delete pTabOut; 444 delete pMRight; 445 throw(AipsError("left and right Tsys data are not conformant")); 446 } 447 if (!shpRight.isEqual(tSysRightArr.shape())) { 448 delete pTabOut; 449 delete pMRight; 450 throw(AipsError("left and right scan tables are not conformant")); 451 } 452 } 453 454 // Make container 480 455 SDContainer sc = left->getSDContainer(i); 481 456 482 // Operate on data and TSys 483 457 // Operate on data and TSys 484 458 if (what==0) { 485 459 MaskedArray<Float> tmp = mLeft + *pMRight; … … 511 485 } 512 486 513 // Put new row in output Table 514 487 // Put new row in output Table 515 488 pTabOut->putSDContainer(sc); 516 489 } 517 490 if (pMRight) delete pMRight; 518 491 pTabOut->resetCursor(); 519 492 520 493 return CountedPtr<SDMemTable>(pTabOut); 521 494 } 522 523 495 524 496 … … 828 800 shapeOut(asap::PolAxis) = 1; // Average all polarizations 829 801 if (shapeIn(asap::PolAxis)==1) { 830 throw(AipsError("The input has only one polarisation")); 802 delete pTabOut; 803 throw(AipsError("The input has only one polarisation")); 831 804 } 832 805 // … … 876 849 ReadOnlyArrayIterator<Bool> itMaskPlane(barr, axes); 877 850 ReadOnlyArrayIterator<Float>* pItTsysPlane = 0; 878 if (wtType==TSYS) pItTsysPlane = new ReadOnlyArrayIterator<Float>(tSys, axes); 851 if (wtType==TSYS) 852 pItTsysPlane = new ReadOnlyArrayIterator<Float>(tSys, axes); 879 853 880 854 // Accumulations … … 900 874 ReadOnlyVectorIterator<Float>* pItTsysVec = 0; 901 875 if (wtType==TSYS) { 902 pItTsysVec = new ReadOnlyVectorIterator<Float>(pItTsysPlane->array(), 1); 876 pItTsysVec = 877 new ReadOnlyVectorIterator<Float>(pItTsysPlane->array(), 1); 903 878 } 904 879 // … … 908 883 909 884 if (useMask) { 910 const MaskedArray<Float> spec(itDataVec.vector(),mask&&itMaskVec.vector()); 885 const MaskedArray<Float> spec(itDataVec.vector(), 886 mask&&itMaskVec.vector()); 911 887 if (wtType==VAR) { 912 888 fac = 1.0 / variance(spec); … … 916 892 } 917 893 } else { 918 const MaskedArray<Float> spec(itDataVec.vector(),itMaskVec.vector()); 894 const MaskedArray<Float> spec(itDataVec.vector(), 895 itMaskVec.vector()); 919 896 if (wtType==VAR) { 920 897 fac = 1.0 / variance(spec); … … 927 904 // Normalize spectrum (without OTF mask) and accumulate 928 905 929 const MaskedArray<Float> spec(fac*itDataVec.vector(), itMaskVec.vector()); 906 const MaskedArray<Float> spec(fac*itDataVec.vector(), 907 itMaskVec.vector()); 930 908 vecSum += spec; 931 909 norm += fac; … … 1001 979 { 1002 980 1003 // Number of channels 1004 981 // Number of channels 1005 982 const uInt nChan = in.nChan(); 1006 983 1007 // Generate Kernel 1008 984 // Generate Kernel 1009 985 VectorKernel::KernelTypes type = VectorKernel::toKernelType(kernelType); 1010 986 Vector<Float> kernel = VectorKernel::make(type, width, nChan, True, False); 1011 987 1012 // Generate Convolver 1013 988 // Generate Convolver 1014 989 IPosition shape(1,nChan); 1015 990 Convolver<Float> conv(kernel, shape); 1016 991 1017 // New Table 1018 992 // New Table 1019 993 SDMemTable* pTabOut = new SDMemTable(in,True); 1020 994 1021 // Output Vectors 1022 995 // Output Vectors 1023 996 Vector<Float> valuesOut(nChan); 1024 997 Vector<Bool> maskOut(nChan); 1025 998 1026 // Get data slice bounds 1027 999 // Get data slice bounds 1028 1000 IPosition start, end; 1029 1001 setCursorSlice (start, end, doAll, in); 1030 1002 1031 // Loop over rows in Table 1032 1003 // Loop over rows in Table 1033 1004 for (uInt ri=0; ri < in.nRow(); ++ri) { 1034 1005 1035 // Get slice of data 1036 1006 // Get slice of data 1037 1007 MaskedArray<Float> dataIn = in.rowAsMaskedArray(ri); 1038 1008 1039 // Deconstruct and get slices which reference these arrays 1040 1009 // Deconstruct and get slices which reference these arrays 1041 1010 Array<Float> valuesIn = dataIn.getArray(); 1042 1011 Array<Bool> maskIn = dataIn.getMask(); 1043 // 1012 1044 1013 Array<Float> valuesIn2 = valuesIn(start,end); // ref to valuesIn 1045 1014 Array<Bool> maskIn2 = maskIn(start,end); 1046 1015 1047 // Iterate through by spectra 1048 1016 // Iterate through by spectra 1049 1017 VectorIterator<Float> itValues(valuesIn2, asap::ChanAxis); 1050 1018 VectorIterator<Bool> itMask(maskIn2, asap::ChanAxis); 1051 1019 while (!itValues.pastEnd()) { 1052 1053 // Smooth 1054 1055 if (kernelType==VectorKernel::HANNING) { 1056 mathutil::hanning(valuesOut, maskOut, itValues.vector(), itMask.vector()); 1057 itMask.vector() = maskOut; 1058 } else { 1059 mathutil::replaceMaskByZero(itValues.vector(), itMask.vector()); 1060 conv.linearConv(valuesOut, itValues.vector()); 1061 } 1062 // 1063 itValues.vector() = valuesOut; 1064 // 1065 itValues.next(); 1066 itMask.next(); 1020 1021 // Smooth 1022 if (kernelType==VectorKernel::HANNING) { 1023 mathutil::hanning(valuesOut, maskOut, itValues.vector(), 1024 itMask.vector()); 1025 itMask.vector() = maskOut; 1026 } else { 1027 mathutil::replaceMaskByZero(itValues.vector(), itMask.vector()); 1028 conv.linearConv(valuesOut, itValues.vector()); 1029 } 1030 1031 itValues.vector() = valuesOut; 1032 itValues.next(); 1033 itMask.next(); 1067 1034 } 1068 1069 // Create and put back 1070 1035 1036 // Create and put back 1071 1037 SDContainer sc = in.getSDContainer(ri); 1072 1038 putDataInSDC(sc, valuesIn, maskIn); 1073 // 1039 1074 1040 pTabOut->putSDContainer(sc); 1075 1041 } 1076 // 1042 1077 1043 return pTabOut; 1078 1044 } … … 1091 1057 SDMemTable* pTabOut = new SDMemTable(in, True); 1092 1058 1093 // Find out how to convert values into Jy and K (e.g. units might be mJy or mK) 1094 // Also automatically find out what we are converting to according to the 1095 // flux unit 1096 1059 // Find out how to convert values into Jy and K (e.g. units might be 1060 // mJy or mK) Also automatically find out what we are converting to 1061 // according to the flux unit 1097 1062 Unit fluxUnit(sh.fluxunit); 1098 1063 Unit K(String("K")); 1099 1064 Unit JY(String("Jy")); 1100 // 1065 1101 1066 Bool toKelvin = True; 1102 1067 Double cFac = 1.0; 1103 1068 if (fluxUnit==JY) { 1104 1105 // 1106 1107 1108 1109 // 1110 1111 1069 cout << "Converting to K" << endl; 1070 1071 Quantum<Double> t(1.0,fluxUnit); 1072 Quantum<Double> t2 = t.get(JY); 1073 cFac = (t2 / t).getValue(); // value to Jy 1074 1075 toKelvin = True; 1076 sh.fluxunit = "K"; 1112 1077 } else if (fluxUnit==K) { 1113 1114 // 1115 1116 1117 1118 // 1119 1120 1078 cout << "Converting to Jy" << endl; 1079 1080 Quantum<Double> t(1.0,fluxUnit); 1081 Quantum<Double> t2 = t.get(K); 1082 cFac = (t2 / t).getValue(); // value to K 1083 1084 toKelvin = False; 1085 sh.fluxunit = "Jy"; 1121 1086 } else { 1122 1087 throw(AipsError("Unrecognized brightness units in Table - must be consistent with Jy or K")); 1123 1088 } 1124 1089 pTabOut->putSDHeader(sh); 1125 1126 // Make sure input values are converted to either Jy or K first... 1127 1090 1091 // Make sure input values are converted to either Jy or K first... 1128 1092 Float factor = cFac; 1129 1093 1130 // Select method 1131 1094 // Select method 1132 1095 if (JyPerK>0.0) { 1133 factor *= JyPerK; 1134 if (toKelvin) factor = 1.0 / JyPerK; 1135 // 1136 cout << "Jy/K = " << JyPerK << endl; 1137 Vector<Float> factors(in.nRow(), factor); 1138 scaleByVector(pTabOut, in, doAll, factors, False); 1096 factor *= JyPerK; 1097 if (toKelvin) factor = 1.0 / JyPerK; 1098 cout << "Jy/K = " << JyPerK << endl; 1099 Vector<Float> factors(in.nRow(), factor); 1100 scaleByVector(pTabOut, in, doAll, factors, False); 1139 1101 } else if (etaAp>0.0) { 1140 1141 1142 1143 1144 1145 1146 1147 1148 1149 1150 // 1151 1152 1102 Bool throwIt = True; 1103 Instrument inst = SDAttr::convertInstrument (sh.antennaname, throwIt); 1104 SDAttr sda; 1105 if (D < 0) D = sda.diameter(inst); 1106 Float JyPerK = SDAttr::findJyPerK (etaAp,D); 1107 cout << "Jy/K = " << JyPerK << endl; 1108 factor *= JyPerK; 1109 if (toKelvin) { 1110 factor = 1.0 / factor; 1111 } 1112 1113 Vector<Float> factors(in.nRow(), factor); 1114 scaleByVector(pTabOut, in, doAll, factors, False); 1153 1115 } else { 1154 1155 // OK now we must deal with automatic look up of values.1156 // We must also deal with the fact that the factors need1157 // to be computed per IF and may be different and may1158 // change per integration.1159 1160 1161 1162 } 1163 // 1116 1117 // OK now we must deal with automatic look up of values. 1118 // We must also deal with the fact that the factors need 1119 // to be computed per IF and may be different and may 1120 // change per integration. 1121 1122 cout << "Looking up conversion factors" << endl; 1123 convertBrightnessUnits (pTabOut, in, toKelvin, cFac, doAll); 1124 } 1125 1164 1126 return pTabOut; 1165 1127 } 1166 1167 1168 1169 1128 1170 1129 … … 1175 1134 { 1176 1135 1177 // Get header and clone output table 1178 1136 // Get header and clone output table 1179 1137 SDHeader sh = in.getSDHeader(); 1180 1138 SDMemTable* pTabOut = new SDMemTable(in, True); 1181 1139 1182 // Get elevation data from SDMemTable and convert to degrees 1183 1140 // Get elevation data from SDMemTable and convert to degrees 1184 1141 const Table& tab = in.table(); 1185 1142 ROScalarColumn<Float> elev(tab, "ELEVATION"); 1186 1143 Vector<Float> x = elev.getColumn(); 1187 1144 x *= Float(180 / C::pi); // Degrees 1188 // 1145 1189 1146 const uInt nC = coeffs.nelements(); 1190 1147 if (fileName.length()>0 && nC>0) { 1191 throw(AipsError("You must choose either polynomial coefficients or an ascii file, not both")); 1192 } 1193 1194 // Correct 1195 1148 throw(AipsError("You must choose either polynomial coefficients or an ascii file, not both")); 1149 } 1150 1151 // Correct 1196 1152 if (nC>0 || fileName.length()==0) { 1197 1198 // Find instrument 1199 1153 // Find instrument 1200 1154 Bool throwIt = True; 1201 1155 Instrument inst = SDAttr::convertInstrument (sh.antennaname, throwIt); 1202 1156 1203 // Set polynomial 1204 1157 // Set polynomial 1205 1158 Polynomial<Float>* pPoly = 0; 1206 1159 Vector<Float> coeff; 1207 1160 String msg; 1208 1161 if (nC>0) { 1209 1210 1211 1162 pPoly = new Polynomial<Float>(nC); 1163 coeff = coeffs; 1164 msg = String("user"); 1212 1165 } else { 1213 1214 1215 1216 1166 SDAttr sdAttr; 1167 coeff = sdAttr.gainElevationPoly(inst); 1168 pPoly = new Polynomial<Float>(3); 1169 msg = String("built in"); 1217 1170 } 1218 // 1171 1219 1172 if (coeff.nelements()>0) { 1220 1173 pPoly->setCoefficients(coeff); 1221 1174 } else { 1222 throw(AipsError("There is no known gain-elevation polynomial known for this instrument")); 1175 delete pPoly; 1176 throw(AipsError("There is no known gain-elevation polynomial known for this instrument")); 1223 1177 } 1224 //1225 1178 cout << "Making polynomial correction with " << msg << " coefficients" << endl; 1226 1179 const uInt nRow = in.nRow(); 1227 1180 Vector<Float> factor(nRow); 1228 1181 for (uInt i=0; i<nRow; i++) { 1229 1182 factor[i] = 1.0 / (*pPoly)(x[i]); 1230 1183 } 1231 1184 delete pPoly; 1232 //1233 1185 scaleByVector (pTabOut, in, doAll, factor, True); 1186 1234 1187 } else { 1235 1236 // Indicate which columns to read from ascii file 1237 1238 String col0("ELEVATION"); 1239 String col1("FACTOR"); 1240 1241 // Read and correct 1242 1243 cout << "Making correction from ascii Table" << endl; 1244 scaleFromAsciiTable (pTabOut, in, fileName, col0, col1, 1245 methodStr, doAll, x, True); 1246 } 1247 // 1248 return pTabOut; 1249 } 1250 1188 1189 // Indicate which columns to read from ascii file 1190 String col0("ELEVATION"); 1191 String col1("FACTOR"); 1192 1193 // Read and correct 1194 1195 cout << "Making correction from ascii Table" << endl; 1196 scaleFromAsciiTable (pTabOut, in, fileName, col0, col1, 1197 methodStr, doAll, x, True); 1198 } 1199 1200 return pTabOut; 1201 } 1251 1202 1252 1203 … … 1254 1205 { 1255 1206 1256 // Get header and clone output table1207 // Get header and clone output table 1257 1208 1258 1209 SDHeader sh = in.getSDHeader(); … … 1278 1229 1279 1230 scaleByVector (pTabOut, in, doAll, factor, True); 1280 // 1231 1281 1232 return pTabOut; 1282 1233 } … … 1289 1240 // 1290 1241 { 1291 1292 1293 1294 // 1242 if (in.nPol() != 4) { 1243 throw(AipsError("You must have 4 polarizations to run this function")); 1244 } 1245 1295 1246 SDHeader sh = in.getSDHeader(); 1296 1247 Instrument inst = SDAttr::convertInstrument (sh.antennaname, False); … … 1324 1275 IPosition shape = data.shape(); 1325 1276 1326 // Get polarization slice references 1327 1277 // Get polarization slice references 1328 1278 Array<Float> C3 = data(start3,end3); 1329 1279 Array<Float> C4 = data(start4,end4); 1330 1280 1331 // Rotate 1332 1281 // Rotate 1333 1282 SDPolUtil::rotatePhase(C3, C4, value); 1334 1283 1335 // Put 1336 1284 // Put 1337 1285 specCol.put(i,data); 1338 1286 } 1339 1287 } 1340 1341 1288 1342 1289 … … 1523 1470 1524 1471 1525 SDMemTable* SDMath::frequencyAlign 1472 SDMemTable* SDMath::frequencyAlign(const SDMemTable& in, 1526 1473 MFrequency::Types freqSystem, 1527 1474 const String& refTime, … … 1551 1498 Matrix<uInt> ddIdx; 1552 1499 SDDataDesc dDesc; 1553 generateDataDescTable (ddIdx, dDesc, nIF, in, tabIn, srcCol, fqIDCol, perFreqID); 1500 generateDataDescTable(ddIdx, dDesc, nIF, in, tabIn, srcCol, 1501 fqIDCol, perFreqID); 1554 1502 1555 1503 // Get reference Epoch to time of first row or given String 1556 1504 1557 1505 Unit DAY(String("d")); 1558 1506 MEpoch::Ref epochRef(in.getTimeReference()); 1559 1507 MEpoch refEpoch; 1560 1508 if (refTime.length()>0) { 1561 1509 refEpoch = epochFromString(refTime, in.getTimeReference()); 1562 1510 } else { 1563 1511 refEpoch = in.getEpoch(0); 1564 1512 } 1565 1513 cout << "Aligning at reference Epoch " << formatEpoch(refEpoch) 1566 1514 << " in frame " << MFrequency::showType(freqSystem) << endl; 1567 1515 1568 // Get Reference Position1569 1516 // Get Reference Position 1517 1570 1518 MPosition refPos = in.getAntennaPosition(); 1571 1572 // Create FrequencyAligner Block. One FA for each possible 1573 // source/freqID (perFreqID=True) or source/IF (perFreqID=False) combination 1574 1519 1520 // Create FrequencyAligner Block. One FA for each possible 1521 // source/freqID (perFreqID=True) or source/IF (perFreqID=False) 1522 // combination 1523 1575 1524 PtrBlock<FrequencyAligner<Float>* > a(dDesc.length()); 1576 1525 generateFrequencyAligners (a, dDesc, in, nChan, freqSystem, refPos, 1577 1526 refEpoch, perFreqID); 1578 1579 // Generate and fill output Frequency Table. WHen perFreqID=True, there is one output FreqID 1580 // for each entry in the SDDataDesc table. However, in perFreqID=False mode, there may be 1581 // some degeneracy, so we need a little translation map 1582 1527 1528 // Generate and fill output Frequency Table. WHen perFreqID=True, 1529 // there is one output FreqID for each entry in the SDDataDesc 1530 // table. However, in perFreqID=False mode, there may be some 1531 // degeneracy, so we need a little translation map 1532 1583 1533 SDFrequencyTable freqTabOut = in.getSDFreqTable(); 1584 1534 freqTabOut.setLength(0); … … 1586 1536 units = String("Hz"); 1587 1537 Bool linear=True; 1588 //1538 // 1589 1539 Vector<uInt> ddFQTrans(dDesc.length(),0); 1590 1540 for (uInt i=0; i<dDesc.length(); i++) { 1591 1541 1592 // Get Aligned SC in Hz 1593 1594 SpectralCoordinate sC = a[i]->alignedSpectralCoordinate(linear); 1595 sC.setWorldAxisUnits(units); 1596 1597 // Add FreqID 1598 1599 uInt idx = freqTabOut.addFrequency(sC.referencePixel()[0], 1600 sC.referenceValue()[0], 1601 sC.increment()[0]); 1602 ddFQTrans(i) = idx; // output FreqID = ddFQTrans(ddIdx) 1542 // Get Aligned SC in Hz 1543 1544 SpectralCoordinate sC = a[i]->alignedSpectralCoordinate(linear); 1545 sC.setWorldAxisUnits(units); 1546 1547 // Add FreqID 1548 1549 uInt idx = freqTabOut.addFrequency(sC.referencePixel()[0], 1550 sC.referenceValue()[0], 1551 sC.increment()[0]); 1552 // output FreqID = ddFQTrans(ddIdx) 1553 ddFQTrans(i) = idx; 1603 1554 } 1604 1605 // Interpolation method1606 1555 1556 // Interpolation method 1557 1607 1558 InterpolateArray1D<Double,Float>::InterpolationMethod interp; 1608 1559 convertInterpString(interp, methodStr); 1609 1610 // New output Table1611 1560 1561 // New output Table 1562 1612 1563 cout << "Create output table" << endl; 1613 1564 SDMemTable* pTabOut = new SDMemTable(in,True); 1614 1565 pTabOut->putSDFreqTable(freqTabOut); 1615 1616 // Loop over rows in Table1617 1566 1567 // Loop over rows in Table 1568 1618 1569 Bool extrapolate=False; 1619 1570 const IPosition polChanAxes(2, asap::PolAxis, asap::ChanAxis); … … 1626 1577 uInt ifIdx, faIdx; 1627 1578 Vector<Double> xIn; 1628 //1579 // 1629 1580 for (uInt iRow=0; iRow<nRows; ++iRow) { 1630 1631 1632 1633 1634 // Get EPoch1635 1581 if (iRow%10==0) { 1582 cout << "Processing row " << iRow << endl; 1583 } 1584 1585 // Get EPoch 1586 1636 1587 Quantum<Double> tQ2(times[iRow],DAY); 1637 1588 MVEpoch mv2(tQ2); 1638 1589 MEpoch epoch(mv2, epochRef); 1639 1640 // Get copy of data1641 1590 1591 // Get copy of data 1592 1642 1593 const MaskedArray<Float>& mArrIn(in.rowAsMaskedArray(iRow)); 1643 1594 Array<Float> values = mArrIn.getArray(); 1644 1595 Array<Bool> mask = mArrIn.getMask(); 1645 1646 // For each row, the Frequency abcissa will be the same regardless 1647 // of polarization. For all other axes (IF and BEAM) the abcissa 1648 // will change. So we iterate through the data by pol-chan planes 1649 // to mimimize the work. Probably won't work for multiple beams 1650 //at this point.1651 1596 1597 // For each row, the Frequency abcissa will be the same 1598 // regardless of polarization. For all other axes (IF and BEAM) 1599 // the abcissa will change. So we iterate through the data by 1600 // pol-chan planes to mimimize the work. Probably won't work for 1601 // multiple beams at this point. 1602 1652 1603 ArrayIterator<Float> itValuesPlane(values, polChanAxes); 1653 1604 ArrayIterator<Bool> itMaskPlane(mask, polChanAxes); 1654 1605 while (!itValuesPlane.pastEnd()) { 1655 1606 1656 // Find the IF index and then the FA PtrBlock index1657 1658 1659 1660 1661 1662 // Generate abcissa for perIF. Could cache this in a Matrix 1663 // on a per scan basis.Pretty expensive doing it for every row.1664 1665 1666 1667 1668 1607 // Find the IF index and then the FA PtrBlock index 1608 1609 const IPosition& pos = itValuesPlane.pos(); 1610 ifIdx = pos(asap::IFAxis); 1611 faIdx = ddIdx(iRow,ifIdx); 1612 1613 // Generate abcissa for perIF. Could cache this in a Matrix on 1614 // a per scan basis. Pretty expensive doing it for every row. 1615 1616 if (!perFreqID) { 1617 xIn.resize(nChan); 1618 uInt fqID = dDesc.secID(ddIdx(iRow,ifIdx)); 1619 SpectralCoordinate sC = in.getSpectralCoordinate(fqID); 1669 1620 Double w; 1670 1621 for (uInt i=0; i<nChan; i++) { 1671 1622 sC.toWorld(w,Double(i)); 1672 1623 xIn[i] = w; 1673 1624 } 1674 1675 // 1676 1677 1678 1679 // Iterate through the plane by vector and align1680 1625 } 1626 1627 VectorIterator<Float> itValuesVec(itValuesPlane.array(), 1); 1628 VectorIterator<Bool> itMaskVec(itMaskPlane.array(), 1); 1629 1630 // Iterate through the plane by vector and align 1631 1681 1632 first = True; 1682 1633 useCachedAbcissa=False; 1683 1634 while (!itValuesVec.pastEnd()) { 1684 1685 1686 1687 1688 1689 1690 1691 1692 1693 //1694 1695 1696 //1697 1698 1699 //1700 1701 1702 1703 1635 if (perFreqID) { 1636 ok = a[faIdx]->align (yOut, maskOut, itValuesVec.vector(), 1637 itMaskVec.vector(), epoch, useCachedAbcissa, 1638 interp, extrapolate); 1639 } else { 1640 ok = a[faIdx]->align (yOut, maskOut, xIn, itValuesVec.vector(), 1641 itMaskVec.vector(), epoch, useCachedAbcissa, 1642 interp, extrapolate); 1643 } 1644 // 1645 itValuesVec.vector() = yOut; 1646 itMaskVec.vector() = maskOut; 1647 // 1648 itValuesVec.next(); 1649 itMaskVec.next(); 1650 // 1651 if (first) { 1652 useCachedAbcissa = True; 1653 first = False; 1654 } 1704 1655 } 1705 //1706 1707 1656 // 1657 itValuesPlane.next(); 1658 itMaskPlane.next(); 1708 1659 } 1709 1710 // Create SDContainer and put back1711 1712 SDContainer sc = in.getSDContainer(iRow);1713 putDataInSDC(sc, values, mask);1714 1715 // Set output FreqIDs1716 1717 for (uInt i=0; i<nIF; i++) {1660 1661 // Create SDContainer and put back 1662 1663 SDContainer sc = in.getSDContainer(iRow); 1664 putDataInSDC(sc, values, mask); 1665 1666 // Set output FreqIDs 1667 1668 for (uInt i=0; i<nIF; i++) { 1718 1669 uInt idx = ddIdx(iRow,i); // Index into SDDataDesc table 1719 1670 freqID(i) = ddFQTrans(idx); // FreqID in output FQ table 1720 }1721 sc.putFreqMap(freqID);1722 //1723 pTabOut->putSDContainer(sc);1671 } 1672 sc.putFreqMap(freqID); 1673 // 1674 pTabOut->putSDContainer(sc); 1724 1675 } 1725 1726 // Now we must set the base and extra frames to the 1727 // input frame 1728 1676 1677 // Now we must set the base and extra frames to the input frame 1729 1678 std::vector<string> info = pTabOut->getCoordInfo(); 1730 1679 info[1] = MFrequency::showType(freqSystem); // Conversion frame … … 1732 1681 pTabOut->setCoordInfo(info); 1733 1682 1734 // Clean up PointerBlock 1735 1683 // Clean up PointerBlock 1736 1684 for (uInt i=0; i<a.nelements(); i++) delete a[i]; 1737 // 1685 1738 1686 return pTabOut; 1739 1687 } 1740 1688 1689 1690 SDMemTable* SDMath::frequencySwitch(const SDMemTable& in) const 1691 { 1692 if (in.nIF() != 2) { 1693 throw(AipsError("nIF != 2 ")); 1694 } 1695 Bool clear = True; 1696 SDMemTable* pTabOut = new SDMemTable(in, clear); 1697 const Table& tabIn = in.table(); 1698 1699 // Shape of input and output data 1700 const IPosition& shapeIn = in.rowAsMaskedArray(0).shape(); 1701 1702 const uInt nRows = in.nRow(); 1703 AlwaysAssert(asap::nAxes==4,AipsError); 1704 1705 ROArrayColumn<Float> tSysCol; 1706 Array<Float> tsys; 1707 tSysCol.attach(tabIn,"TSYS"); 1708 1709 for (uInt iRow=0; iRow<nRows; iRow++) { 1710 // Get data for this row 1711 MaskedArray<Float> marr(in.rowAsMaskedArray(iRow)); 1712 tSysCol.get(iRow, tsys); 1713 1714 // whole Array for IF 0 1715 IPosition start(asap::nAxes,0); 1716 IPosition end = shapeIn-1; 1717 end(asap::IFAxis) = 0; 1718 1719 MaskedArray<Float> on = marr(start,end); 1720 Array<Float> ton = tsys(start,end); 1721 // Make a copy as "src" is a refrence which is manipulated. 1722 // oncopy is needed for the inverse quotient 1723 MaskedArray<Float> oncopy = on.copy(); 1724 1725 // whole Array for IF 1 1726 start(asap::IFAxis) = 1; 1727 end(asap::IFAxis) = 1; 1728 1729 MaskedArray<Float> off = marr(start,end); 1730 Array<Float> toff = tsys(start,end); 1731 1732 on /= off; on -= 1.0f; 1733 on *= ton; 1734 off /= oncopy; off -= 1.0f; 1735 off *= toff; 1736 1737 SDContainer sc = in.getSDContainer(iRow); 1738 putDataInSDC(sc, marr.getArray(), marr.getMask()); 1739 pTabOut->putSDContainer(sc); 1740 } 1741 return pTabOut; 1742 } 1741 1743 1742 1744 void SDMath::fillSDC(SDContainer& sc, … … 2163 2165 2164 2166 2165 MEpoch SDMath::epochFromString 2167 MEpoch SDMath::epochFromString(const String& str, MEpoch::Types timeRef) const 2166 2168 { 2167 2169 Quantum<Double> qt; … … 2184 2186 2185 2187 2186 void SDMath::generateFrequencyAligners 2187 2188 2189 2190 2191 2192 2188 void SDMath::generateFrequencyAligners(PtrBlock<FrequencyAligner<Float>* >& a, 2189 const SDDataDesc& dDesc, 2190 const SDMemTable& in, uInt nChan, 2191 MFrequency::Types system, 2192 const MPosition& refPos, 2193 const MEpoch& refEpoch, 2194 Bool perFreqID) const 2193 2195 { 2194 2196 for (uInt i=0; i<dDesc.length(); i++) { … … 2214 2216 } 2215 2217 2216 Vector<uInt> SDMath::getRowRange 2218 Vector<uInt> SDMath::getRowRange(const SDMemTable& in) const 2217 2219 { 2218 2220 Vector<uInt> range(2); … … 2223 2225 2224 2226 2225 Bool SDMath::rowInRange 2227 Bool SDMath::rowInRange(uInt i, const Vector<uInt>& range) const 2226 2228 { 2227 2229 return (i>=range[0] && i<=range[1]); 2228 2230 } 2231 -
branches/Release-2-fixes/src/SDMath.h
r653 r699 63 63 64 64 // Copy Constructor (copy semantics) 65 SDMath 65 SDMath(const SDMath& other); 66 66 67 67 // Assignment (copy semantics) … … 72 72 73 73 // Binary Table operators. op=ADD, SUB, MUL, DIV, QUOTIENT 74 casa::CountedPtr<SDMemTable> binaryOperate 74 casa::CountedPtr<SDMemTable> binaryOperate(const casa::CountedPtr<SDMemTable>& left, 75 75 const casa::CountedPtr<SDMemTable>& right, 76 76 const casa::String& op, casa::Bool preserve, … … 99 99 100 100 // Smooth 101 SDMemTable* smooth 101 SDMemTable* smooth(const SDMemTable& in, const casa::String& kernel, 102 102 casa::Float width, casa::Bool doAll) const; 103 103 104 104 // Flux conversion between Jansky and Kelvin 105 SDMemTable* convertFlux 105 SDMemTable* convertFlux(const SDMemTable& in, casa::Float D, casa::Float etaAp, 106 106 casa::Float JyPerK, casa::Bool doAll) const; 107 107 108 108 // Gain-elevation correction 109 SDMemTable* gainElevation 109 SDMemTable* gainElevation(const SDMemTable& in, const casa::Vector<casa::Float>& coeffs, 110 110 const casa::String& fileName, 111 111 const casa::String& method, casa::Bool doAll) const; 112 112 113 113 // Frequency Alignment 114 SDMemTable* frequencyAlignment 114 SDMemTable* frequencyAlignment(const SDMemTable& in, const casa::String& refTime, 115 115 const casa::String& method, casa::Bool perFreqID) const; 116 116 117 117 // Opacity correction 118 SDMemTable* opacity 118 SDMemTable* opacity(const SDMemTable& in, casa::Float tau, casa::Bool doAll) const; 119 119 120 120 // Simple unary mathematical operations. what=0 (mul) or 1 (add) … … 126 126 const casa::String& wtStr) const; 127 127 128 SDMemTable* frequencySwitch(const SDMemTable& in) const; 129 130 128 131 // Rotate XY phase. Value in degrees. 129 void rotateXYPhase 132 void rotateXYPhase(SDMemTable& in, casa::Float value, casa::Bool doAll); 130 133 131 134 // Rotate Q & U by operating on the raw correlations.Value in degrees. 132 void rotateLinPolPhase 135 void rotateLinPolPhase(SDMemTable& in, casa::Float value, casa::Bool doAll); 133 136 134 137 … … 158 161 159 162 // Work out conversion factor for converting Jy<->K per IF per row and apply 160 void convertBrightnessUnits 163 void convertBrightnessUnits(SDMemTable* pTabOut, const SDMemTable& in, 161 164 casa::Bool toKelvin, casa::Float sFac, casa::Bool doAll) const; 162 165 163 166 // Convert weight string to enum value 164 167 165 void convertWeightString 168 void convertWeightString(WeightType& wt, const casa::String& weightStr, 166 169 casa::Bool listType) const; 167 170 … … 185 188 186 189 // Scale data and optionally TSys by values in a Vector 187 void scaleByVector 190 void scaleByVector(SDMemTable* pTabOut, const SDMemTable& in, 188 191 casa::Bool doAll, const casa::Vector<casa::Float>& factor, 189 192 casa::Bool doTSys) const; 190 193 191 194 // Convert time String to Epoch 192 casa::MEpoch epochFromString 195 casa::MEpoch epochFromString(const casa::String& str, casa::MEpoch::Types timeRef) const; 193 196 194 197 // Function to fill Scan Container when averaging in time 195 198 196 void fillSDC 199 void fillSDC(SDContainer& sc, const casa::Array<casa::Bool>& mask, 197 200 const casa::Array<casa::Float>& data, 198 201 const casa::Array<casa::Float>& tSys, … … 205 208 206 209 // Align in Frequency 207 SDMemTable* frequencyAlign 210 SDMemTable* frequencyAlign(const SDMemTable& in, 208 211 casa::MFrequency::Types system, 209 212 const casa::String& timeRef, … … 212 215 213 216 // Generate frequency aligners 214 void generateFrequencyAligners 217 void generateFrequencyAligners(casa::PtrBlock<casa::FrequencyAligner<casa::Float>* >& a, 215 218 const SDDataDesc& dDesc, 216 219 const SDMemTable& in, casa::uInt nChan, … … 220 223 casa::Bool perFreqID) const; 221 224 222 // Generate data description table 223 void generateDataDescTable 225 // Generate data description table(combines source and freqID) 226 void generateDataDescTable(casa::Matrix<casa::uInt>& ddIdx, 224 227 SDDataDesc& dDesc, 225 228 casa::uInt nIF, … … 231 234 232 235 // Get row range from SDMemTable state 233 casa::Vector<casa::uInt> getRowRange 236 casa::Vector<casa::uInt> getRowRange(const SDMemTable& in) const; 234 237 235 238 // Is row in the row range ? 236 casa::Bool rowInRange 239 casa::Bool rowInRange(casa::uInt i, const casa::Vector<casa::uInt>& range) const; 237 240 238 241 // Set slice to cursor or all axes 239 void setCursorSlice 242 void setCursorSlice(casa::IPosition& start, casa::IPosition& end, 240 243 casa::Bool doAll, const SDMemTable& in) const; 241 244 242 245 // Function to normalize data when averaging in time 243 246 244 void normalize 247 void normalize(casa::MaskedArray<casa::Float>& data, 245 248 const casa::Array<casa::Float>& sumSq, 246 249 const casa::Array<casa::Float>& tSysSumSq, … … 250 253 251 254 // Put the data and mask into the SDContainer 252 void putDataInSDC 255 void putDataInSDC(SDContainer& sc, const casa::Array<casa::Float>& data, 253 256 const casa::Array<casa::Bool>& mask) const; 254 257 255 258 // Read ascii file into a Table 256 259 257 casa::Table readAsciiFile 260 casa::Table readAsciiFile(const casa::String& fileName) const; 258 261 }; 259 262
Note:
See TracChangeset
for help on using the changeset viewer.