[ido,l1,cc,ch,wa1,wa2,wa3,wa4]=passb5(ido,l1,cc,ch,wa1,wa2,wa3,wa4);

function [ido,l1,cc,ch,wa1,wa2,wa3,wa4]=passb5(ido,l1,cc,ch,wa1,wa2,wa3,wa4);
persistent ci2 ci3 ci4 ci5 cr2 cr3 cr4 cr5 di2 di3 di4 di5 dr2 dr3 dr4 dr5 i k pi ti11 ti12 ti2 ti3 ti4 ti5 tr11 tr12 tr2 tr3 tr4 tr5 ; 

if isempty(ci2), ci2=0; end;
if isempty(ci3), ci3=0; end;
if isempty(ci4), ci4=0; end;
if isempty(ci5), ci5=0; end;
if isempty(cr2), cr2=0; end;
if isempty(cr3), cr3=0; end;
if isempty(cr4), cr4=0; end;
if isempty(cr5), cr5=0; end;
if isempty(di2), di2=0; end;
if isempty(di3), di3=0; end;
if isempty(di4), di4=0; end;
if isempty(di5), di5=0; end;
if isempty(dr2), dr2=0; end;
if isempty(dr3), dr3=0; end;
if isempty(dr4), dr4=0; end;
if isempty(dr5), dr5=0; end;
if isempty(pi), pi=0; end;
if isempty(ti11), ti11=0; end;
if isempty(ti12), ti12=0; end;
if isempty(ti2), ti2=0; end;
if isempty(ti3), ti3=0; end;
if isempty(ti4), ti4=0; end;
if isempty(ti5), ti5=0; end;
if isempty(tr11), tr11=0; end;
if isempty(tr12), tr12=0; end;
if isempty(tr2), tr2=0; end;
if isempty(tr3), tr3=0; end;
if isempty(tr4), tr4=0; end;
if isempty(tr5), tr5=0; end;
if isempty(i), i=0; end;
if isempty(k), k=0; end;
%***BEGIN PROLOGUE  PASSB5
%***SUBSIDIARY
%***PURPOSE  Calculate the fast Fourier transform of subvectors of
%            length five.
%***LIBRARY   SLATEC (FFTPACK)
%***TYPE      SINGLE PRECISION (PASSB5-S)
%***AUTHOR  Swarztrauber, P. N., (NCAR)
%***ROUTINES CALLED  (NONE)
%***REVISION HISTORY  (YYMMDD)
%   790601  DATE WRITTEN
%   830401  Modified to use SLATEC library source file format.
%   860115  Modified by Ron Boisvert to adhere to Fortran 77 by
%           (a) changing dummy array size declarations (1) to (*),
%           (b) changing definition of variables PI, TI11, TI12,
%               TR11, TR12 by using FORTRAN intrinsic functions ATAN
%               and SIN instead of DATA statements.
%   881128  Modified by Dick Valent to meet prologue standards.
%   890831  Modified array declarations.  (WRB)
%   891214  Prologue converted to Version 4.0 format.  (BAB)
%   900402  Added TYPE section.  (WRB)
%***end PROLOGUE  PASSB5
cc_shape=size(cc);cc=reshape([cc(:).',zeros(1,ceil(numel(cc)./prod([ido,5])).*prod([ido,5])-numel(cc))],ido,5,[]);
ch_orig=ch;ch_shape=[ido,l1,5];ch=reshape([ch_orig(1:min(prod(ch_shape),numel(ch_orig))),zeros(1,max(0,prod(ch_shape)-numel(ch_orig)))],ch_shape);
wa1_shape=size(wa1);wa1=reshape(wa1,1,[]);
wa2_shape=size(wa2);wa2=reshape(wa2,1,[]);
wa3_shape=size(wa3);wa3=reshape(wa3,1,[]);
wa4_shape=size(wa4);wa4=reshape(wa4,1,[]);
%***FIRST EXECUTABLE STATEMENT  PASSB5
pi = 4..*atan(1.);
tr11 = sin(.1.*pi);
ti11 = sin(.4.*pi);
tr12 = -sin(.3.*pi);
ti12 = sin(.2.*pi);
if( ido==2 )
for k = 1 : l1;
ti5 = cc(2,2,k) - cc(2,5,k);
ti2 = cc(2,2,k) + cc(2,5,k);
ti4 = cc(2,3,k) - cc(2,4,k);
ti3 = cc(2,3,k) + cc(2,4,k);
tr5 = cc(1,2,k) - cc(1,5,k);
tr2 = cc(1,2,k) + cc(1,5,k);
tr4 = cc(1,3,k) - cc(1,4,k);
tr3 = cc(1,3,k) + cc(1,4,k);
ch(1,k,1) = cc(1,1,k) + tr2 + tr3;
ch(2,k,1) = cc(2,1,k) + ti2 + ti3;
cr2 = cc(1,1,k) + tr11.*tr2 + tr12.*tr3;
ci2 = cc(2,1,k) + tr11.*ti2 + tr12.*ti3;
cr3 = cc(1,1,k) + tr12.*tr2 + tr11.*tr3;
ci3 = cc(2,1,k) + tr12.*ti2 + tr11.*ti3;
cr5 = ti11.*tr5 + ti12.*tr4;
ci5 = ti11.*ti5 + ti12.*ti4;
cr4 = ti12.*tr5 - ti11.*tr4;
ci4 = ti12.*ti5 - ti11.*ti4;
ch(1,k,2) = cr2 - ci5;
ch(1,k,5) = cr2 + ci5;
ch(2,k,2) = ci2 + cr5;
ch(2,k,3) = ci3 + cr4;
ch(1,k,3) = cr3 - ci4;
ch(1,k,4) = cr3 + ci4;
ch(2,k,4) = ci3 - cr4;
ch(2,k,5) = ci2 - cr5;
end; k = fix(l1+1);
cc_shape=zeros(cc_shape);cc_shape(:)=cc(1:numel(cc_shape));cc=cc_shape;
ch_orig(1:prod(ch_shape))=ch;ch=ch_orig;
wa1_shape=zeros(wa1_shape);wa1_shape(:)=wa1(1:numel(wa1_shape));wa1=wa1_shape;
wa2_shape=zeros(wa2_shape);wa2_shape(:)=wa2(1:numel(wa2_shape));wa2=wa2_shape;
wa3_shape=zeros(wa3_shape);wa3_shape(:)=wa3(1:numel(wa3_shape));wa3=wa3_shape;
wa4_shape=zeros(wa4_shape);wa4_shape(:)=wa4(1:numel(wa4_shape));wa4=wa4_shape;
return;
elseif( fix(ido./2)<l1 ) ;
for i = 2 : 2: ido ;
%DIR$ IVDEP
for k = 1 : l1;
ti5 = cc(i,2,k) - cc(i,5,k);
ti2 = cc(i,2,k) + cc(i,5,k);
ti4 = cc(i,3,k) - cc(i,4,k);
ti3 = cc(i,3,k) + cc(i,4,k);
tr5 = cc(i-1,2,k) - cc(i-1,5,k);
tr2 = cc(i-1,2,k) + cc(i-1,5,k);
tr4 = cc(i-1,3,k) - cc(i-1,4,k);
tr3 = cc(i-1,3,k) + cc(i-1,4,k);
ch(i-1,k,1) = cc(i-1,1,k) + tr2 + tr3;
ch(i,k,1) = cc(i,1,k) + ti2 + ti3;
cr2 = cc(i-1,1,k) + tr11.*tr2 + tr12.*tr3;
ci2 = cc(i,1,k) + tr11.*ti2 + tr12.*ti3;
cr3 = cc(i-1,1,k) + tr12.*tr2 + tr11.*tr3;
ci3 = cc(i,1,k) + tr12.*ti2 + tr11.*ti3;
cr5 = ti11.*tr5 + ti12.*tr4;
ci5 = ti11.*ti5 + ti12.*ti4;
cr4 = ti12.*tr5 - ti11.*tr4;
ci4 = ti12.*ti5 - ti11.*ti4;
dr3 = cr3 - ci4;
dr4 = cr3 + ci4;
di3 = ci3 + cr4;
di4 = ci3 - cr4;
dr5 = cr2 + ci5;
dr2 = cr2 - ci5;
di5 = ci2 - cr5;
di2 = ci2 + cr5;
ch(i-1,k,2) = wa1(i-1).*dr2 - wa1(i).*di2;
ch(i,k,2) = wa1(i-1).*di2 + wa1(i).*dr2;
ch(i-1,k,3) = wa2(i-1).*dr3 - wa2(i).*di3;
ch(i,k,3) = wa2(i-1).*di3 + wa2(i).*dr3;
ch(i-1,k,4) = wa3(i-1).*dr4 - wa3(i).*di4;
ch(i,k,4) = wa3(i-1).*di4 + wa3(i).*dr4;
ch(i-1,k,5) = wa4(i-1).*dr5 - wa4(i).*di5;
ch(i,k,5) = wa4(i-1).*di5 + wa4(i).*dr5;
end; k = fix(l1+1);
end; i = fix(ido +1);
else;
for k = 1 : l1;
%DIR$ IVDEP
for i = 2 : 2: ido ;
ti5 = cc(i,2,k) - cc(i,5,k);
ti2 = cc(i,2,k) + cc(i,5,k);
ti4 = cc(i,3,k) - cc(i,4,k);
ti3 = cc(i,3,k) + cc(i,4,k);
tr5 = cc(i-1,2,k) - cc(i-1,5,k);
tr2 = cc(i-1,2,k) + cc(i-1,5,k);
tr4 = cc(i-1,3,k) - cc(i-1,4,k);
tr3 = cc(i-1,3,k) + cc(i-1,4,k);
ch(i-1,k,1) = cc(i-1,1,k) + tr2 + tr3;
ch(i,k,1) = cc(i,1,k) + ti2 + ti3;
cr2 = cc(i-1,1,k) + tr11.*tr2 + tr12.*tr3;
ci2 = cc(i,1,k) + tr11.*ti2 + tr12.*ti3;
cr3 = cc(i-1,1,k) + tr12.*tr2 + tr11.*tr3;
ci3 = cc(i,1,k) + tr12.*ti2 + tr11.*ti3;
cr5 = ti11.*tr5 + ti12.*tr4;
ci5 = ti11.*ti5 + ti12.*ti4;
cr4 = ti12.*tr5 - ti11.*tr4;
ci4 = ti12.*ti5 - ti11.*ti4;
dr3 = cr3 - ci4;
dr4 = cr3 + ci4;
di3 = ci3 + cr4;
di4 = ci3 - cr4;
dr5 = cr2 + ci5;
dr2 = cr2 - ci5;
di5 = ci2 - cr5;
di2 = ci2 + cr5;
ch(i-1,k,2) = wa1(i-1).*dr2 - wa1(i).*di2;
ch(i,k,2) = wa1(i-1).*di2 + wa1(i).*dr2;
ch(i-1,k,3) = wa2(i-1).*dr3 - wa2(i).*di3;
ch(i,k,3) = wa2(i-1).*di3 + wa2(i).*dr3;
ch(i-1,k,4) = wa3(i-1).*dr4 - wa3(i).*di4;
ch(i,k,4) = wa3(i-1).*di4 + wa3(i).*dr4;
ch(i-1,k,5) = wa4(i-1).*dr5 - wa4(i).*di5;
ch(i,k,5) = wa4(i-1).*di5 + wa4(i).*dr5;
end; i = fix(ido +1);
end; k = fix(l1+1);
cc_shape=zeros(cc_shape);cc_shape(:)=cc(1:numel(cc_shape));cc=cc_shape;
ch_orig(1:prod(ch_shape))=ch;ch=ch_orig;
wa1_shape=zeros(wa1_shape);wa1_shape(:)=wa1(1:numel(wa1_shape));wa1=wa1_shape;
wa2_shape=zeros(wa2_shape);wa2_shape(:)=wa2(1:numel(wa2_shape));wa2=wa2_shape;
wa3_shape=zeros(wa3_shape);wa3_shape(:)=wa3(1:numel(wa3_shape));wa3=wa3_shape;
wa4_shape=zeros(wa4_shape);wa4_shape(:)=wa4(1:numel(wa4_shape));wa4=wa4_shape;
return;
end;
cc_shape=zeros(cc_shape);cc_shape(:)=cc(1:numel(cc_shape));cc=cc_shape;
ch_orig(1:prod(ch_shape))=ch;ch=ch_orig;
wa1_shape=zeros(wa1_shape);wa1_shape(:)=wa1(1:numel(wa1_shape));wa1=wa1_shape;
wa2_shape=zeros(wa2_shape);wa2_shape(:)=wa2(1:numel(wa2_shape));wa2=wa2_shape;
wa3_shape=zeros(wa3_shape);wa3_shape(:)=wa3(1:numel(wa3_shape));wa3=wa3_shape;
wa4_shape=zeros(wa4_shape);wa4_shape(:)=wa4(1:numel(wa4_shape));wa4=wa4_shape;
end
%DECK PASSB