- Timestamp:
- 11/19/08 14:14:26 (16 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/src/STWriter.cpp
r1443 r1449 39 39 #include <casa/Utilities/Assert.h> 40 40 41 #include <atnf/PKSIO/PKSrecord.h> 41 42 #include <atnf/PKSIO/PKSMS2writer.h> 42 43 #include <atnf/PKSIO/PKSSDwriter.h> … … 155 156 } 156 157 157 Double srcVel; 158 159 String fieldName, srcName, tcalTime; 160 Vector<Float> calFctr, sigma, tcal, tsys; 161 Vector<Double> direction(2), scanRate(2), srcDir(2), srcPM(2); 162 Matrix<Float> spectra; 163 Matrix<uChar> flagtra; 164 Complex xCalFctr; 158 165 159 Int count = 0; 166 Int scanno = 1; 160 PKSrecord pksrec; 161 pksrec.scanNo = 1; 167 162 // use spearate iterators to ensure renumbering of all numbers 168 163 TableIterator scanit(table, "SCANNO"); … … 170 165 Table stable = scanit.table(); 171 166 TableIterator beamit(stable, "BEAMNO"); 172 Int beamno = 1;167 pksrec.beamNo = 1; 173 168 while (!beamit.pastEnd() ) { 174 169 Table btable = beamit.table(); 175 170 TableIterator cycit(btable, "CYCLENO"); 176 171 ROArrayColumn<Double> srateCol(btable, "SCANRATE"); 177 srateCol.get(0, scanRate);172 srateCol.get(0, pksrec.scanRate); 178 173 ROArrayColumn<Double> spmCol(btable, "SRCPROPERMOTION"); 179 spmCol.get(0, srcPM);174 spmCol.get(0, pksrec.srcPM); 180 175 ROArrayColumn <Double> sdirCol(btable, "SRCDIRECTION"); 181 sdirCol.get(0, srcDir);176 sdirCol.get(0, pksrec.srcDir); 182 177 ROScalarColumn<Double> svelCol(btable, "SRCVELOCITY"); 183 svelCol.get(0, srcVel);178 svelCol.get(0, pksrec.srcVel); 184 179 ROScalarColumn<uInt> bCol(btable, "BEAMNO"); 185 beamno = bCol(0)+1;186 Int cycno = 1;180 pksrec.beamNo = bCol(0)+1; 181 pksrec.cycleNo = 1; 187 182 while (!cycit.pastEnd() ) { 188 183 Table ctable = cycit.table(); 189 184 TableIterator ifit(ctable, "IFNO"); 190 185 MDirection::ScalarColumn dirCol(ctable, "DIRECTION"); 191 Vector<Double>direction = dirCol(0).getAngle("rad").getValue();192 Int ifno = 1;186 pksrec.direction = dirCol(0).getAngle("rad").getValue(); 187 pksrec.IFno = 1; 193 188 while (!ifit.pastEnd() ) { 194 189 Table itable = ifit.table(); … … 197 192 const TableRecord& rec = row.get(0); 198 193 ROArrayColumn<Float> specCol(itable, "SPECTRA"); 199 ifno = rec.asuInt("IFNO")+1;194 pksrec.IFno = rec.asuInt("IFNO")+1; 200 195 uInt nchan = specCol(0).nelements(); 201 Double cdelt,crval,crpix, restfreq; 202 Float focusAxi, focusTan, focusRot, 203 temperature, pressure, humidity, windSpeed, windAz; 196 Double crval,crpix; 204 197 Float tmp0,tmp1,tmp2,tmp3,tmp4; 205 Vector<Float> tcalval;206 String stmp0,stmp1, tcalt;207 in->frequencies().getEntry(crpix,crval,cdelt,rec.asuInt("FREQ_ID"));208 in->focus().getEntry( focusAxi, focusTan, focusRot,209 tmp0,tmp1,tmp2,tmp3,tmp4,198 String stmp0,stmp1; 199 in->frequencies().getEntry(crpix,crval, pksrec.freqInc, 200 rec.asuInt("FREQ_ID")); 201 in->focus().getEntry(pksrec.focusAxi, pksrec.focusTan, 202 pksrec.focusRot, tmp0,tmp1,tmp2,tmp3,tmp4, 210 203 rec.asuInt("FOCUS_ID")); 211 in->molecules().getEntry(restfreq,stmp0,stmp1,rec.asuInt("MOLECULE_ID")); 212 in->tcal().getEntry(tcalt,tcalval,rec.asuInt("TCAL_ID")); 213 in->weather().getEntry(temperature, pressure, humidity, 214 windSpeed, windAz, 215 rec.asuInt("WEATHER_ID")); 204 in->molecules().getEntry(pksrec.restFreq,stmp0,stmp1, 205 rec.asuInt("MOLECULE_ID")); 206 in->tcal().getEntry(pksrec.tcalTime, pksrec.tcal, 207 rec.asuInt("TCAL_ID")); 208 in->weather().getEntry(pksrec.temperature, pksrec.pressure, 209 pksrec.humidity, pksrec.windSpeed, 210 pksrec.windAz, rec.asuInt("WEATHER_ID")); 216 211 Double pixel = Double(nchan/2); 217 Double refFreqNew = (pixel-crpix)*cdelt+ crval;212 pksrec.refFreq = (pixel-crpix)*pksrec.freqInc + crval; 218 213 // ok, now we have nrows for the n polarizations in this table 219 Matrix<Float> specs; 220 Matrix<uChar> flags; 221 Vector<Complex> xpol; 222 polConversion(specs, flags, xpol, itable); 223 Vector<Float> tsys = tsysFromTable(itable); 214 polConversion(pksrec.spectra, pksrec.flagged, pksrec.xpol, itable); 215 pksrec.tsys = tsysFromTable(itable); 224 216 // dummy data 225 uInt npol = specs.ncolumn(); 226 227 Matrix<Float> baseLin(npol,2, 0.0f); 228 Matrix<Float> baseSub(npol,9, 0.0f); 229 Complex xCalFctr; 230 Vector<Double> scanRate(2, 0.0); 231 Vector<Float> sigma(npol, 0.0f); 232 Vector<Float> calFctr(npol, 0.0f); 233 status = writer_->write(scanno, cycno, rec.asDouble("TIME"), 234 rec.asDouble("INTERVAL"), 235 rec.asString("FIELDNAME"), 236 rec.asString("SRCNAME"), 237 srcDir, srcPM, srcVel,hdr.obstype, 238 ifno, 239 refFreqNew, nchan*abs(cdelt), cdelt, 240 restfreq, 241 tcal, 242 tcalt, 243 rec.asFloat("AZIMUTH"), 244 rec.asFloat("ELEVATION"), 245 rec.asFloat("PARANGLE"), 246 focusAxi, focusTan, focusRot, 247 temperature, 248 pressure, humidity, windSpeed, windAz, 249 rec.asInt("REFBEAMNO")+1, beamno, 250 direction, 251 scanRate, 252 tsys, 253 sigma, calFctr,// not in scantable 254 baseLin, baseSub,// not in scantable 255 specs, flags, 256 xCalFctr,// 257 xpol); 217 uInt npol = spectra.ncolumn(); 218 219 pksrec.mjd = rec.asDouble("TIME"); 220 pksrec.interval = rec.asDouble("INTERVAL"); 221 pksrec.fieldName = rec.asString("FIELDNAME"); 222 pksrec.srcName = rec.asString("SRCNAME"); 223 pksrec.obsType = hdr.obstype; 224 pksrec.bandwidth = nchan * abs(pksrec.freqInc); 225 pksrec.azimuth = rec.asFloat("AZIMUTH"); 226 pksrec.elevation = rec.asFloat("ELEVATION"); 227 pksrec.parAngle = rec.asFloat("PARANGLE"); 228 pksrec.refBeam = rec.asInt("REFBEAMNO") + 1; 229 pksrec.sigma.resize(npol); 230 pksrec.sigma = 0.0f; 231 pksrec.calFctr.resize(npol); 232 pksrec.calFctr = 0.0f; 233 pksrec.baseLin.resize(npol,2); 234 pksrec.baseLin = 0.0f; 235 pksrec.baseSub.resize(npol,9); 236 pksrec.baseSub = 0.0f; 237 pksrec.xCalFctr = 0.0; 238 239 status = writer_->write(pksrec); 258 240 if ( status ) { 259 241 writer_->close(); … … 261 243 } 262 244 ++count; 263 //++ ifno;245 //++pksrec.IFno; 264 246 ++ifit; 265 247 } 266 ++ cycno;248 ++pksrec.cycleNo; 267 249 ++cycit; 268 250 } 269 //++ beamno;251 //++pksrec.beamNo; 270 252 ++beamit; 271 253 } 272 ++ scanno;254 ++pksrec.scanNo; 273 255 ++scanit; 274 256 }
Note:
See TracChangeset
for help on using the changeset viewer.