1//
2//  PRDFTHighPass.m
3//  PRICE
4//
5//  Created by Riccardo Mottola on Fri Oct 24 2003.
6//  Copyright (c) 2003-2009 Carduus. All rights reserved.
7//
8// This application is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version.
9// This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
10
11#import "PRDFTHighPass.h"
12#import "PRDFTFilter.h"
13#import "FFT.h"
14#include "math.h"
15
16extern double *cosinus;
17extern double *sinus;
18extern unsigned int *bitRevIndex;
19
20@implementation PRDFTHighPass
21
22- (PRImage *)filterImage:(PRImage *)image with:(NSArray *)parameters progressPanel:(PRCProgress *)progressPanel
23{
24    BOOL autoRange;
25    float bandPassFreq;
26    float bandStopFreq;
27
28    autoRange = [[parameters objectAtIndex:0] boolValue];
29    bandPassFreq = [[parameters objectAtIndex:1] floatValue];
30    bandStopFreq = [[parameters objectAtIndex:2] floatValue];
31    return [self transformImage:image :autoRange :bandPassFreq :bandStopFreq :progressPanel];
32}
33
34- (NSString *)actionName
35{
36    return @"DFT High Pass";
37}
38
39
40- (PRImage *)transformImage:(PRImage *)srcImage :(BOOL)autoRange :(float) BPFreq :(float) BSFreq :(PRCProgress *)prPanel
41{
42    PRDFTFilter  *filter;
43    double       **Fa, **Fb;
44    unsigned int i, j;
45    unsigned int w, h;
46    unsigned int bitNum;
47    unsigned int num;
48    unsigned int halfNum;
49    float        radius;
50    PRImage      *destImage;
51
52    progressSteps = 0;
53    totalProgressSteps = 5;
54    progPanel = prPanel;
55
56
57    /* get source image representation and associated information */
58    if (progPanel != nil)
59    {
60        [self setActivity:@"Get image size"];
61        [self advanceProgress];
62    }
63    w = [srcImage width];
64    h = [srcImage height];
65
66    /* we set the image square to the nearest power of two on the longer side */
67    if (w > h)
68        bitNum = binaryLog(w);
69    else
70        bitNum = binaryLog(h);
71    num = binpow(bitNum);
72
73    /* calculate the half comprising the zero */
74    halfNum = (num >> 1);
75    NSLog(@"halfNum %d", halfNum);
76
77    /* allocate filter  matrix */
78    Fa = (double**)calloc(num, sizeof(double*));
79    for (i = 0; i < num; i++)
80    {
81        Fa[i] = (double*)calloc(num, sizeof(double));
82        if (!Fa[i])
83        {
84            printf("out of memory allocating float matrix?\n");
85            return srcImage;
86        }
87        memset(Fa[i], '\000', num);
88    }
89    Fb = (double**)calloc(num, sizeof(double*));
90    for (i = 0; i < num; i++)
91    {
92        Fb[i] = (double*)calloc(num, sizeof(double));
93        if (!Fb[i])
94        {
95            printf("out of memory allocating float matrix?\n");
96            return srcImage;
97        }
98        memset(Fb[i], '\000', num);
99    }
100
101    /* calculate the filter */
102    NSLog(@"BP %f, BS %f", BPFreq, BSFreq);
103    /* top left quadrant */
104    for (i = 0; i < halfNum; i++)
105        for (j = 0; j < halfNum; j++)
106        {
107            radius = sqrt(i*i + j*j) / halfNum;
108            if (radius > BPFreq)
109            {
110                Fa[i][j] = 1;
111            } else if (radius > BSFreq)
112            {
113                Fa[i][j] = 0.5*(cos(M_PI*(BPFreq-radius)/(BPFreq-BSFreq))+1);
114            }
115        }
116
117    /* top right quadrant */
118    for (i = 0; i < halfNum; i++)
119        for (j = 0; j < halfNum; j++)
120            Fa[i][num - j - 1] = Fa[i][j];
121    /* bottom left quadrant */
122    for (i = 0; i < halfNum; i++)
123        for (j = 0; j < halfNum; j++)
124            Fa[num - i - 1][j] = Fa[i][j];
125    /* bottom right quadrant */
126    for (i = 0; i < halfNum; i++)
127        for (j = 0; j < halfNum; j++)
128            Fa[num - i - 1][num - j - 1] = Fa[i][j];
129
130    if (progPanel != nil)
131    {
132        [self setActivity:@"run the filter"];
133        [self advanceProgress];
134    }
135
136    /* instantiate the filter */
137    filter = [[PRDFTFilter alloc] init];
138
139    /* run it */
140    destImage = [[filter filterImage :srcImage :Fa :Fb :num :autoRange :NULL] retain];
141
142    /* release filter */
143    [filter release];
144
145    /* now we free our filter matrix */
146    for (i = 0; i < num; i++)
147        free(Fa[i]);
148    free(Fa);
149    for (i = 0; i < num; i++)
150        free(Fb[i]);
151    free(Fb);
152
153    if (progPanel != nil)
154    {
155        [self setActivity:@"Done"];
156        [self showProgress];
157    }
158    [destImage autorelease];
159    return destImage;
160}
161
162
163@end
164