fits-parse 0.0.1 → 0.3.1
raw patch · 10 files changed
+972/−362 lines, 10 filesdep +containersdep +microlensdep +microlens-thdep −data-defaultdep −parser-combinatorsdep ~JuicyPixelsdep ~binarydep ~bytestringsetup-changedPVP ok
version bump matches the API change (PVP)
Dependencies added: containers, microlens, microlens-th, mtl, tasty, tasty-hunit
Dependencies removed: data-default, parser-combinators
Dependency ranges changed: JuicyPixels, binary, bytestring, fast-logger, megaparsec, optparse-applicative, statistics, text, text-latin1, vector
API changes (from Hackage documentation)
- Data.Fits: Axis :: Int -> Int -> Axis
- Data.Fits: ComplexType :: NumberType
- Data.Fits: DataString :: StringType
- Data.Fits: EmptyString :: StringType
- Data.Fits: HeaderData :: SimpleFormat -> BitPixFormat -> [Axis] -> StringValue -> StringValue -> StringValue -> StringValue -> StringValue -> StringValue -> StringValue -> StringValue -> HeaderData
- Data.Fits: IntegerType :: NumberType
- Data.Fits: Negative :: NumberModifier
- Data.Fits: NullString :: StringType
- Data.Fits: NumberValue :: NumberType -> NumberModifier -> Text -> Maybe NumberModifier -> Maybe Text -> Maybe NumberModifier -> Maybe Int -> NumberValue
- Data.Fits: Positive :: NumberModifier
- Data.Fits: RealType :: NumberType
- Data.Fits: StringValue :: StringType -> Maybe Text -> StringValue
- Data.Fits: Zero :: NumberModifier
- Data.Fits: [authorIdentifier] :: HeaderData -> StringValue
- Data.Fits: [axes] :: HeaderData -> [Axis]
- Data.Fits: [axisElementCount] :: Axis -> Int
- Data.Fits: [axisNumber] :: Axis -> Int
- Data.Fits: [bitPixFormat] :: HeaderData -> BitPixFormat
- Data.Fits: [exponentModifier] :: NumberValue -> Maybe NumberModifier
- Data.Fits: [exponent] :: NumberValue -> Maybe Int
- Data.Fits: [headerData] :: HeaderDataUnit -> HeaderData
- Data.Fits: [imaginaryModifier] :: NumberValue -> Maybe NumberModifier
- Data.Fits: [imaginaryPart] :: NumberValue -> Maybe Text
- Data.Fits: [instrumentIdentifier] :: HeaderData -> StringValue
- Data.Fits: [numberType] :: NumberValue -> NumberType
- Data.Fits: [objectIdentifier] :: HeaderData -> StringValue
- Data.Fits: [observationDate] :: HeaderData -> StringValue
- Data.Fits: [observerIdentifier] :: HeaderData -> StringValue
- Data.Fits: [originIdentifier] :: HeaderData -> StringValue
- Data.Fits: [payloadData] :: HeaderDataUnit -> ByteString
- Data.Fits: [realModifier] :: NumberValue -> NumberModifier
- Data.Fits: [realPart] :: NumberValue -> Text
- Data.Fits: [referenceString] :: HeaderData -> StringValue
- Data.Fits: [simpleFormat] :: HeaderData -> SimpleFormat
- Data.Fits: [stringType] :: StringValue -> StringType
- Data.Fits: [stringValue] :: StringValue -> Maybe Text
- Data.Fits: [telescopeIdentifier] :: HeaderData -> StringValue
- Data.Fits: data Axis
- Data.Fits: data HeaderData
- Data.Fits: data NumberModifier
- Data.Fits: data NumberType
- Data.Fits: data NumberValue
- Data.Fits: data StringType
- Data.Fits: data StringValue
- Data.Fits: instance Data.Default.Class.Default Data.Fits.Axis
- Data.Fits: instance Data.Default.Class.Default Data.Fits.HeaderData
- Data.Fits: instance Data.Default.Class.Default Data.Fits.NumberValue
- Data.Fits: instance Data.Default.Class.Default Data.Fits.StringValue
- Data.Fits: instance GHC.Show.Show Data.Fits.StringType
- Data.Fits: instance GHC.Show.Show Data.Fits.StringValue
+ Data.Fits: BinTable :: Int -> ByteString -> Extension
+ Data.Fits: Comment :: Text -> Comment
+ Data.Fits: Dimensions :: BitPixFormat -> Axes -> Dimensions
+ Data.Fits: F :: LogicalConstant
+ Data.Fits: Float :: Float -> Value
+ Data.Fits: Header :: Map Keyword Value -> Header
+ Data.Fits: Image :: Extension
+ Data.Fits: Integer :: Int -> Value
+ Data.Fits: Keyword :: Text -> Keyword
+ Data.Fits: Logic :: LogicalConstant -> Value
+ Data.Fits: Primary :: Extension
+ Data.Fits: String :: Text -> Value
+ Data.Fits: T :: LogicalConstant
+ Data.Fits: [_axes] :: Dimensions -> Axes
+ Data.Fits: [_bitpix] :: Dimensions -> BitPixFormat
+ Data.Fits: [_dimensions] :: HeaderDataUnit -> Dimensions
+ Data.Fits: [_extension] :: HeaderDataUnit -> Extension
+ Data.Fits: [_header] :: HeaderDataUnit -> Header
+ Data.Fits: [_keywords] :: Header -> Map Keyword Value
+ Data.Fits: [_mainData] :: HeaderDataUnit -> ByteString
+ Data.Fits: [heap] :: Extension -> ByteString
+ Data.Fits: [pCount] :: Extension -> Int
+ Data.Fits: axes :: Lens' Dimensions Axes
+ Data.Fits: bitPixToByteSize :: BitPixFormat -> Int
+ Data.Fits: bitpix :: Lens' Dimensions BitPixFormat
+ Data.Fits: data Dimensions
+ Data.Fits: data Extension
+ Data.Fits: data LogicalConstant
+ Data.Fits: data Value
+ Data.Fits: dimensions :: Lens' HeaderDataUnit Dimensions
+ Data.Fits: extension :: Lens' HeaderDataUnit Extension
+ Data.Fits: header :: Lens' HeaderDataUnit Header
+ Data.Fits: instance Data.String.IsString Data.Fits.Comment
+ Data.Fits: instance Data.String.IsString Data.Fits.Keyword
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.BitPixFormat
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.Comment
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.Dimensions
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.Extension
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.Header
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.Keyword
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.LogicalConstant
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.SimpleFormat
+ Data.Fits: instance GHC.Classes.Eq Data.Fits.Value
+ Data.Fits: instance GHC.Classes.Ord Data.Fits.Comment
+ Data.Fits: instance GHC.Classes.Ord Data.Fits.Keyword
+ Data.Fits: instance GHC.Show.Show Data.Fits.Comment
+ Data.Fits: instance GHC.Show.Show Data.Fits.Dimensions
+ Data.Fits: instance GHC.Show.Show Data.Fits.Extension
+ Data.Fits: instance GHC.Show.Show Data.Fits.Header
+ Data.Fits: instance GHC.Show.Show Data.Fits.HeaderDataUnit
+ Data.Fits: instance GHC.Show.Show Data.Fits.Keyword
+ Data.Fits: instance GHC.Show.Show Data.Fits.LogicalConstant
+ Data.Fits: instance GHC.Show.Show Data.Fits.SimpleFormat
+ Data.Fits: instance GHC.Show.Show Data.Fits.Value
+ Data.Fits: keywords :: Lens' Header (Map Keyword Value)
+ Data.Fits: lookup :: Keyword -> Header -> Maybe Value
+ Data.Fits: mainData :: Lens' HeaderDataUnit ByteString
+ Data.Fits: newtype Comment
+ Data.Fits: newtype Header
+ Data.Fits: newtype Keyword
+ Data.Fits: pixDimsByCol :: Axes -> [Int]
+ Data.Fits: pixDimsByRow :: Axes -> [Int]
+ Data.Fits: toFloat :: Value -> Maybe Float
+ Data.Fits: toInt :: Value -> Maybe Int
+ Data.Fits: toText :: Value -> Maybe Text
+ Data.Fits: type Axes = [Int]
+ Data.Fits.Read: InvalidData :: String -> FitsError
+ Data.Fits.Read: InvalidKey :: Keyword -> Value -> FitsError
+ Data.Fits.Read: MissingHDU :: String -> Int -> FitsError
+ Data.Fits.Read: MissingKey :: Keyword -> FitsError
+ Data.Fits.Read: ParseError :: ParseErr -> FitsError
+ Data.Fits.Read: data FitsError
+ Data.Fits.Read: eitherFail :: MonadFail m => Either String a -> m a
+ Data.Fits.Read: getHDU :: String -> Int -> [HeaderDataUnit] -> Either String HeaderDataUnit
+ Data.Fits.Read: getKeyword :: Text -> (Value -> Maybe a) -> HeaderDataUnit -> Either String a
+ Data.Fits.Read: instance GHC.Classes.Eq Data.Fits.Read.FitsError
+ Data.Fits.Read: instance GHC.Show.Show Data.Fits.Read.FitsError
+ Data.Fits.Read: maybeError :: FitsError -> Maybe a -> Either String a
+ Data.Fits.Read: readHDUs :: ByteString -> Either String [HeaderDataUnit]
+ Data.Fits.Read: readPrimaryHDU :: ByteString -> Either String HeaderDataUnit
- Data.Fits: HeaderDataUnit :: HeaderData -> ByteString -> HeaderDataUnit
+ Data.Fits: HeaderDataUnit :: Header -> Dimensions -> Extension -> ByteString -> HeaderDataUnit
- Data.Fits: bitPixToWordSize :: BitPixFormat -> Natural
+ Data.Fits: bitPixToWordSize :: BitPixFormat -> Int
Files
- Setup.hs +0/−2
- examples/omnibus/Main.hs +26/−19
- fits-parse.cabal +67/−48
- fits_files/Spiral_2_30_0_300_10_0_NoGrad.fits binary
- fits_files/nonconformant/testkeys.fits binary
- fits_files/nonconformant/testkeys2.fits binary
- src/Data/Fits.hs +159/−144
- src/Data/Fits/MegaParser.hs +240/−149
- src/Data/Fits/Read.hs +94/−0
- test/Spec.hs +386/−0
− Setup.hs
@@ -1,2 +0,0 @@-import Distribution.Simple-main = defaultMain
examples/omnibus/Main.hs view
@@ -19,15 +19,20 @@ ---- statistics import Statistics.Sample ( range, mean, variance, stdDev, stdErrMean ) -import System.Log.FastLogger( TimedFastLogger, ToLogStr, LogType( LogStderr )+---- fastlogger+import System.Log.FastLogger( TimedFastLogger, ToLogStr, LogType'( LogStderr ) , defaultBufSize, newTimeCache, simpleTimeFormat , toLogStr, newTimedFastLogger, withTimedFastLogger ) +---- microlens+import Lens.Micro ((^.)) -- local library imports-import Data.Fits ( HeaderDataUnit(..), HeaderData(..), Axis(..)- , parsePix, isBitPixInt, isBitPixFloat, pixsUnwrapI, pixsUnwrapD )-import Data.Fits.MegaParser ( getAllHDUs )+import qualified Data.Fits as FITS+import Data.Fits ( HeaderDataUnit(..), dimensions, extension, mainData+ , Axes, Dimensions(..), axes, bitpix+ , pixDimsByCol, parsePix, isBitPixInt, isBitPixFloat, pixsUnwrapI, pixsUnwrapD )+import Data.Fits.Read ( readHDUs ) -- | Paramaterized input type for files or standard input. data Input = FileInput FilePath -- ^ The 'FileInput' constructor needs a path name@@ -86,9 +91,11 @@ withTimedFastLogger timeCache (LogStderr defaultBufSize) $ \logger -> do fits <- bs i myLog logger $ "[DEBUG] input file size " ++ show (BS.length fits) ++ " bytes\n"- hdus <- getAllHDUs fits- myLog logger ("[DEBUG] found " ++ show (length hdus) ++ " hdu record(s)\n")- mapM_ (processHDU logger) hdus+ case readHDUs fits of+ Left err -> myLog logger ("[ERROR] cannot parse HDUs " ++ show err)+ Right hdus -> do+ myLog logger ("[DEBUG] found " ++ show (length hdus) ++ " hdu record(s)\n")+ mapM_ (processHDU logger) hdus where bs (FileInput f) = BS.readFile f bs StdInput = BS.hGetContents stdin@@ -97,9 +104,9 @@ processHDU logger hdu = do myLog logger $ "[DEBUG] Bit Format " ++ show bpf ++ "\n" myLog logger $ "[DEBUG] data block size " ++ show (BS.length pd) ++ " bytes\n"- myLog logger $ "[DEBUG] " ++ show (length ax) ++ " Axes\n"- unless (null ax) (mapM_ logAxes ax)- let pixCount = foldr (\x acc -> axisElementCount x * acc) 1 ax+ myLog logger $ "[DEBUG] " ++ show (length (hdu ^. dimensions . axes)) ++ " Axes\n"+ unless (null ax) (mapM_ logAxes (hdu ^. dimensions . axes))+ let pixCount = sum (pixDimsByCol $ (hdu ^. dimensions . axes)) pixs <- parsePix pixCount bpf (LBS.fromStrict pd) let pxsI = if isBitPixInt bpf then pixsUnwrapI bpf pixs else [] pxsD = if isBitPixFloat bpf then pixsUnwrapD bpf pixs else []@@ -117,24 +124,24 @@ myLog logger "[DEBUG] skipping bitmap analysis.\n" where- hd = headerData hdu- pd = payloadData hdu- ax = axes hd- bpf = bitPixFormat hd+ hd = hdu ^. FITS.header+ pd = hdu ^. mainData+ ax = hdu ^. dimensions . axes+ bpf = hdu ^. dimensions . bitpix logAxes a = myLog logger $ "[DEBUG] Axis: "- ++ show (axisNumber a)+ ++ show a ++ " count: "- ++ show (axisElementCount a) ++ "\n"+ ++ show a ++ "\n" {- | If we happen to be working in floating point 2D, let's try the following. -} bitMapProcess :: TimedFastLogger- -> [Axis] -- ^ Metadata about the column oriented axes+ -> Axes -- ^ Metadata about the column oriented axes -> V.Vector Double -- ^ Data is stored in column-row major order -> IO ()-bitMapProcess logger [] _ = myLog logger "[ERROR] BitMap processing run with no axes.\n"-bitMapProcess logger [_] _ = myLog logger "[ERROR] BitMap processing run with only one axis.\n"+bitMapProcess logger [] _ = myLog logger "[ERROR] BitMap processing run with no axes.\n"+bitMapProcess logger [_] _ = myLog logger "[ERROR] BitMap processing run with only one axis.\n" bitMapProcess logger (y:x) v = do myLog logger $ "[DEBUG] Mean: " ++ show (mean v) ++ "\n" myLog logger $ "[DEBUG] Range: " ++ show (range v) ++ "\n"
fits-parse.cabal view
@@ -1,67 +1,86 @@-cabal-version: >=1.10-name: fits-parse-version: 0.0.1-license: BSD2-license-file: LICENSE-copyright: Copyright (c) 2019 Zac Slade-maintainer: krakrjak@gmail.com-author: Zac Slade-homepage: https://github.com/krakrjak/fits-parse#readme-synopsis: Parse FITS files-description:- Parse and manipulate FITS data natively in Haskell-category: Science-build-type: Simple+cabal-version: 1.12+name: fits-parse+version: 0.3.1+license: BSD2+license-file: LICENSE+copyright: Copyright (c) 2023 Zac Slade+maintainer:+ Zac Slade <krakrjak@gmail.com>,+ Sean Hess++author: Zac Slade+homepage: https://github.com/krakrjak/fits-parse#readme+synopsis: Parse FITS files+description: Parse and manipulate FITS data natively in Haskell+category: Science+build-type: Simple extra-source-files:+ README.md+ docs/Makefile+ docs/make.bat docs/conf.py- docs/examples/omnibus.rst docs/index.rst- docs/make.bat- docs/Makefile- fits_files/nonconformant/testkeys.fits- fits_files/nonconformant/testkeys2.fits- fits_files/Spiral_2_30_0_300_10_0_NoGrad.fits- README.md+ docs/examples/omnibus.rst flag examples- description:- Do you want to build the examples?- default: False- manual: True+ description: Do you want to build the examples?+ default: False+ manual: True library exposed-modules: Data.Fits Data.Fits.MegaParser- hs-source-dirs: src- other-modules:- Paths_fits_parse+ Data.Fits.Read++ hs-source-dirs: src+ other-modules: Paths_fits_parse default-language: Haskell2010 build-depends: base >=4.7 && <5,- binary >=0.8.6.0 && <0.9,- bytestring >=0.10.8.2 && <0.11,- data-default >=0.7.1.1 && <0.8,- megaparsec >=7.0.4 && <7.1,- parser-combinators >=1.0.1 && <1.1,- text >=1.2.3.1 && <1.3,- text-latin1 >=0.3.1 && <0.4+ binary <0.9,+ bytestring <0.12,+ containers <0.7,+ megaparsec <9.4,+ microlens <0.5,+ microlens-th <0.5,+ text <2.1,+ text-latin1 <0.4 executable omnibus- main-is: Main.hs- hs-source-dirs: examples/omnibus- other-modules:- Paths_fits_parse+ main-is: Main.hs+ hs-source-dirs: examples/omnibus+ other-modules: Paths_fits_parse default-language: Haskell2010 build-depends:- JuicyPixels >=3.3.3 && <3.4,+ JuicyPixels <3.4, base >=4.7 && <5,- bytestring >=0.10.8.2 && <0.11,- fast-logger >=2.4.13 && <2.5,- fits-parse -any,- optparse-applicative >=0.14.3.0 && <0.15,- statistics >=0.15.0.0 && <0.16,- vector >=0.12.0.2 && <0.13- + bytestring <0.12,+ fast-logger <3.3,+ fits-parse,+ microlens <0.5,+ microlens-th <0.5,+ optparse-applicative <0.18,+ statistics <0.17,+ vector <0.14+ if flag(examples)- buildable: True++test-suite fits-tests+ type: exitcode-stdio-1.0+ main-is: Spec.hs+ hs-source-dirs: test+ other-modules: Paths_fits_parse+ default-language: Haskell2010+ build-depends:+ base >=4.7 && <5,+ bytestring <0.12,+ containers <0.7,+ fits-parse,+ megaparsec <9.4,+ microlens <0.5,+ microlens-th <0.5,+ mtl <2.3,+ tasty <1.5,+ tasty-hunit <0.11,+ text <2.1
− fits_files/Spiral_2_30_0_300_10_0_NoGrad.fits
binary file changed (1545444 → absent bytes)
− fits_files/nonconformant/testkeys.fits
binary file changed (195840 → absent bytes)
− fits_files/nonconformant/testkeys2.fits
binary file changed (195840 → absent bytes)
src/Data/Fits.hs view
@@ -1,7 +1,7 @@ {-| Module : Data.Fits Description : Types for FITS Data Units-Copyright : (c) Zac Slade, 2018+Copyright : (c) Zac Slade, 2023 License : BSD2 Maintainer : krakrjak@gmail.com Stability : experimental@@ -9,9 +9,11 @@ Definitions for the data types needed to parse an HDU in a FITS block. -} -{-# LANGUAGE PartialTypeSignatures, DataKinds, ExistentialQuantification- , ScopedTypeVariables, GADTs- , OverloadedStrings, TypeOperators, TypeFamilies #-}+{-# LANGUAGE+ GeneralizedNewtypeDeriving+ , OverloadedStrings+ , TemplateHaskell+#-} module Data.Fits ( -- * Data payload functions parsePix@@ -20,23 +22,37 @@ -- * Main data types , HeaderDataUnit(..)- , HeaderData(..)- , BitPixFormat(..)+ -- ^ lens exports+ , dimensions+ , header+ , extension+ , mainData , Pix(..) -- ** Header Data Types+ , Header(..)+ , keywords -- ^ lens for Keyword Map in Header+ , Extension(..)+ , Data.Fits.lookup+ , Keyword(..)+ , Value(..)+ , toInt, toFloat, toText+ , LogicalConstant(..)+ , Dimensions(..)+ , axes+ , bitpix+ , Comment(..) , SimpleFormat(..)- , Axis(..)- , StringType(..)- , StringValue(..)- , NumberType(..)- , NumberModifier(..)- , NumberValue(..)+ , BitPixFormat(..)+ , Axes -- * Utility , isBitPixInt , isBitPixFloat , bitPixToWordSize+ , bitPixToByteSize+ , pixDimsByCol+ , pixDimsByRow -- ** Constants , hduRecordLength@@ -50,23 +66,27 @@ ---- bytestring import qualified Data.ByteString as BS import qualified Data.ByteString.Lazy as BL+import qualified Data.Map as Map ----- base-import Numeric.Natural ( Natural )+import Data.String (IsString) ---- ghc import GHC.TypeNats (KnownNat, Nat) ---- text import Data.Text ( Text )+import Data.Map ( Map )+import Data.List ( intercalate ) ---- bytestring import Data.ByteString ( ByteString ) ----- default-import Data.Default( Default, def )+---- microlens+import Lens.Micro ((^.))+---- microlens-th+import Lens.Micro.TH ( makeLenses ) ----- binaryi+ import Data.Binary import Data.Binary.Get @@ -93,99 +113,40 @@ hduBlockSize :: Int hduBlockSize = hduRecordLength * hduMaxRecords -{-| There are many types of strings defined in the FITS documentation.- Refered to as "character string(s)" in the documentation, they can be- null, empty, undefined, or contain characters (printable ASCII only).--} -data StringType = NullString -- ^ The NULL character string (e.g. AUTHOR=)- | EmptyString -- ^ An empty string (e.g. AUTHOR="")- | DataString -- ^ Plain ASCII data -instance Show StringType where- show NullString = "Null String"- show EmptyString = "Empty quoted String"- show DataString = "String"---- | A 'StringValue' is a type paired with a possible value.-data StringValue = StringValue- { stringType :: StringType -- ^ Which 'StringType' is this value?- , stringValue :: Maybe Text -- ^ The payload of a character string- }--{-| The default instance of 'StringValue' is an undefined string with no- text.--}-instance Default StringValue where- def = StringValue NullString Nothing--instance Show StringValue where- show (StringValue NullString _) = show NullString- show (StringValue EmptyString _) = show EmptyString- show (StringValue DataString Nothing) = "No good " ++ show DataString- show (StringValue DataString (Just s)) = show DataString ++ T.unpack s--{-| The FITS standard allows for the encoding of unsigned integers, signed- integers, real numbers, and complex numbers. They are always ASCII- encoded. See 5.2 of the standard for more details.--}-data NumberType =- IntegerType -- ^ HDU ASCII encoded integer number- | RealType -- ^ HDU ASCII encoded real number- | ComplexType -- ^ HDU ASCII encoded complex number---- | Utility data type to help with the ASCII representation of numbers-data NumberModifier =- Positive -- ^ HDU positive number value- | Negative -- ^ HDU negative number value- | Zero -- ^ HDU numeric value is zero, could be positive or negative--{-| 'NumberValue' contains an encoded numeric record from a data field.- This data type still needs to be converted into more useful Haskell data- types.--}-data NumberValue = NumberValue- { numberType :: NumberType- -- ^ Key to decoding the structure- , realModifier :: NumberModifier- -- ^ Encoding the sign of the real part- , realPart :: Text- -- ^ All 'NumberType' have a real part (sign stripped)- , imaginaryModifier :: Maybe NumberModifier- -- ^ Encoding the sign of the imaginary part- , imaginaryPart :: Maybe Text- -- ^ Only 'ComplexType' have an imaginary part (sign stripped)- , exponentModifier :: Maybe NumberModifier- -- ^ 'Positive', 'Negative', or 'Zero'- , exponent :: Maybe Int- -- ^ All 'NumberType' may have an exponent- }--{-| The default instance for 'NumberValue' is a 32-bit integer of zero with- no imaginary or exponent parts.--}-instance Default NumberValue where- def = NumberValue IntegerType Zero "0" Nothing Nothing Nothing Nothing- {-| The standard defines two possible values for the SIMPLE keyword, T and F. The T refers to a 'Conformant' format while F refers to a 'NonConformant' format. At this time only the 'Conformant', T, format is supported. -}-data SimpleFormat = Conformant+data SimpleFormat = Conformant | NonConformant+ deriving (Eq, Show) -- ^ Value of SIMPLE=T in the header. /supported/- | NonConformant+ -- NonConformat -- ^ Value of SIMPLE=F in the header. /unsupported/ --- | 'Axis' represents a single NAXIS record.-data Axis = Axis- { axisNumber :: Int -- ^ The axis number under consideration- , axisElementCount :: Int -- ^ The number of elements in this axis- }+{-| Direct encoding of a `Bool` for parsing `Value` -}+data LogicalConstant = T | F+ deriving (Show, Eq) --- | The default instance for 'Axis' is NAXIS=0 with zero elements.-instance Default Axis where- def = Axis 0 0+{-| The `Text` wrapper for HDU the keyword data for lines of the form:+ KEYWORD=VALUE+-}+newtype Keyword = Keyword Text+ deriving (Show, Eq, Ord, IsString) +{-| `Value` datatype for discriminating valid FITS KEYWORD=VALUE types in an HDU. -}+data Value+ = Integer Int+ | Float Float+ | String Text+ | Logic LogicalConstant+ deriving (Show, Eq)+++{-| 'Axes' represents the combination of NAXIS + NAXISn. The spec supports up to 999 axes -}+type Axes = [Int]+ {-| The 'BitPixFormat' is the nitty gritty of how the 'Axis' data is layed out in the file. The standard recognizes six formats: unsigned 8 bit integer, two's complement binary integers at 16, 32, and 64 bits along@@ -198,19 +159,20 @@ | SixtyFourBitInt -- ^ BITPIX = 64; two's complement binary integer of 64 bits | ThirtyTwoBitFloat -- ^ BITPIX = -32; IEEE single precision floating point of 32 bits | SixtyFourBitFloat -- ^ BITPIX = -64; IEEE double precision floating point of 64 bits+ deriving (Eq) instance Show BitPixFormat where- show EightBitInt = "8 bit unsigned integer"- show SixteenBitInt = "16 bit signed integer"- show ThirtyTwoBitInt = "32 bit signed integer"- show SixtyFourBitInt = "64 bit signed interger"- show ThirtyTwoBitFloat = "32 bit IEEE single precision float"- show SixtyFourBitFloat = "64 bit IEEE double precision float"+ show EightBitInt = "8 bit unsigned integer"+ show SixteenBitInt = "16 bit signed integer"+ show ThirtyTwoBitInt = "32 bit signed integer"+ show SixtyFourBitInt = "64 bit signed interger"+ show ThirtyTwoBitFloat = "32 bit IEEE single precision float"+ show SixtyFourBitFloat = "64 bit IEEE double precision float" {-| This utility function can be used to get the word count for data in an HDU. -}-bitPixToWordSize :: BitPixFormat -> Natural+bitPixToWordSize :: BitPixFormat -> Int bitPixToWordSize EightBitInt = 8 bitPixToWordSize SixteenBitInt = 16 bitPixToWordSize ThirtyTwoBitInt = 32@@ -221,7 +183,7 @@ {-| This utility function can be used to get the size in bytes of the - format. -}-bitPixToByteSize :: BitPixFormat -> Natural+bitPixToByteSize :: BitPixFormat -> Int bitPixToByteSize EightBitInt = 1 bitPixToByteSize SixteenBitInt = 2 bitPixToByteSize ThirtyTwoBitInt = 4@@ -308,57 +270,110 @@ {- `pixDimsByCol` takes a list of Axis and gives a column-row major list of axes dimensions. -}-pixDimsByCol :: [Axis] -> [Int]-pixDimsByCol = map axisElementCount+pixDimsByCol :: Axes -> [Int]+pixDimsByCol = id {- `pixDimsByRow` takes a list of Axis and gives a row-column major list of axes dimensions. -}-pixDimsByRow :: [Axis] -> [Int]+pixDimsByRow :: Axes -> [Int] pixDimsByRow = reverse . pixDimsByCol {-| The header part of the HDU is vital carrying not only authorship metadata, but also specifying how to make sense of the binary payload that starts 2,880 bytes after the start of the 'HeaderData'. -}+newtype Header = Header { _keywords :: Map Keyword Value }+ deriving (Eq)+$(makeLenses ''Header) -data HeaderData = HeaderData- { simpleFormat :: SimpleFormat- -- ^ SIMPLE- , bitPixFormat :: BitPixFormat- -- ^ BITPIX- , axes :: [Axis]- -- ^ Axes metadata- , objectIdentifier :: StringValue- -- ^ OBJECT- , observationDate :: StringValue- -- ^ DATE- , originIdentifier :: StringValue- -- ^ OBJECT- , telescopeIdentifier :: StringValue- -- ^ TELESCOP- , instrumentIdentifier :: StringValue- -- ^ INSTRUME- , observerIdentifier :: StringValue- -- ^ OBSERVER- , authorIdentifier :: StringValue- -- ^ CREATOR- , referenceString :: StringValue- -- ^ REFERENC- }+instance Show Header where+ show h =+ let kvs = Map.toList (h ^. keywords) :: [(Keyword, Value)]+ in T.unpack $ T.intercalate "\n" $ fmap line kvs+ where+ --+ -- init :: [Text]+ -- init = map T.pack+ -- [ "BITPIX =" <> show h.size.bitpix+ -- , "NAXES =" <> show h.size.naxes+ -- ] -instance Default HeaderData where- def = HeaderData NonConformant EightBitInt []- (def :: StringValue) (def :: StringValue) (def :: StringValue)- (def :: StringValue) (def :: StringValue) (def :: StringValue)- (def :: StringValue) (def :: StringValue)+ line :: (Keyword, Value) -> Text+ line (Keyword k, v) =+ T.justifyLeft 8 ' ' k+ <> "="+ <> T.justifyLeft (hduRecordLength - 10) ' ' (T.pack $ val v) + val (Integer n) = show n+ val (Float f) = show f+ val (Logic T) = " T"+ val (String t) = T.unpack t++lookup :: Keyword -> Header -> Maybe Value+lookup k h = Map.lookup k (h ^. keywords)+++data Extension+ -- | Any header data unit can use the primary format. The first MUST be+ -- Primary. This is equivalent to having no extension+ = Primary++ -- | An encoded image. PCOUNT and GCOUNT are required but irrelevant+ | Image++ -- | A Binary table. PCOUNT is the number of bytes that follow the data+ -- in the 'heap'+ | BinTable { pCount :: Int, heap :: ByteString }+ deriving (Eq)++instance Show Extension where+ show Primary = "Primary"+ show Image = "Image"+ show (BinTable p _) = "BinTable: heap = " <> show p <> " Bytes"+++toInt :: Value -> Maybe Int+toInt (Integer i) = Just i+toInt _ = Nothing++toFloat :: Value -> Maybe Float+toFloat (Float n) = Just n+toFloat _ = Nothing++toText :: Value -> Maybe Text+toText (String s) = Just s+toText _ = Nothing++{-| When we load a header, we parse the BITPIX and NAXIS(N) keywords so we+ - can know how long the data array is+-}+data Dimensions = Dimensions+ { _bitpix :: BitPixFormat+ , _axes :: Axes+ } deriving (Show, Eq)+$(makeLenses ''Dimensions)++newtype Comment = Comment Text+ deriving (Show, Eq, Ord, IsString)++ {-| The 'HeaderDataUnit' is the full HDU. Both the header information is- encoded alongside the 'Axis' payload.+ encoded alongside the data payload. -} data HeaderDataUnit = HeaderDataUnit- { headerData :: HeaderData- -- ^ Just the header part of the HDU- , payloadData :: ByteString- -- ^ The actual data payload+ { _header :: Header -- ^ The heeader contains metadata about the payload+ , _dimensions :: Dimensions -- ^ This dimensions of the main data array+ , _extension :: Extension -- ^ Extensions may vary the data format+ , _mainData :: ByteString -- ^ The main data array }+ +$(makeLenses ''HeaderDataUnit)++instance Show HeaderDataUnit where+ show hdu = intercalate "\n" + [ "HeaderDataUnit:"+ , " headers = " <> show (Map.size (hdu ^. header . keywords))+ , " extension = " <> show (hdu ^. extension)+ , " mainData = " <> show (BS.length (hdu ^. mainData)) <> " Bytes"+ ]
src/Data/Fits/MegaParser.hs view
@@ -1,7 +1,7 @@ {-| Module : Data.Fits.MegaParser Description : MegaParsec based parser for an HDU.-Copyright : (c) Zac Slade, 2018+Copyright : (c) Zac Slade, 2023 License : BSD2 Maintainer : krakrjak@gmail.com Stability : experimental@@ -15,39 +15,66 @@ module Data.Fits.MegaParser where -- qualified imports+---- bytestring import qualified Data.ByteString as BS import qualified Data.ByteString.Internal as BS ( c2w )+---- containers+import qualified Data.Map.Lazy as Map+---- megaparsec import qualified Text.Megaparsec as M import qualified Text.Megaparsec.Stream as M-import qualified Text.Megaparsec.Char as M-import qualified Text.Megaparsec.Char.Lexer as MCL+import qualified Text.Megaparsec.Pos as MP+import qualified Text.Megaparsec.Byte as M+import qualified Text.Megaparsec.Byte.Lexer as MBL+---- text import qualified Data.Text as T import qualified Data.Text.Encoding as TE+---- local imports+import qualified Data.Fits as Fits+import qualified Data.Text.Encoding as C8+import qualified Data.Binary as C8 --- explicit imports-import Control.Applicative ( (<$>) )-import Control.Monad ( void )-import Numeric.Natural ( Natural )++-- symbol imports+---- bytestring import Data.ByteString ( ByteString )+---- containers+import Data.Map ( Map )+---- text import Data.Text ( Text )-import Data.Word ( Word8, Word16, Word32, Word64 )-import Data.Void ( Void )-import Data.Default ( def )+---- megaparsec import Text.Ascii ( isAscii ) import Text.Megaparsec ( Parsec, ParseErrorBundle, (<|>), (<?>))---- full module imports----- parser-combinators-import Control.Applicative.Permutations+---- microlens+import Lens.Micro ((^.)) ---- base-import Data.Proxy---- local imports+import Control.Applicative ( (<$>) )+import Control.Exception ( Exception(displayException) )+import Control.Monad ( void, foldM )+import Data.Bifunctor ( first )+import Data.Char ( ord )+import Data.Maybe ( catMaybes, fromMaybe )+import Data.Word ( Word8, Word16, Word32, Word64 )+import Data.Void ( Void )+---- local imports import Data.Fits+ ( Axes+ , Comment(Comment)+ , Dimensions(Dimensions)+ , Header(Header)+ , HeaderDataUnit(HeaderDataUnit)+ , Keyword(Keyword)+ , BitPixFormat(..)+ , Extension(..)+ , LogicalConstant(..)+ , Value(..)+ , bitPixToByteSize, hduRecordLength+ ) -type Parser = Parsec Void Text-type ParseErr = ParseErrorBundle Text Void +type Parser = Parsec Void ByteString+type ParseErr = ParseErrorBundle ByteString Void+ data DataUnitValues = FITSUInt8 Word8 | FITSInt16 Word16@@ -56,164 +83,228 @@ | FITSFloat32 Float | FITSFloat64 Double -headerBlockParse :: Parser HeaderData-headerBlockParse = do- -- TODO: Parse other sections of the neader for god's sake- -- HISTORY and COMMENT along with CONTINUE handling is missing- (simple, bitpix, axesDesc, bZero, bScale, ref, obs, instr, tele, object, pCreator, pDate, pEnd) <-- runPermutation $- (,,,,,,,,,,,,) <$> toPermutation (parseSimple <?> "simple")- <*> toPermutation (parseBitPix <?> "bitpix")- <*> toPermutation ((parseAxisCount >>= parseNaxes) <?> "axis parsing")- <*> toPermutationWithDefault 0 (parseBzero <?> "bzero")- <*> toPermutationWithDefault 0 (parseBscale <?> "bscale")- <*> toPermutationWithDefault (StringValue NullString Nothing) (parseReference <?> "reference")- <*> toPermutationWithDefault (StringValue NullString Nothing) (parseObserver <?> "observer")- <*> toPermutationWithDefault (StringValue NullString Nothing) (parseInstrument <?> "instrument")- <*> toPermutationWithDefault (StringValue NullString Nothing) (parseTelescope <?> "telescope")- <*> toPermutationWithDefault (StringValue NullString Nothing) (parseObject <?> "object")- <*> toPermutationWithDefault (StringValue NullString Nothing) (parseCreator <?> "creator")- <*> toPermutationWithDefault (StringValue NullString Nothing) (parseDate <?> "date")- <*> toPermutation (parseEnd <?> "end")- return defHeader { simpleFormat = simple- , bitPixFormat = bitpix- , axes = axesDesc- , referenceString = ref- , observerIdentifier = obs- , instrumentIdentifier = instr- , telescopeIdentifier = tele- , objectIdentifier = object- , observationDate = pDate- , authorIdentifier = pCreator }- where- defHeader = def :: HeaderData -parseSimple :: Parser SimpleFormat-parseSimple = M.string' "simple" >> parseEquals- >> ((conformParse >> consumeDead >> return Conformant)- <|> (nonConformParse >> consumeDead >> return NonConformant))- where- conformParse = M.char' 't'- nonConformParse = M.anySingle+toWord :: Char -> Word8+toWord = fromIntegral . ord -parseEquals :: Parser ()-parseEquals = M.space >> M.char '=' >> M.space+wordsText :: [Word8] -> Text+wordsText = TE.decodeUtf8 . BS.pack -parseBitPix :: Parser BitPixFormat-parseBitPix = M.string' "bitpix" >> parseEquals- >> ((M.chunk "8" >> consumeDead >> return EightBitInt)- <|> (M.chunk "16" >> consumeDead >> return SixteenBitInt)- <|> (M.chunk "32" >> consumeDead >> return ThirtyTwoBitInt)- <|> (M.chunk "64" >> consumeDead >> return SixtyFourBitInt)- <|> (M.chunk "-32" >> consumeDead >> return ThirtyTwoBitFloat)- <|> (M.chunk "-64" >> consumeDead >> return SixtyFourBitFloat)) -parseAxisCount :: Parser Natural-parseAxisCount = M.string' "naxis" >> parseEquals >> parseNatural+-- | Consumes ALL header blocks until end, then all remaining space+parseHeader :: Parser (Map Keyword Value)+parseHeader = do+ pairs <- M.manyTill parseRecordLine (M.string' "end")+ M.space -- consume space padding all the way to the end of the next 2880 bytes header block+ return $ Map.fromList $ catMaybes pairs -parseNaxes :: Natural -> Parser [Axis]-parseNaxes n | n == 0 = return []-parseNaxes n = do- axisNum <- M.string' "naxis" >> parseNatural- elemCount <- parseEquals >> parseNatural- ([buildAxis axisNum elemCount] ++) <$> parseNaxes (n - 1)- where- defAxis = def :: Axis- buildAxis an ec = defAxis { axisNumber = fromIntegral an- , axisElementCount = fromIntegral ec }+parseRecordLine :: Parser (Maybe (Keyword, Value))+parseRecordLine = do+ M.try $ Just <$> parseKeywordRecord+ <|> Nothing <$ parseLineComment+ <|> Nothing <$ M.string' (BS.replicate hduRecordLength (toWord ' ')) -parseBzero :: Parser Int-parseBzero = M.string' "bzero" >> parseEquals >> parseInteger+-- | Combinator to allow for parsing a record with inline comments+withComments :: Parser a -> Parser a+withComments parse = do+ start <- parsePos+ a <- parse+ M.space+ M.optional $ parseInlineComment start+ M.space+ return a -parseBscale :: Parser Int-parseBscale = M.string' "bscale" >> parseEquals >> parseInteger+parseKeywordRecord :: Parser (Keyword, Value)+parseKeywordRecord = withComments parseKeywordValue -parseReference :: Parser StringValue-parseReference = M.string' "referenc" >> parseEquals >> parseStringValue+-- | Parses the specified keyword+parseKeywordRecord' :: ByteString -> Parser a -> Parser a+parseKeywordRecord' k pval = withComments $ do+ M.string' k+ parseEquals+ pval -parseObserver :: Parser StringValue-parseObserver = M.string' "observer" >> parseEquals >> parseStringValue+parseKeywordValue :: Parser (Keyword, Value)+parseKeywordValue = do+ key <- parseKeyword+ parseEquals+ val <- parseValue+ return (key, val) -parseInstrument :: Parser StringValue-parseInstrument = M.string' "instrume" >> parseEquals >> parseStringValue+parseInlineComment :: Int -> Parser Comment+parseInlineComment start = do+ M.char $ toWord '/'+ M.space+ com <- parsePos+ let end = start + hduRecordLength+ let rem = end - com+ c <- M.count rem M.anySingle+ return $ Comment (wordsText c) -parseTelescope :: Parser StringValue-parseTelescope = M.string' "telescop" >> parseEquals >> parseStringValue+parseLineComment :: Parser Comment+parseLineComment = do+ let keyword = "COMMENT " :: ByteString+ M.string' keyword+ c <- M.count (hduRecordLength - BS.length keyword) M.anySingle+ return $ Comment (wordsText c) -parseObject :: Parser StringValue-parseObject = M.string' "object" >> parseEquals >> parseStringValue+-- | Anything but a space or equals+parseKeyword :: Parser Keyword+parseKeyword = Keyword . wordsText <$> M.some (M.noneOf $ fmap toWord [' ', '=']) -parseCreator :: Parser StringValue-parseCreator = M.string' "creator" >> parseEquals >> parseStringValue+parseValue :: Parser Value+parseValue =+ -- try is required here because Megaparsec doesn't automatically backtrack if the parser consumes anything+ M.try (Float <$> parseFloat)+ <|> M.try (Integer <$> parseInt)+ <|> (Logic <$> parseLogic)+ <|> (String <$> parseStringContinue) -parseDate :: Parser StringValue-parseDate = M.string' "date" >> parseEquals >> parseStringValue+parseInt :: Num a => Parser a+parseInt = MBL.signed M.space MBL.decimal -skipEmpty :: Parser ()-skipEmpty = void (M.many $ M.satisfy ('\0' ==))+parseFloat :: Parser Float+parseFloat = MBL.signed M.space MBL.float -consumeDead :: Parser ()-consumeDead = M.space >> skipEmpty+parseLogic :: Parser LogicalConstant+parseLogic = do+ M.string' "T"+ return T -parseEnd :: Parser ()-parseEnd = M.string' "end" >> M.space <* M.eof+parseStringContinue :: Parser Text+parseStringContinue = do+ t <- parseStringValue -parseNatural :: Parser Natural-parseNatural = do- v <- MCL.decimal- consumeDead- return $ fromIntegral v+ mc <- M.optional $ do+ M.string' "CONTINUE"+ M.space+ parseStringContinue -parseInteger :: Parser Int-parseInteger = do- v <- MCL.decimal- consumeDead- return v+ case mc of+ Nothing -> return t+ Just tc -> return $ T.dropWhileEnd (=='&') t <> tc -parseStringValue :: Parser StringValue-parseStringValue = do +parseStringValue :: Parser Text+parseStringValue = do -- The rules are weird, NULL means a NULL string, '' is an empty -- string, a ' followed by a bunch of spaces and a close ' is -- considered an empty string, and trailing whitespace is ignored -- within the quotes, but not leading spaces.- ls <- M.between (M.char '\'') (M.char '\'') $ M.many $ M.anySingleBut '\''+ ls <- M.between (M.char quote) (M.char quote) $ M.many $ M.anySingleBut quote consumeDead- let v = M.tokensToChunk (Proxy :: Proxy Text) ls- if T.length v < 1- then return (StringValue EmptyString Nothing)- else return (StringValue DataString (Just v))+ return (T.stripEnd $ wordsText ls)+ where quote = toWord '\'' -countHeaderDataUnits :: ByteString -> IO Natural-countHeaderDataUnits bs = fromIntegral . length <$> getAllHDUs bs+requireKeyword :: Keyword -> Header -> Parser Value+requireKeyword k kvs = do+ case Fits.lookup k kvs of+ Nothing -> fail $ "Missing: " <> show k+ Just v -> return v --- TODO: make the recursive case work. Currently limited to one HDU.--- The current issue is that when the parser fails on an HDU parse, it--- blows them all up instead of accepting the valid parsings.-getAllHDUs :: ByteString -> IO [HeaderDataUnit]-getAllHDUs bs = do- (hdu, rest) <- getOneHDU bs- return [hdu]--- if BS.length rest < hduBlockSize then return [hdu] else return [hdu]+requireNaxis :: Header -> Parser Int+requireNaxis kvs = do+ v <- requireKeyword "NAXIS" kvs+ case v of+ Integer n -> return n+ _ -> fail "Invalid NAXIS header" -getOneHDU :: ByteString -> IO (HeaderDataUnit, ByteString)-getOneHDU bs =- if isAscii header- then- case M.runParser headerBlockParse "FITS" (TE.decodeUtf8 header) of- Right mainHeader -> do- let (dataUnit, remainder) = BS.splitAt (fromIntegral $ dataSize mainHeader) rest- return (HeaderDataUnit mainHeader dataUnit, remainder)- Left e -> let err = M.errorBundlePretty e in error err- else error "Header data is not ASCII. Please Check your input file and try again"+skipEmpty :: Parser ()+skipEmpty = void (M.many $ M.satisfy (toWord '\0' ==))++consumeDead :: Parser ()+consumeDead = M.space >> skipEmpty++parseEnd :: Parser ()+parseEnd = M.string' "end" >> M.space <* M.eof++parseEquals :: Parser ()+parseEquals = M.space >> M.char (toWord '=') >> M.space++parsePos :: Parser Int+parsePos = MP.unPos . MP.sourceColumn <$> M.getSourcePos+++parseBitPix :: Parser BitPixFormat+parseBitPix = do+ v <- parseKeywordRecord' "BITPIX" parseValue+ toBitpix v+ where+ toBitpix (Integer 8) = return EightBitInt+ toBitpix (Integer 16) = return SixteenBitInt+ toBitpix (Integer 32) = return ThirtyTwoBitInt+ toBitpix (Integer 64) = return SixtyFourBitInt+ toBitpix (Integer (-32)) = return ThirtyTwoBitFloat+ toBitpix (Integer (-64)) = return SixtyFourBitFloat+ toBitpix _ = fail "Invalid BITPIX header"++parseNaxes :: Parser Axes+parseNaxes = do+ n <- parseKeywordRecord' "NAXIS" parseInt+ mapM parseN [1..n] where- (header, rest) = BS.splitAt hduBlockSize bs+ parseN :: Int -> Parser Int+ parseN n = withComments $ do+ M.string' "NAXIS"+ M.string' $ BS.pack $ map toWord (show n)+ parseEquals+ parseInt -dataSize :: HeaderData -> Natural-dataSize h = paddedsize+-- | We don't parse simple here, because it isn't required on all HDUs+parseDimensions :: Parser Dimensions+parseDimensions = do+ bp <- parseBitPix+ Dimensions bp <$> parseNaxes++parsePrimary :: Parser HeaderDataUnit+parsePrimary = do+ parseKeywordRecord' "SIMPLE" parseLogic+ dm <- M.lookAhead parseDimensions+ hd <- parseHeader+ dt <- parseMainData dm+ return $ HeaderDataUnit (Header hd) dm Primary dt++parseImage :: Parser HeaderDataUnit+parseImage = do+ withComments $ M.string' "XTENSION= 'IMAGE '"+ dm <- M.lookAhead parseDimensions+ hd <- parseHeader+ dt <- parseMainData dm+ return $ HeaderDataUnit (Header hd) dm Image dt++parseBinTable :: Parser HeaderDataUnit+parseBinTable = do+ (dm, pc) <- M.lookAhead parseBinTableKeywords+ hd <- parseHeader+ dt <- parseMainData dm+ hp <- parseBinTableHeap+ let tab = BinTable pc hp+ return $ HeaderDataUnit (Header hd) dm tab dt+ where+ parseBinTableHeap = return ""++parseBinTableKeywords :: Parser (Dimensions, Int)+parseBinTableKeywords = do + withComments $ M.string' "XTENSION= 'BINTABLE'"+ sz <- parseDimensions+ pc <- parseKeywordRecord' "PCOUNT" parseInt+ return (sz, pc)++parseMainData :: Dimensions -> Parser ByteString+parseMainData size = do+ let len = dataSize size+ M.takeP (Just ("Data Array of " <> show len <> " Bytes")) (fromIntegral len)++parseHDU :: Parser HeaderDataUnit+parseHDU =+ parsePrimary <|> parseImage <|> parseBinTable++parseHDUs :: Parser [HeaderDataUnit]+parseHDUs = do+ M.many parseHDU++dataSize :: Dimensions -> Int+dataSize (Dimensions bitpix axes) = size bitpix * count axes where- axesCount = length $ axes h- wordCount = product $ map axisElementCount $ axes h- wordsize = fromIntegral . bitPixToWordSize $ bitPixFormat h- datasize = wordsize * wordCount- padding = if axesCount == 0 then 0 else fromIntegral hduBlockSize - datasize `mod` fromIntegral hduBlockSize- paddedsize = fromIntegral (datasize + padding)+ count [] = 0+ count ax = fromIntegral $ product ax+ size = fromIntegral . bitPixToByteSize
+ src/Data/Fits/Read.hs view
@@ -0,0 +1,94 @@+{-# LANGUAGE OverloadedRecordDot #-}+{-# LANGUAGE OverloadedStrings #-}+module Data.Fits.Read where++import Control.Exception ( displayException )+import Data.Bifunctor ( first )+import Data.ByteString ( ByteString )+import Data.Maybe ( listToMaybe )+import Data.Text ( Text, unpack )+import qualified Data.ByteString as BS+import qualified Data.Map.Lazy as Map+import qualified Text.Megaparsec as M+import Data.List ( find )++---- local imports+import Data.Fits as Fits+import Data.Fits.MegaParser (ParseErr(..), parseHDU, parseHDUs)+import Data.Fits (HeaderDataUnit(..))+++-- | Parse and read all HDUs in the input string+readHDUs :: ByteString -> Either String [HeaderDataUnit]+readHDUs bs = do+ first (show . ParseError) $ M.runParser parseHDUs "FITS" bs++-- | Parse and read only the Primary HDU from the input string+readPrimaryHDU :: ByteString -> Either String HeaderDataUnit+readPrimaryHDU bs = do+ first (show . ParseError) $ M.runParser parseHDU "FITS" bs++-- | Look up a keyword and parse it into the expected format+getKeyword :: Text -> (Value -> Maybe a) -> HeaderDataUnit -> Either String a+getKeyword k fromVal hdu = do+ let key = Keyword k+ v <- maybeError (MissingKey key) $ Map.lookup key (_keywords . _header $ hdu)+ maybeError (InvalidKey key v) $ fromVal v+ where+ findKey :: Keyword -> Header -> Maybe Value+ findKey key h = Map.lookup key (_keywords h)++-- | Get the HDU at an index and fail with a readable error+getHDU :: String -> Int -> [HeaderDataUnit] -> Either String HeaderDataUnit+getHDU name n hdus = do+ maybeError (MissingHDU name n) $ listToMaybe $ drop n hdus++maybeError :: FitsError -> Maybe a -> Either String a+maybeError e Nothing = Left (show e)+maybeError _ (Just a) = Right a++eitherFail :: MonadFail m => Either String a -> m a +eitherFail (Left e) = fail e+eitherFail (Right a) = return a++data FitsError+ = ParseError ParseErr+ | MissingKey Keyword+ | InvalidKey Keyword Value+ | MissingHDU String Int + | InvalidData String+ deriving (Eq)++instance Show FitsError where+ show (ParseError e) = displayException e+ show (MissingKey (Keyword k)) = "Keyword Missing: " <> unpack k+ show (InvalidKey (Keyword k) val) = "Keyword: " <> unpack k <> " was invalid. Got " <> show val+ show (MissingHDU name n) = "HDU Missing: " <> name <> " at index " <> show n+ show (InvalidData err) = "Data Invalid: " <> err++-- -- | An example of how to use the library+-- example :: IO ()+-- example = do+-- bs <- BS.readFile "./fits_files/nso_dkist.fits"+--+-- (tel, obs, dm) <- throwLeft $ exampleReadMyData bs+--+-- putStrLn $ "TELESCOPE: " <> unpack tel+-- putStrLn $ "OBSERVATORY: " <> unpack obs+-- putStrLn $ "DATAMIN: " <> show dm+--+-- where+-- throwLeft :: Show e => Either e a -> IO a+-- throwLeft (Left e) = fail $ show e+-- throwLeft (Right a) = return a+--+-- -- You can parse the file and lookup relevant data in the same function+-- exampleReadMyData :: ByteString -> Either String (Text, Text, Float)+-- exampleReadMyData bs = do+-- hdus <- readHDUs bs+-- hdu <- getHDU "Main Binary Table" 1 hdus+-- tel <- getKeyword "TELESCOP" toText hdu+-- obs <- getKeyword "OBSRVTRY" toText hdu+-- dm <- getKeyword "DATAMIN" toFloat hdu+-- return (tel, obs, dm)+--
+ test/Spec.hs view
@@ -0,0 +1,386 @@+{-# LANGUAGE GeneralizedNewtypeDeriving #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE OverloadedStrings #-}++module Main where+-- qualified imports+---- base+import qualified Data.Map.Strict as Map+---- bytestring+import qualified Data.ByteString as BS+import qualified Data.ByteString.Char8 as C8+---- megaparsec+import qualified Text.Megaparsec as M+import qualified Text.Megaparsec.Char as M+---- text+import qualified Data.Text as T+import qualified Data.Text.Encoding as TE+---- local imports+import qualified Data.Fits as Fits+-- symbol based imports+---- base+import Control.Exception (Exception(displayException), throwIO)+import Control.Monad.Writer+import Data.List ( unfoldr )+---- bytestring+import Data.ByteString ( ByteString )+---- ghc+import GHC.RTS.Flags (MiscFlags(numIoWorkerThreads))+---- microlens+import Lens.Micro ((^.))+---- tasty+import Test.Tasty+import Test.Tasty.HUnit+---- text+import Data.Text ( Text )+---- local-imports+import Data.Fits ( Axes+ , BitPixFormat(..)+ , bitPixToByteSize+ , Comment(Comment)+ , Dimensions+ , axes+ , bitpix+ , Extension(..)+ , Header+ , keywords+ , HeaderDataUnit+ , header+ , dimensions+ , extension+ , mainData+ , Keyword(Keyword)+ , LogicalConstant(..)+ , Value(..)+ , hduBlockSize+ )+import Data.Fits.MegaParser+import Data.Fits.Read++main :: IO ()+main =+ testMain $ runTests "Tests" $ do+ basicParsing+ keywordValueLines+ comments+ continue+ fullRecord+ fullRecordLine+ headerMap+ requiredHeaders+ dataArray+ sampleSpiral+ sampleNSOHeaders+ sampleNSO++parse :: Parser a -> ByteString -> IO a+parse p inp =+ case M.parse p "Test" inp of+ Left e -> fail $ displayException e+ Right v -> pure v++flattenKeywords :: [ByteString] -> ByteString+flattenKeywords ts = mconcat (map pad ts) <> "END"++pad :: ByteString -> ByteString+pad m =+ let n = 80 - BS.length m+ in m <> C8.replicate n ' '+++basicParsing :: Test ()+basicParsing = describe "Basic Parsing" $ do+ it "should parse a string" $ do+ res <- parse parseStringValue "'hello there' "+ res @?= "hello there"++ it "should parse a number value" $ do+ res <- parse parseValue "42 "+ res @?= Integer 42++ it "should parse a keyword" $ do+ res <- parse parseKeyword "WS_TEMP ="+ res @?= Keyword "WS_TEMP"++ it "should handle keyword symbols" $ do+ res <- parse parseKeyword "OBSGEO-X= -5466045.256954942 / [m]"+ res @?= Keyword "OBSGEO-X"++++keywordValueLines :: Test ()+keywordValueLines = describe "parse keyword=value" $ do+ it "should parse an integer" $ do+ res <- parse parseKeywordValue "KEY=42 "+ res @?= ("KEY", Integer 42)++ it "should parse a string" $ do+ res <- parse parseKeywordValue "KEY='value'" + res @?= ("KEY", String "value")++ it "should absorb spaces" $ do+ res <- parse parseKeywordValue "KEY = 'value' "+ res @?= ("KEY", String "value")++ it "should parse a float" $ do+ res <- parse parseKeywordValue "KEY = 44.88 "+ res @?= ("KEY", Float 44.88)++ it "should parse a negative number" $ do+ res <- parse parseKeywordValue "KEY = -44.88"+ res @?= ("KEY", Float ( -44.88 ))++ it "should parse a logical constant" $ do+ res <- parse parseKeywordValue "KEY= T " + res @?= ("KEY", Logic T)++ it "should ignore comments" $ do+ res <- parse parseKeywordValue "SIMPLE = T / conforms to FITS standard" + res @?= ("SIMPLE", Logic T)++ it "should strip trailing spaces from strings" $ do+ res <- parse parseKeywordValue "INSTRUME= 'VISP '"+ res @?= ("INSTRUME", String "VISP")++fullRecord :: Test ()+fullRecord = describe "parseKeywordRecord" $ do+ it "should parse an 80 character record" $ do+ res <- parse parseKeywordRecord (flattenKeywords ["KEYWORD = 12345"])+ res @?= ("KEYWORD", Integer 12345)++ it "should parse an a record and comment" $ do+ res <- parse parseKeywordRecord (flattenKeywords ["KEYWORD = 12345 / this is a comment"])+ res @?= ("KEYWORD", Integer 12345)++ it "should parse a record, comment, followed by next keyword" $ do+ res <- parse parseKeywordRecord $ flattenKeywords ["SIMPLE = T / conforms to FITS standard"]+ res @?= ("SIMPLE", Logic T)++ it "should handle keyword symbols" $ do+ res <- parse parseKeywordRecord $ flattenKeywords ["OBSGEO-X= -5466045.256954942 / [m]"]+ res @?= ("OBSGEO-X", Float (-5466045.256954942))++ it "should handle extension" $ do+ res <- parse parseKeywordRecord $ flattenKeywords ["XTENSION= 'IMAGE '"]+ res @?= ("XTENSION", String "IMAGE")+++fullRecordLine :: Test ()+fullRecordLine = describe "parseRecordLine" $ do+ it "should parse a normal line" $ do+ res <- parse parseRecordLine "NAXIS1 = 100 / [pix] END"+ res @?= Just ("NAXIS1", Integer 100)++ it "should parse a comment line" $ do+ res <- parse parseRecordLine "COMMENT ------------------------------ Telescope -------------------------------END"+ res @?= Nothing++ it "should parse a blank line" $ do+ res <- parse parseRecordLine $ flattenKeywords [" "]+ res @?= Nothing+++comments :: Test ()+comments = describe "Full-line comments" $ do+ it "should parse full-line comments" $ do+ res <- parse parseLineComment $ flattenKeywords ["COMMENT --------------------------- VISP Instrument ----------------------------"]+ res @?= Comment "--------------------------- VISP Instrument ----------------------------"++ it "should parse comments with text" $ do+ res <- parse parseLineComment $ flattenKeywords ["COMMENT Keys describing the pointing and operation of the telescope. Including "]+ res @?= Comment " Keys describing the pointing and operation of the telescope. Including "++ it "should parse blank comments" $ do+ res <- parse parseLineComment $ flattenKeywords ["COMMENT "]+ res @?= Comment " "+++continue :: Test ()+continue = describe "Continue Keyword" $ do++ it "should be picked up in parseValue" $ do+ res <- parse parseValue $ flattenKeywords ["'hello&'CONTINUE '!'"]+ res @?= String "hello!"++ it "should combine continue into previous keyword" $ do+ let h = [ "CAL_URL = 'https://docs.dkist.nso.edu/projects/visp/en/v2.0.1/l0_to_l1_visp.ht&'"+ , "CONTINUE 'ml' "+ ]++ m <- parse parseHeader $ flattenKeywords h+ Map.lookup "CAL_URL" m @?= Just (String "https://docs.dkist.nso.edu/projects/visp/en/v2.0.1/l0_to_l1_visp.html")+++headerMap :: Test ()+headerMap = describe "full header" $ do+ it "should parse single header" $ do+ res <- parse parseHeader $ flattenKeywords ["KEY1='value'"]+ Map.size res @?= 1+ Map.lookup "KEY1" res @?= Just (String "value")++ it "should parse multiple headers " $ do+ res <- parse parseHeader $ flattenKeywords ["KEY1='value'", "KEY2= 23"]+ Map.size res @?= 2+ Map.lookup "KEY2" res @?= Just (Integer 23)++ it "should ignore comments" $ do+ res <- parse parseHeader $ flattenKeywords ["KEY1='value' / this is a comment"]+ Map.size res @?= 1+ Map.lookup "KEY1" res @?= Just (String "value")++ it "should handle xtension" $ do+ res <- parse parseHeader $ flattenKeywords ["XTENSION= 'IMAGE '"]+ Map.size res @?= 1+ Map.lookup "XTENSION" res @?= Just (String "IMAGE")+++requiredHeaders :: Test ()+requiredHeaders = describe "required headers" $ do++ -- res <- parse parseSimple $ keywords ["SIMPLE= T"]+ -- res @?= Primary++ it "should parse bitpix" $ do+ res <- parse parseBitPix $ flattenKeywords ["BITPIX = 16"]+ res @?= SixteenBitInt++ it "should parse NAxes" $ do+ res <- parse parseNaxes $ flattenKeywords ["NAXIS = 3", "NAXIS1=1", "NAXIS2=2", "NAXIS3=3"]+ res @?= [1,2,3]++ it "should parse size" $ do+ res <- parse parseDimensions $ flattenKeywords ["BITPIX = -32", "NAXIS=2", "NAXIS1=10", "NAXIS2=20"]+ (res ^. bitpix) @?= ThirtyTwoBitFloat+ (res ^. axes) @?= [10,20]++ it "should include required headers in the keywords" $ do+ let fakeData = "1234" -- Related to NAXIS!+ h <- parse parsePrimary $ flattenKeywords ["SIMPLE = T", "BITPIX = 8", "NAXIS=2", "NAXIS1=2", "NAXIS2=2", "TEST='hi'"] <> fakeData+ (h ^. extension) @?= Primary+ Map.size (h ^. header . keywords) @?= 5+ Fits.lookup "NAXIS" (h ^. header) @?= Just (Integer 2) ++ it "should parse full extension" $ do+ h <- parse parseBinTable $ flattenKeywords ["XTENSION= 'BINTABLE'", "BITPIX = -32", "NAXIS=0", "PCOUNT=0", "GCOUNT=1"]+ (h ^. extension) @?= BinTable 0 ""++dataArray :: Test ()+dataArray = describe "data array" $ do+ it "should grab correct data array" $ do+ let fakeData = "1234" -- Related to NAXIS!+ h <- parse parsePrimary $ flattenKeywords ["SIMPLE = T", "BITPIX = 8", "NAXIS=2", "NAXIS1=2", "NAXIS2=2", "TEST='hi'"] <> " " <> fakeData+ (h ^. mainData) @?= fakeData+++sampleSpiral :: Test ()+sampleSpiral =+ describe "Spiral Sample FITS Parse" $ do+ it "should parse" $ do+ let fileSizeOnDisk = 1545444+ bs <- BS.readFile "./fits_files/Spiral_2_30_0_300_10_0_NoGrad.fits"+ hdu <- eitherFail $ readPrimaryHDU bs+ -- hdu.header.size.bitpix @?= ThirtyTwoBitFloat+ -- hdu.header.size.naxes @?= NAxes [621, 621]+ --+ Fits.lookup "NAXIS" (hdu ^. header) @?= Just (Integer 2) ++ let payloadSize = BS.length (hdu ^. mainData)++ -- Make sure we took the right number of bytes out of the file+ payloadSize @?= fileSizeOnDisk - hduBlockSize+++sampleNSOHeaders :: Test ()+sampleNSOHeaders = do+ describe "NSO Stripped Headers" $ do++ it "should parse comment block" $ do+ let h = [ "DATASUM = '550335088' / data unit checksum updated 2023-04-22T04:10:59 " + , " "+ , "COMMENT ------------------------------ Telescope -------------------------------"+ , "COMMENT Keys describing the pointing and operation of the telescope. Including "+ , "COMMENT the FITS WCS keys describing the world coordinates of the array. "+ ]+ m <- parse parseHeader $ flattenKeywords h+ Map.size m @?= 1+++ describe "sample header file" $ do+ bs <- liftIO $ BS.readFile "./fits_files/nso_dkist_headers.txt"+ let ts = filter (not . ignore) $ T.lines $ TE.decodeUtf8 bs++ it "should parse all keywords individually" $ do+ forM_ ( zip [1..] ts ) $ \(n, t) -> do+ _ <- parse parseRecordLine $ flattenKeywords [TE.encodeUtf8 t]+ pure ()++ it "should parse NAxes correctly" $ do+ (sz, _) <- parse parseBinTableKeywords $ mconcat $ C8.lines bs+ (sz ^. axes) @?= [32, 998]++ where+ ignore t = T.isPrefixOf "CONTINUE" t || T.isPrefixOf "END" t++++sampleNSO :: Test ()+sampleNSO = do+ describe "NSO Sample FITS Parse" $ do+ bs <- liftIO $ BS.readFile "./fits_files/nso_dkist.fits"+ it "should parse empty primary header" $ do+ h0 <- eitherFail $ readPrimaryHDU bs+ -- first header doesn't have any data+ BS.length (h0 ^. mainData) @?= 0++ it "should parse both HDUs" $ do+ hdus <- eitherFail $ readHDUs bs+ length hdus @?= 2+ [_, h2] <- pure hdus+ Fits.lookup "INSTRUME" (h2 ^. header) @?= Just (String "VISP")+ Fits.lookup "NAXIS" (h2 ^. header) @?= Just (Integer 2)++ let sizeOnDisk = 161280+ countedHeaderBlocks = 11 -- this was manually counted... until end of all headers+ payloadLength = BS.length (h2 ^. mainData)+ headerLength = countedHeaderBlocks * hduBlockSize+ heapLength = pCount (h2 ^. extension)++ -- Payload size is as expected+ payloadLength @?= 32 * 998 * fromIntegral (bitPixToByteSize EightBitInt)+ pCount (h2 ^. extension) @?= 95968+ assertBool "The end of the heap has some null data" $ do+ C8.all (/='\0') $ C8.take 100 $ C8.drop (headerLength + payloadLength + heapLength - 100) bs+ assertBool "The remainder of the file contains real data" $ do+ C8.all (=='\0') $ C8.drop (headerLength + payloadLength + heapLength) bs+ where+ pCount :: Extension -> Int+ pCount (BinTable p h) = p+ pCount _ = 0++-- Test monad with describe/it+newtype Test a = Test {runTest :: WriterT [TestTree] IO a}+ deriving (Functor, Applicative, Monad, MonadIO, MonadWriter [TestTree])++runTests :: TestName -> Test () -> IO TestTree+runTests n (Test t) = do+ tests <- execWriterT t :: IO [TestTree]+ return $ testGroup n tests+++describe :: TestName -> Test () -> Test ()+describe n t = do+ ts <- liftIO $ runTests n t+ tell [ts]++it :: TestName -> IO () -> Test ()+it n a = do+ tell [testCase n a]+++testMain :: IO TestTree -> IO ()+testMain mtt = do+ tt <- mtt+ defaultMain tt+++