From 4a098e64a76cc71116bc1fa14d7f794123ae5c14 Mon Sep 17 00:00:00 2001 From: Francis Russell Date: Sat, 5 Jan 2013 13:16:43 +0000 Subject: [PATCH] Initial untested implementation of spectral reasoning. --- OFC/Spectral.hs | 86 +++++++++++++++++++++++++++++++++++++++++++------ 1 file changed, 77 insertions(+), 9 deletions(-) diff --git a/OFC/Spectral.hs b/OFC/Spectral.hs index 9c07fd4..972d43a 100644 --- a/OFC/Spectral.hs +++ b/OFC/Spectral.hs @@ -5,7 +5,8 @@ module OFC.Spectral , frequencyMultiple , Band(..) , negateBand - , isEmptyBand + , nullBand + , isNullBand , containsFrequency , containsBand , Spectrum @@ -68,6 +69,22 @@ constantFrequency = Frequency 0 frequencyMultiple :: Integer -> Frequency frequencyMultiple = flip Frequency $ 0 +multiplyFrequency :: Frequency -> Integer -> Frequency +multiplyFrequency (Frequency m o) i = Frequency (m*i) (o*i) + +divideFrequency :: Frequency -> Frequency -> Integer +divideFrequency f1 f2 = fromIntegral $ sign * absDiv f1 f2 + where + absDiv a b = length . takeWhile (>= 0) . tail . iterate (subtract $ abs b) $ abs a + sign = if (f1 >= 0 && f2 >= 0) || (f1 < 0 && f2 < 0) + then 1 else (-1) + +halveFrequency :: Frequency -> Frequency +halveFrequency freq@(Frequency m o) = + if m `mod` 2 == 0 + then Frequency (m `div` 2) (o `div` 2) + else error $ "Cannot halve frequency " ++ show freq + data Band = Band Frequency Frequency deriving (Eq, Show) @@ -78,8 +95,11 @@ instance Ord Band where negateBand :: Band -> Band negateBand (Band low high) = Band (-high + 1) (-low + 1) -isEmptyBand :: Band -> Bool -isEmptyBand (Band low high) = low == high +nullBand :: Band +nullBand = Band 0 0 + +isNullBand :: Band -> Bool +isNullBand (Band low high) = low == high restrictBandToPositive :: Band -> Band restrictBandToPositive (Band low high) = Band (max low 0) (max high 0) @@ -95,14 +115,19 @@ overlaps :: Band -> Band -> Bool overlaps b1@(Band low1 _) b2@(Band low2 _) = containsFrequency b1 low2 || containsFrequency b2 low1 +bandIntersection :: Band -> Band -> Band +bandIntersection b1@(Band l1 h1) b2@(Band l2 h2) = if overlaps b1 b2 + then Band (max l1 l2) (min h1 h2) + else nullBand + data Spectrum = Spectrum { - numSamples :: Frequency, + samplingRate :: Frequency, content :: [(Band, ContentType)] } constructEmptySpectrum :: Integer -> Bool -> Spectrum constructEmptySpectrum multiple oddSamples = Spectrum { - numSamples = (frequencyMultiple $ 2*multiple) + + samplingRate = (frequencyMultiple $ 2*multiple) + (constantFrequency $ if oddSamples then 1 else 0), content = [(Band (constantFrequency 0) (frequencyMultiple multiple), Zero)] } @@ -110,7 +135,7 @@ constructEmptySpectrum multiple oddSamples = Spectrum { mergeContent :: [(Band, ContentType)] -> [(Band, ContentType)] mergeContent content = zip bands $ findContent <$> bands where - bands = removeOverlaps . map fst $ content + bands = splitOverlaps . map fst $ content findContent band = mconcat $ snd <$> filter (overlaps band . fst) content addContent :: [(Band, ContentType)] -> [(Band, ContentType)] -> [(Band, ContentType)] @@ -118,13 +143,48 @@ addContent content1 content2 = mergeContent $ content1 ++ content2 multiplyBands :: Band -> Band -> [Band] multiplyBands (Band low1 high1) (Band low2 high2) = - sort . nub $ filter (not . isEmptyBand) (restrictBandToPositive <$> allResultBands) + sort . nub $ filter (not . isNullBand) (restrictBandToPositive <$> allResultBands) where addedAndSubtracted = [ Band (low1 + low2) (high1 + high2 - 1) , Band (low1 - high2 + 1) (high1 - low2) ] allResultBands = addedAndSubtracted ++ (negateBand <$> addedAndSubtracted) +wrapContent :: Frequency -> (Band, ContentType) -> [(Band, ContentType)] +wrapContent samplingRate (Band low high, content) = + mergeContent $ undamaged:(damagedPositive ++ damagedNegative) + where + nyquist = halveFrequency samplingRate + numFreqs = nyquist + 1 + negNyquist = - nyquist + undamaged = (bandIntersection (Band low high) (Band negNyquist numFreqs), content) + damagedPositive = damagedBand <$> wrapBand samplingRate (Band (max numFreqs low) (max numFreqs high)) + damagedNegative = damagedBand <$> wrapBand samplingRate (Band (min negNyquist low) (min negNyquist high)) + damagedBand band = (band, Damaged) + +wrapBand :: Frequency -> Band -> [Band] +wrapBand samplingRate (Band low high) = + if samplingRate <= 0 + then error "wrapBand takes a sampling rate greater than 0" + else mergeOverlaps wrappedBands + where + wrappedBands = (wrapBand' positiveBand) ++ (negateBand <$> (wrapBand' (negateBand negativeBand))) + positiveBand = Band (max 0 low) (max 0 high) + negativeBand = Band (min 0 low) (min 0 high) + wrapBand' (Band low high) = unwrapped : wrapped : laterBands + where + laterBands = if localHigh == newHigh + then [] + else wrapBand' (Band 0 (newHigh - samplingRate)) + nyquist = halveFrequency samplingRate + numFreqs = nyquist + 1 + offsetFrequency = samplingRate `multiplyFrequency` (low `divideFrequency` samplingRate) + newLow = low - offsetFrequency + newHigh = high - offsetFrequency + localHigh = min newHigh samplingRate + unwrapped = Band (min newLow numFreqs) (min newHigh numFreqs) + wrapped = Band (min (samplingRate - localHigh + 1) numFreqs) (min (samplingRate - newLow + 1) numFreqs) + multiplyContent :: [(Band, ContentType)] -> [(Band, ContentType)] -> [(Band, ContentType)] multiplyContent bandContent1 bandContent2 = mergeContent $ do (band1, content1) <- bandContent1 @@ -132,8 +192,16 @@ multiplyContent bandContent1 bandContent2 = mergeContent $ do resultBand <- multiplyBands band1 band2 return (resultBand, mappend content1 content2) -removeOverlaps :: [Band] -> [Band] -removeOverlaps = sort . foldl' addBand [] +mergeOverlaps :: [Band] -> [Band] +mergeOverlaps = mergeContiguous . sort . splitOverlaps + where + mergeContiguous ((Band l1 h1):(Band l2 h2):bs) = if h1==l2 + then (Band l1 h2):(mergeContiguous bs) + else (Band l1 h2):(Band l2 h2):(mergeContiguous bs) + mergeContiguous bs = bs + +splitOverlaps :: [Band] -> [Band] +splitOverlaps = sort . foldl' addBand [] where addBand bands band = nonoverlapping ++ (createBands edges) where -- 2.47.3