(* Ft.m    copyright 1989  by Victor W. Sparrow *)

(* Permission is granted to copy the file provided
	that the copyright notice is left intact and the
	copies are not distributed commercially. *)

(* The author wishes to thank Roman Maeder, George
	Swenson, and Mike White, who helped in the
	development of this package.  Dr. Maeder for his
	help in the treatment of complex numbers in pattern
	matching, and Drs. Swenson and White for their
	encouragement and suggestions. *)

BeginPackage["Ft`"]

(* Version 9:  added Fta[], Ftb[], iFta[], iFtb[],
   conversion functions for alternate transform
   definitions. 04/10/89 *)
(* Version 8:  fixed delay with optional scaling for
   forward transforms. Also some other rearrangements.
   04/04/89 *)
(* Version 7:  fixed differentiation transforms
   and some other things on 03/27/89 *)
(* Version 6: more inverse transforms and many
   improvements to forward ones 03/23/89 *)
(* Version 5: started adding inverse transforms
   on 03/22/89 *)
(* Version 4: more usage statements, and additional
   improvements to many functions. 03/20/89 *)
(* Version 3: more forward transforms 03/03/89 *)
(* In version 2 Vic Sparrow, the author of this
   package, had a little help on complex numbers from
   Roman Maeder. 02/28/89 *)
   
Ft::usage = "Ft[expr,t,omega] gives the symbolic
 Fourier transform of expr in t space to omega space."
 
iFt::usage = "iFt[expr,t,omega] gives the symbolic
 inverse Fourier transform of expr in omega space to
 t space."
 
Fta::usage = "Fta[expr,t,omega] gives the symbolic
 Fourier transform of expr in t space to omega space.
 In this transform the 2 Pi appears in the exponent
 of the exponential function." 

iFta::usage = "iFta[expr,t,omega] gives the symbolic
 inverse Fourier transform of expr in omega space to
 t space.  In this transform the 2 Pi appears in the
 exponent of the exponential function." 

Ftb::usage = "Ftb[expr,t,omega] gives the symbolic
 Fourier transform of expr in t space to omega space.
 This transform includes Sqrt[2 Pi] as a scaling
 factor"
 
iFtb::usage = "iFtb[expr,t,omega] gives the symbolic
 inverse Fourier transform of expr in omega space to
 t space. This transform includes Sqrt[2 Pi] as a
 scaling factor."
 
convolution::usage = "convolution[a,b] is the convolution
 integral of a and b."
 
delta::usage = "delta[t] is the dirac delta function
 of t."

u::usage = "u[t] is the unit step function of t,
 u[t]=0 for t<0 and u[t]=1 for t>0."

rect::usage = "rect[t] is the rectagular function of t,
 rect[t]=1 for Abs[t]<1/2 and rect[t]=0 otherwise."

sinc::usage = "sinc[t] is the function Sin[t]/t."

triangle::usage = "triangle[t] is the triangular
 function of t, triangle[t]=1-Abs[t] for 
 Abs[t]<1 and triangle[t]=0 otherwise."  
  
Begin["`private`"]

(* Fourier Transform Section: *)  
(* constants *)
Ft[c_,t_,om_]:= 2 Pi c delta[om] /; FreeQ[c,t]

(* linearity *)
Ft[a_ + b_,t_,om_]:= Ft[a,t,om] + Ft[b,t,om]

(* pick off constants *)
Ft[c_ a_,t_,om_]:= c Ft[a,t,om] /; FreeQ[c,t]

(* frequency translation *)
Ft[a_ Exp[c_Complex omn_. t_],t_,om_]:=
 Ft[a,t,om-Im[c]omn] /;
	(FreeQ[omn,t] && Re[c]==0)

(* amplitude modulation *)
Ft[a_ Cos[omn_. t_],t_,om_]:= Ft[a,t,om+omn]/2 +
	Ft[a,t,om-omn]/2 /; FreeQ[omn,t]

(* time differentiation *)
Unprotect[Derivative]
Derivative/: Ft[Derivative[n_Integer?Positive][a_][t_],t_,om_]:=
 (I om)^n Ft[a[t],t,om] 
Protect[Derivative]
 
(* time convolution *)
Ft[convolution[a_,b_],t_,om_]:=
	Ft[a,t,om] Ft[b,t,om]

(* complex conjugate *)
Ft[Conjugate[a_],t_,om_]:= Conjugate[Ft[a,t,-om]]

(* sign function *)
Ft[Sign[t_],t_,om_]:= 2 / (I om)

(* delta function *)
Ft[delta[t_],t_,om_]:= 1

(* step function *)
Ft[u[t_],t_,om_]:= Pi delta[om] - I/om

(* complex exponential *)
Ft[Exp[c_Complex omn_. t_],t_,om_]:=
 2 Pi delta[om - Im[c]omn] /;
	(FreeQ[omn,t] && Re[c]==0)

(* cosine *)
Ft[Cos[omn_. t_],t_,om_]:= Pi(delta[om -omn] +
	delta[om+omn]) /; FreeQ[omn,t]

(* sine *)
Ft[Sin[omn_. t_],t_,om_]:= -I Pi(delta[om - omn] -
	delta[om+omn]) /; FreeQ[omn,t]

(* step exponential *)
Ft[Exp[c_?NumberQ con_. t_] u[t_],t_,om_]:=
 1/(-c con + I om) /;
	(FreeQ[con,t] && c<0 && Im[c]==0)

(* step exponential 2 *)
Ft[t_ Exp[c_?NumberQ con_. t_] u[t_],t_,om_]:=
 1/(-c con + I om)^2 /;
	(FreeQ[con,t] && c<0 && Im[c]==0)

(* double sided exponential *)
Ft[Exp[c_?NumberQ con_. Abs[t_]],t_,om_]:=
 -2 c con/((- c con)^2 + om^2) /;
	(FreeQ[con,t] && c<0 && Im[c]==0)

(* Gaussian:  c=-1/(2 sigma^2)  sigma=Sqrt[-1/2c]
	assume c<0 and sigma>0 *)
Ft[Exp[c_?NumberQ ccon_. t_^2],t_,om_]:=
 Sqrt[- Pi/(c ccon)] Exp[om^2/( 4 c ccon)] /;
	(FreeQ[ccon,t] && c<0 && Im[c]==0)

(* Sign in freq. domain. *)
Ft[Power[t_,-1],t_,om_]:= -I Pi Sign[om] /;
	FreeQ[c,t]

(* Rectangular fn. *)
Ft[rect[t_],t_,om_]:= sinc[om/2]

(* sinc fn. *)
Ft[sinc[t_],t_,om_]:= Pi rect[om/2]

(* Triangle fn. *)
Ft[triangle[t_],t_,om_]:=(sinc[om/2])^2

(* Very general rules should be last. *)
(* frequency differentiation *)
Ft[t_^n_. a_,t_,om_]:=
 I^n Derivative[n][ Ft[a,t,om] ][om] /;
 FreeQ[n,t] && n>0

(* multiplication by delta function *)
Ft[a_ delta[t_ + tnot_.],t_,om_]:=
 Limit[a Exp[-I om t],t->-tnot] /; FreeQ[tnot,t]
 
(* scaling *)
Ft[a_[alph_ t_],t_,om_]:= (1/Abs[alph]
 Ft[a[t],t,om/alph] /; FreeQ[alph,t] )

(* delay with optional scaling *)
Ft[a_[alph_. t_ + tnot_],t_,om_]:= (1/Abs[alph]
	Ft[a[t],t,om/alph] Exp[I tnot om / alph] /;
 FreeQ[alph,t] && FreeQ[tnot,t])
	 
(* frequency convolution for 2 *)
Ft[a_[t_] b_[t_],t_,om_]:=
	convolution[Ft[a[t],t,om],Ft[b[t],t,om]]/(2 Pi)
	
		
(* Inverse Fourier Transform Section: *)
(* constants *)
iFt[c_,t_,om_]:= c delta[t] /; FreeQ[c,om]

(* linearity *)
iFt[a_ + b_,t_,om_]:= iFt[a,t,om] + iFt[b,t,om]

(* pick off constants *)
iFt[c_ a_,t_,om_]:= c iFt[a,t,om] /; FreeQ[c,om]

(* one sided exponential *)
iFt[1/(c_ + I om_),t_,om_]:= Exp[-c t] u[t] /;
 FreeQ[c,om]

(* t and one sided exponential *)
iFt[1/(c_ + I om_)^2,t_,om_]:= t Exp[-c t] u[t] /;
 FreeQ[c,om]

(* double sided exponential *)
iFt[1/(c_ + om_^2),t_,om_]:= ( Exp[-Sqrt[c] Abs[t]]
 /(2 Sqrt[c]) /; FreeQ[c,om] )
 
(* Gaussian *)
(* c=-sigma^2/2)  sigma=Sqrt[-2c]
	assume c<0 and sigma>0 *)
iFt[Exp[c_?NumberQ ccon_. om_^2],t_,om_]:=
 Exp[t^2/(4 c ccon)]/(2 Sqrt[-Pi c ccon]) /;
	(FreeQ[ccon,om] && c<0 && Im[c]==0)

(* Sign[] in time *)
iFt[1/om_,t_,om_]:= I Sign[t]/2

(* Sign[] for frequency *)
iFt[Sign[om_],t_,om_]:= I/(Pi t)

(* delta function *)
iFt[delta[om_],t_,om_]:= 1/( 2 Pi)

(* sine: this could be improved *)
iFt[delta[om_-omnot_]-delta[om_+omnot_],t_,om_]:=
 I Sin[omnot t]/Pi /; FreeQ[omnot,om]
 
(* cosine: this could be improved *)
iFt[delta[om_-omnot_]+delta[om_+omnot_],t_,om_]:=
 Cos[omnot t]/Pi /; FreeQ[omnot,om]

(* sinc function *)
iFt[sinc[om_],t_,om_]:= rect[t/2]/2 
 
(* rect function *)
iFt[rect[om_],t_,om_]:= sinc[t/2]/(2 Pi) 

(* sinc squared *)
iFt[sinc[om_ tau_.]^2,t_,om_]:=
 triangle[t/(2 tau)]/(2 tau) /; FreeQ[tau,om]
 
(* time differentiation *)
iFt[om_^n_. a_,t_,om_]:=
 (-I)^n Derivative[n][ iFt[a,t,om] ][t] /;
 FreeQ[n,om] && n>0
 
(* frequency differentiation *)
Unprotect[Derivative]
Derivative/: iFt[Derivative[n_Integer?Positive][a_][om_],t_,om_]:=
 (t/I)^n iFt[a[om],t,om]
Protect[Derivative]
  
(* general rules last *)
(* multiplication by delta function *)
iFt[a_ delta[om_ + omnot_.],t_,om_]:=
 Limit[a Exp[I om t]/(2 Pi),om->-omnot] /; FreeQ[omnot,om]
 
(* scaling *)
iFt[a_[alph_ om_],t_,om_]:=
 ( iFt[a[om],t/alph,om]/Abs[alph] /; FreeQ[alph,om] )
 
(* frequency translation *)
iFt[a_[om_ + omnot_],t_,om_]:=
 iFt[a[om],t,om] Exp[-I omnot t]/;
 FreeQ[omnot,om]  

(* time translation *)
iFt[a_[om_] Exp[c_Complex tnot_. om_],t_,om_]:=
 iFt[a[om],t + Im[c] tnot, om] /;
 FreeQ[tnot,om] && Re[c]==0
 
(* frequency multiplication *)
iFt[a_[om_] b_[om_],t_,om_]:=
 convolution[iFt[a[om],t,om],iFt[b[om],t,om]]
  
(* Alternative definition transform section. *)
(* Alternative a:  2Pi only in exponent of expoentials. *)
Fta[a_,t_,om_]:=Ft[a,t,2 Pi om]
iFta[a_,t_,om_]:=iFt[2 Pi a, 2 Pi t, om]

(* Alternative b:  Sqrt[2 Pi] appears in both forward and
  inverse trasforms. *)
Ftb[a_,t_,om_]:=Ft[a,t,om]/Sqrt[2 Pi]
iFtb[a_,t_,om_]:=Sqrt[2 Pi] iFt[a,t,om]

End[]
EndPackage[]
