import glob import os.path import sys from itertools import izip import astropy.io.fits import lsst.afw.table as afwTable DesCatType = "AMPINFO" IsSubaru = True def massageHeader(inPath, outPath): """Massage headers""" print "Massaging header keywords in %r" % (inPath,) blob = astropy.io.fits.open(inPath) catType = blob[1].header["AFW_TYPE"] if catType != DesCatType: raise RuntimeError("catalog %r is of the wrong type: %s != %s" % (inPath, catType, DesCatType)) blob[1].header["AFW_TYPE"] = "BASE" blob.writeto(outPath, clobber=True) print "wrote stripped table to %r" % (outPath,) def convertAmpInfoCatalog(inPath, outPath): print "Converting amp info from %r" % (inPath,) oldCat = afwTable.BaseCatalog.readFits(inPath) newSchema = afwTable.AmpInfoTable.makeMinimalSchema() newCat = afwTable.AmpInfoCatalog(newSchema) newCat.reserve(len(oldCat)) for i in range(len(oldCat)): newCat.addNew() numGoodFields = 0 numBadFields = 0 oldLinCoeffsName = None for oldField in oldCat.schema.getOrderedNames(): if oldField.startswith("linearity") and oldField.endswith("coeffs"): if oldLinCoeffsName is None: oldLinCoeffsName = oldField else: raise RuntimeError("Too many fields look like linearity coeffs") newField = oldField.replace(".", "_") try: for oldRecord, newRecord in zip(oldCat, newCat): newRecord.set(newField, oldRecord.get(oldField)) numGoodFields += 1 except Exception as e: print "Failed on new field %r: %s" % (newField, e) numBadFields += 1 if numBadFields > 0: raise RuntimeError("Failed! Copied %s fields and failed to copy %s fields" % (numGoodFields, numBadFields)) print "Copied %s fields" % (numGoodFields,) if IsSubaru: print "SUBARU: setting suspect level to 2nd linearity coefficient (or nan if 0) " \ "and setting that coeff to nan" print "Old linearity coeffs=%s" % (oldCat["linearity.coeffs"]) imageTypeMax = 65535 newLinCoeffsList = [] for oldRec, newRec in izip(oldCat, newCat): oldLinCoeffs = oldRec[oldLinCoeffsName] newLinCoeffs = list(oldLinCoeffs[0:2]) + [float("nan"), float("nan")] suspectLevel = oldLinCoeffs[2] if suspectLevel == 0 or suspectLevel >= imageTypeMax: suspectLevel = float("nan") newRec.setSuspectLevel(suspectLevel) newRec.setLinearityCoeffs(newLinCoeffs) newLinCoeffsList.append(newRec.getLinearityCoeffs()) print "New linearity coeffs=%s" % (newLinCoeffsList,) print "Suspect levels=%s" % (newCat["suspectLevel"]) else: print "Not subaru: setting suspect level to nan" newCat["suspectLevel"] = float("nan") newCat.writeFits(outPath) print "Wrote new catalog to %r" % (outPath,) if __name__ == "__main__": if len(sys.argv) != 3 or sys.argv[1].startswith("-"): print "Usage: convertAmpInfoCatalog.py indir outdir" sys.exit(1) inDir, outDir = sys.argv[1:3] if os.path.abspath(inDir) == os.path.abspath(outDir): print ("Error: output directory must not match input directory") sys.exit(1) fileNames = [os.path.basename(p) for p in glob.glob(os.path.join(inDir, "*.fits"))] print "Converting %s files" % (len(fileNames),) for fn in fileNames: inPath = os.path.join(inDir, fn) tempPath = os.path.join(outDir, "temp_" + fn) outPath = os.path.join(outDir, fn) massageHeader(inPath, tempPath) convertAmpInfoCatalog(tempPath, outPath)