1#  This program is free software; you can redistribute it and/or modify
2#  it under the terms of the GNU General Public License as published by
3#  the Free Software Foundation; either version 2 of the License, or
4#  (at your option) any later version.
5#
6#  This program is distributed in the hope that it will be useful,
7#  but WITHOUT ANY WARRANTY; without even the implied warranty of
8#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
9#  GNU General Public License for more details.
10#
11#  A copy of the GNU General Public License is available at
12#  http://www.r-project.org/Licenses/
13
14
15charStable <- function(theta, tau, pm = 0)
16	{
17	# pm is the type parametrization as described by Nolan(2009)
18	# It takes the value 0 or 1
19
20	# const can fixe parameters. It is NULL for no constraint or
21	# a matrix in which case the constraint is theta[const[,1]]=const[,2]
22
23	a <- theta[1]
24	b <- theta[2]
25	g <- theta[3]
26	d <- theta[4]
27	if(pm == 0)
28		{
29		if(a == 1)
30			{
31			if(g == 0)
32				{
33				the_car <- exp(complex(imaginary = d*tau))
34				}
35			else
36				{
37				re_p <- -g * abs(tau)
38				im_p <- d * tau
39				im_p[tau!=0] <- im_p[tau != 0] + re_p[tau != 0]*2/pi*b*sign(tau[tau != 0])*log(g*abs(tau[tau != 0]))
40				the_car <- exp(complex(real = re_p, imaginary = im_p))
41				}
42			}
43		else
44			{
45			if(g == 0)
46				{
47				the_car <- exp(complex(imaginary = d*tau))
48				}
49			else
50				{
51				phi <- tan(pi*a/2)
52				re_p <- -g^a*abs(tau)^a
53				im_p <- d*tau*1i
54				im_p[tau != 0] <- im_p[tau != 0] + re_p*( b*phi*sign(tau[tau != 0])*(abs(g*tau[tau != 0])^(1-a) - 1) )
55				the_car <- exp(complex(real = re_p, imaginary = im_p))
56				}
57			}
58		}
59
60	if(pm == 1)
61		{
62		if(a == 1)
63			{
64			re_p <- -g*abs(tau)
65			im_p <- d*tau
66			im_p[tau!=0] <- im_p[tau != 0]+re_p*(b*2/pi*sign(tau[tau != 0])*log(abs(tau[tau!=0])))
67			the_car <- exp(complex(real = re_p, imaginary = im_p))
68			}
69		else
70			{
71			phi <- tan(pi*a/2)
72			re_p <- -g^a*abs(tau)^a
73			im_p <- re_p*(-phi*b*sign(tau)) + d*tau
74			the_car <- exp(complex(real = re_p, imaginary = im_p))
75			}
76		}
77	return(the_car)
78	}
79
80
81
82