MagickCore 7.1.1
Convert, Edit, Or Compose Bitmap Images
Loading...
Searching...
No Matches
gem.c
1/*
2%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3% %
4% %
5% %
6% GGGG EEEEE M M %
7% G E MM MM %
8% G GG EEE M M M %
9% G G E M M %
10% GGGG EEEEE M M %
11% %
12% %
13% Graphic Gems - Graphic Support Methods %
14% %
15% Software Design %
16% Cristy %
17% August 1996 %
18% %
19% %
20% Copyright @ 1999 ImageMagick Studio LLC, a non-profit organization %
21% dedicated to making software imaging solutions freely available. %
22% %
23% You may not use this file except in compliance with the License. You may %
24% obtain a copy of the License at %
25% %
26% https://imagemagick.org/script/license.php %
27% %
28% Unless required by applicable law or agreed to in writing, software %
29% distributed under the License is distributed on an "AS IS" BASIS, %
30% WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. %
31% See the License for the specific language governing permissions and %
32% limitations under the License. %
33% %
34%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
35%
36%
37%
38*/
39
40/*
41 Include declarations.
42*/
43#include "MagickCore/studio.h"
44#include "MagickCore/color-private.h"
45#include "MagickCore/draw.h"
46#include "MagickCore/gem.h"
47#include "MagickCore/gem-private.h"
48#include "MagickCore/image.h"
49#include "MagickCore/image-private.h"
50#include "MagickCore/log.h"
51#include "MagickCore/memory_.h"
52#include "MagickCore/pixel-accessor.h"
53#include "MagickCore/quantum.h"
54#include "MagickCore/quantum-private.h"
55#include "MagickCore/random_.h"
56#include "MagickCore/resize.h"
57#include "MagickCore/transform.h"
58#include "MagickCore/signature-private.h"
59
60/*
61%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
62% %
63% %
64% %
65% E x p a n d A f f i n e %
66% %
67% %
68% %
69%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
70%
71% ExpandAffine() computes the affine's expansion factor, i.e. the square root
72% of the factor by which the affine transform affects area. In an affine
73% transform composed of scaling, rotation, shearing, and translation, returns
74% the amount of scaling.
75%
76% The format of the ExpandAffine method is:
77%
78% double ExpandAffine(const AffineMatrix *affine)
79%
80% A description of each parameter follows:
81%
82% o expansion: ExpandAffine returns the affine's expansion factor.
83%
84% o affine: A pointer the affine transform of type AffineMatrix.
85%
86*/
87MagickExport double ExpandAffine(const AffineMatrix *affine)
88{
89 assert(affine != (const AffineMatrix *) NULL);
90 return(sqrt(fabs(affine->sx*affine->sy-affine->rx*affine->ry)));
91}
92
93/*
94%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
95% %
96% %
97% %
98% G e n e r a t e D i f f e r e n t i a l N o i s e %
99% %
100% %
101% %
102%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
103%
104% GenerateDifferentialNoise() generates differential noise.
105%
106% The format of the GenerateDifferentialNoise method is:
107%
108% double GenerateDifferentialNoise(RandomInfo *random_info,
109% const Quantum pixel,const NoiseType noise_type,const double attenuate)
110%
111% A description of each parameter follows:
112%
113% o random_info: the random info.
114%
115% o pixel: noise is relative to this pixel value.
116%
117% o noise_type: the type of noise.
118%
119% o attenuate: attenuate the noise.
120%
121*/
122MagickPrivate double GenerateDifferentialNoise(RandomInfo *random_info,
123 const Quantum pixel,const NoiseType noise_type,const double attenuate)
124{
125#define SigmaUniform (attenuate*0.015625)
126#define SigmaGaussian (attenuate*0.015625)
127#define SigmaImpulse (attenuate*0.1)
128#define SigmaLaplacian (attenuate*0.0390625)
129#define SigmaMultiplicativeGaussian (attenuate*0.5)
130#define SigmaPoisson (attenuate*12.5)
131#define SigmaRandom (attenuate)
132#define TauGaussian (attenuate*0.078125)
133
134 double
135 alpha,
136 beta,
137 noise,
138 sigma;
139
140 alpha=GetPseudoRandomValue(random_info);
141 switch (noise_type)
142 {
143 case UniformNoise:
144 default:
145 {
146 noise=(double) pixel+(double) QuantumRange*SigmaUniform*(alpha-0.5);
147 break;
148 }
149 case GaussianNoise:
150 {
151 double
152 gamma,
153 tau;
154
155 if (fabs(alpha) < MagickEpsilon)
156 alpha=1.0;
157 beta=GetPseudoRandomValue(random_info);
158 gamma=sqrt(-2.0*log(alpha));
159 sigma=gamma*cos((double) (2.0*MagickPI*beta));
160 tau=gamma*sin((double) (2.0*MagickPI*beta));
161 noise=(double) pixel+sqrt((double) pixel)*SigmaGaussian*sigma+
162 (double) QuantumRange*TauGaussian*tau;
163 break;
164 }
165 case ImpulseNoise:
166 {
167 if (alpha < (SigmaImpulse/2.0))
168 noise=0.0;
169 else
170 if (alpha >= (1.0-(SigmaImpulse/2.0)))
171 noise=(double) QuantumRange;
172 else
173 noise=(double) pixel;
174 break;
175 }
176 case LaplacianNoise:
177 {
178 if (alpha <= 0.5)
179 {
180 if (alpha <= MagickEpsilon)
181 noise=(double) (pixel-QuantumRange);
182 else
183 noise=(double) pixel+(double) QuantumRange*SigmaLaplacian*
184 log(2.0*alpha)+0.5;
185 break;
186 }
187 beta=1.0-alpha;
188 if (beta <= (0.5*MagickEpsilon))
189 noise=(double) (pixel+QuantumRange);
190 else
191 noise=(double) pixel-(double) QuantumRange*SigmaLaplacian*
192 log(2.0*beta)+0.5;
193 break;
194 }
195 case MultiplicativeGaussianNoise:
196 {
197 sigma=1.0;
198 if (alpha > MagickEpsilon)
199 sigma=sqrt(-2.0*log(alpha));
200 beta=GetPseudoRandomValue(random_info);
201 noise=(double) pixel+(double) pixel*SigmaMultiplicativeGaussian*sigma*
202 cos((double) (2.0*MagickPI*beta))/2.0;
203 break;
204 }
205 case PoissonNoise:
206 {
207 double
208 poisson;
209
210 ssize_t
211 i;
212
213 poisson=exp(-SigmaPoisson*QuantumScale*(double) pixel);
214 for (i=0; alpha > poisson; i++)
215 {
216 beta=GetPseudoRandomValue(random_info);
217 alpha*=beta;
218 }
219 noise=(double) QuantumRange*i*PerceptibleReciprocal(SigmaPoisson);
220 break;
221 }
222 case RandomNoise:
223 {
224 noise=(double) QuantumRange*SigmaRandom*alpha;
225 break;
226 }
227 }
228 return(noise);
229}
230
231/*
232%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
233% %
234% %
235% %
236% G e t O p t i m a l K e r n e l W i d t h %
237% %
238% %
239% %
240%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
241%
242% GetOptimalKernelWidth() computes the optimal kernel radius for a convolution
243% filter. Start with the minimum value of 3 pixels and walk out until we drop
244% below the threshold of one pixel numerical accuracy.
245%
246% The format of the GetOptimalKernelWidth method is:
247%
248% size_t GetOptimalKernelWidth(const double radius,
249% const double sigma)
250%
251% A description of each parameter follows:
252%
253% o width: GetOptimalKernelWidth returns the optimal width of a
254% convolution kernel.
255%
256% o radius: the radius of the Gaussian, in pixels, not counting the center
257% pixel.
258%
259% o sigma: the standard deviation of the Gaussian, in pixels.
260%
261*/
262MagickPrivate size_t GetOptimalKernelWidth1D(const double radius,
263 const double sigma)
264{
265 double
266 alpha,
267 beta,
268 gamma,
269 normalize,
270 value;
271
272 size_t
273 width;
274
275 ssize_t
276 i,
277 j;
278
279 if (IsEventLogging() != MagickFalse)
280 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"...");
281 if (radius > MagickEpsilon)
282 return((size_t) (2.0*ceil(radius)+1.0));
283 gamma=fabs(sigma);
284 if (gamma <= MagickEpsilon)
285 return(3UL);
286 alpha=PerceptibleReciprocal(2.0*gamma*gamma);
287 beta=(double) PerceptibleReciprocal((double) MagickSQ2PI*gamma);
288 for (width=5; ; )
289 {
290 normalize=0.0;
291 j=(ssize_t) (width-1)/2;
292 for (i=(-j); i <= j; i++)
293 normalize+=exp(-((double) (i*i))*alpha)*beta;
294 value=exp(-((double) (j*j))*alpha)*beta/normalize;
295 if ((value < QuantumScale) || (value < MagickEpsilon))
296 break;
297 width+=2;
298 }
299 return((size_t) (width-2));
300}
301
302MagickPrivate size_t GetOptimalKernelWidth2D(const double radius,
303 const double sigma)
304{
305 double
306 alpha,
307 beta,
308 gamma,
309 normalize,
310 value;
311
312 size_t
313 width;
314
315 ssize_t
316 j,
317 u,
318 v;
319
320 if (IsEventLogging() != MagickFalse)
321 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"...");
322 if (radius > MagickEpsilon)
323 return((size_t) (2.0*ceil(radius)+1.0));
324 gamma=fabs(sigma);
325 if (gamma <= MagickEpsilon)
326 return(3UL);
327 alpha=PerceptibleReciprocal(2.0*gamma*gamma);
328 beta=(double) PerceptibleReciprocal((double) Magick2PI*gamma*gamma);
329 for (width=5; ; )
330 {
331 normalize=0.0;
332 j=(ssize_t) (width-1)/2;
333 for (v=(-j); v <= j; v++)
334 for (u=(-j); u <= j; u++)
335 normalize+=exp(-((double) (u*u+v*v))*alpha)*beta;
336 value=exp(-((double) (j*j))*alpha)*beta/normalize;
337 if ((value < QuantumScale) || (value < MagickEpsilon))
338 break;
339 width+=2;
340 }
341 return((size_t) (width-2));
342}
343
344MagickPrivate size_t GetOptimalKernelWidth(const double radius,
345 const double sigma)
346{
347 return(GetOptimalKernelWidth1D(radius,sigma));
348}