Season's greetings from Phastar's SAS Art competition winners!

*CHRISTMAS ART 2022;

dm log 'clear; output clear';
proc datasets mt = data lib = work kill nolist; run; quit;

*create "page" and set ideal number of snowflakes;
%let xmax = 8.3;
%let ymax = 11.7;
%let it_max = 14;

%macro centres();

data centres;
num=0;
call streaminit(5678);
run;

%let i = 1;
%let it_current = 1;

%do %while(&i <= &it_max.);

*create first centre and radius;
data new;
x_new = rand("Uniform", 0, &xmax.);
y_new = rand("Uniform", 0, &ymax.);
size = rand("Uniform", 0.7, 2.5);
num = 1;
i = &i.;
run;

*reset centres if size puts them over the edge;
data new2;
set new;
if x_new+size>&xmax. then x_new = &xmax. - size;
if x_new-size<0 then x_new = size;
if y_new+size>&ymax. then y_new = &ymax. - size;
if y_new-size<0 then y_new = size;
run;

*calculate edges of section so that none can overlap, with x1 being the bottom left and x2, x3, x4 going clockwise;
data new3;
set new2;
*create small r=Rcos(1/5pi);
r = size * cos(0.2*constant("pi"));
*rename size to be bigR, the circumradius of our pentagons;
rename size=bigR ;
run;

%if &i=1 %then %do;
%let i = %sysevalf(&i. + 1);

data centres;
set new3(rename=(x_new=x_centre y_new=y_centre))
centres;
if num ne 0;
run;

%end;

%else %do;

*merge with existing sections;
data new4;
merge new3 (keep = x_new y_new bigR num
rename=(bigR = bigR_new))
centres (keep = x_centre y_centre bigR num) ;
by num;
*find the maximum distnace these areas are allowed to be - the two radii from eachother;
d1 = bigR_new + bigR;
*find the length between each centre;
d2 = sqrt( (x_new-x_centre)**2 + (y_new-y_centre)**2 );
*if the distnace between centres is smaller than the two radii then no;
if d2 <= d1 then overlap = "Y";
else overlap = "N";
run;

proc sql noprint;
create table overlap as
select *
from new4
where overlap = "Y";
quit;

%if &sqlobs. = 0 %then %do;

data centres;
set centres new3(rename=(x_new=x_centre y_new=y_centre));
run;

%let i = %sysevalf(&i. + 1);

%end;

*find how many posible centres we have used;
%let it_current = %eval(&it_current. + 1);

%if &it_current >= 700 %then %do;
%let j = %eval(&i.-1);
%put "WARN" "ING: macro has tried 700 possible points without accepting &it_max points. The program will now progress with the &j. accepted points.";
%let i = %eval(&it_max.+1);
%end;

%end;

%end;

%put &it_current.;

%mend;

%centres();

*********************************
PENTAFLAKE
*********************************;

%macro pentaflake(in=, out=, iterations=);

%let n = &iterations.;

*let centre dataset be the centre for the pentaflake we want to make;
data new_centres0;
set &in.;
run;

*if n=0 we want to just record the centre we have and finish;
%if &n. = 0 %then %do;
data &out.;
set new_centres0;
run;
%end;

*if n=1 we want to make 5 new points as well as our original;
%if &n >= 1 %then %do;

data new_centres;
set new_centres0;
*make r2 smaller than the radius of the whole area;
r2 = r/(1+(2*cos(36*(constant("pi")/180))));
*create d based on radius;
d = 2*r2*cos(36*(constant("pi")/180));
*start at 90 and increase for each c;
direction = (90/&it_max.)*&k.;
*make 5 new points off of the original centre point;
do c1=1 to 5;
*update x and y based on the radius and direction;
cos_num = direction*(constant("pi")/180);
sin_num = direction*(constant("pi")/180);
x_new = x_centre+(d*cos(cos_num));
y_new = y_centre+(d*sin(sin_num));
output;
direction = direction + 72;
end;
drop x_centre y_centre;
run;

*add your new centres into the points dataset and set the iteration down to go onto the next round;
*also add new_centres0 so that more points can be made around it;
data new_centres1;
set new_centres0(in=a) new_centres(rename=(y_new=y_centre x_new=x_centre));
if a then do;
c1=0;
r2 = r/(1+(2*cos(36*(constant("pi")/180))));
d = 2*r2*cos(36*(constant("pi")/180));
direction = 90;
end;
n = 1;
run;

%end;

%do n1=2 %to &n.;

%let m = %eval(&n1. - 1);

%do j=0 %to 5;
data new_centres&n1._&j.;
set new_centres&m.;
where c&m. = &j.;
r2 = r2/(1+(2*cos(36*(constant("pi")/180))));
d = 2*r2*cos(36*(constant("pi")/180));
direction = (90/&it_max.)*&k. + (180*&n1.);
do c&n1.=0 to 5;
*update x and y based on the radius and direction;
cos_num = direction*(constant("pi")/180);
sin_num = direction*(constant("pi")/180);
x_new = x_centre+(d*cos(cos_num));
y_new = y_centre+(d*sin(sin_num));
output;
direction = direction + 72;
end;
drop x_centre y_centre;
run;
%end;

data new_centres&n1.;
set new_centres&n1._0 - new_centres&n1._5(rename=( y_new=y_centre x_new=x_centre));
n = &n.;
run;

%end;

data &out.;
set new_centres0-new_centres&n. ;
n=&n.;
it=&k.;
run;

%mend;

****************************************
RECURSIVE ALGORITHM FOR PLOTTING CENTRES FOR EACH PENTAGON
****************************************;

%macro p();

*create empty out dataset;
data points;
x_centre=.; y_centre=.;
run;

%do k=1 %to &it_max.;

*subset centres dataset for specific area;
data centres2;
set centres;
where i=&k.;
r2=r;
run;

*choose number of iterations;
data num_it;
set centres;
if r < 0.7 then num_it=4;
else num_it=5;
run;
proc sql noprint;
select num_it
into :num_it
from num_it;
run;

%put &num_it.;

%pentaflake(in=centres2, out=points&k., iterations=&num_it.);

data points;
set points points&k.;
run;

%end;

%mend;
%p();

data final;
set points;
*set the middle point for each snowflake on its own;
if c1=. then do;
x_c1=x_centre; y_c1=y_centre;
end;
*then set the other 5 n=1 centres in a new group;
else if c2=. and c1 ne 0 then do;
x_c2=x_centre; y_c2=y_centre;
end;
*set the rest of the points in their own group;
else if c3=. then do;
x_c3=x_centre; y_c3=y_centre;
end;
else if c4=. then do;
x_c4=x_centre; y_c4=y_centre;
end;
else if c5=. then do;
x_c5=x_centre; y_c5=y_centre;
end;
run;

*add on text;
data text;
length text $20.;
x_tex=&xmax./2; y_tex=6.9; text="Seasons"; output;
x_tex=&xmax./2; y_tex=6.3; text="Greetings"; output;
run;

*add on boxes for behind text;
data box;
length x_box y_box box_num 8.;
input x_box y_box box_num;
datalines;
2.07 5.65 1
2.07 7.45 1
6.23 7.45 1
6.23 5.65 1
1.99 5.57 2
1.99 7.53 2
6.31 7.53 2
6.31 5.57 2

;
data final2(compress=YES);
set final text box;
run;

*****************************************
PLOT
*****************************************;

title;
filename odsout "P:\Training Area\TrainingCode_UK_2021\Maura Brown";
ods graphics / imagename="christmas2022" height=11.7in width=8.3in antialiasmax=200700;
ods _all_ close;
ods html image_dpi=250 path=odsout file="christmas2022.html";

proc sgplot data=final2 noautolegend noborder ;
styleattrs wallcolor=DEGB ;
xaxis display=none type=linear offsetmin=0 offsetmax=0 min=0 max=&xmax.;
yaxis display=none type=linear offsetmin=0 offsetmax=0 min=0 max=&ymax.;
scatter X=x_c5 Y=y_c5 / markerattrs=(symbol=trianglefilled size=2 color=white);
scatter X=x_c4 Y=y_c4 / markerattrs=(symbol=triangledownfilled size=2.5 color=silver) ;
scatter X=x_c3 Y=y_c3 / markerattrs=(symbol=trianglefilled size=3 color=silver);
scatter X=x_c2 Y=y_c2 / markerattrs=(symbol=triangledownfilled size=8 color=VPAB) ;
scatter X=x_c1 Y=y_c1 / markerattrs=(symbol=trianglefilled size=12 color=VPAB );

polygon x=x_box y=y_box id=box_num / fill fillattrs=(color=cx73002e transparency=0.75) outline lineattrs=(color=cxaf95a6 thickness=3 ) ;
text x=x_tex y=y_tex text=text / textattrs=(family="Georgia" size=35 color=white ) fillattrs=(transparency=0.5);
run;

ods html close;

/*
* KochSnowflakev3.0
*
* Generates fractals and animates them
* We do not talk about versions 1.0 and 2.0 :(
*/

proc datasets library=work nolist;
delete initdata updated poly: mat: vec:;
quit;

%*Register all Windows fonts;
proc fontreg mode=add MSGLEVEL=VERBOSE;
fontpath '!SYSTEMROOT\fonts';
run;

dm log 'clear';

%let out_path = P:\Training Area\TrainingCode_UK_US_2020\James P\Graphics;
%let out_name = snowCard19;
%let out_title = Seasons Greetings;

%*Larger snowflakes take more time to render (keep: large < 5, medium < 10, small < 40);
%let numsmall = 30;
%let nummedium = 6;
%let numlarge = 2;

%*Total number of snowflakes to render;
%let num = %eval(&numsmall. + &nummedium. + &numlarge.);

%*Number of frames in gif;
%let frames = 220;

%*Random seed;
%let rand_seed = 1340;

%*Size of graph;
%let width = 300;
%let height = 300;

/*
* Set up math utilities
*/

%let pi = 3.14159265;
%let err = 2E-6;

/*
* vec2
* creates a two dimensional vector vec2_&name taking in &x and &y as components
* [x]
* [y]
*/
%macro vec2(name = zero, x = 0, y = 0);

%global vec2_&name.;

%let vec2_&name. = &x. &y.;

%mend;

/*
* mat2
* creates a two by two matrix mat2_&name taking in &ii, &ij, &ji, &jj as components
* [ii ij]
* [ji jj]
*/
%macro mat2(name = zero, ii = 0, ij = 0, ji = 0, jj = 0);

%global mat2_&name.;

%let mat2_&name. = &ii. &ij. &ji. &jj.;

%mend;

/*
* mat2_rot
* creates a rotation matrix mat2_rot_&name with a given angle of &theta degrees
*/
%macro mat2_rot(name = zero, theta = 0);

%global mat2_&name.;

%let cos = %sysfunc(cos(&pi. * &theta. / 180));
%let sin = %sysfunc(sin(&pi. * &theta. / 180));

%let mat2_&name. = &cos. -&sin. &sin. &cos.;

%mend mat2_rot;

/*
* mul
* multiplies a matrix &mat and a vector &vec creating a vector &name
*/
%macro mul(name = , mat = , vec = );

%global vec2_&name.;

%let ii = %scan(&&mat2_&mat.., 1, " "); %let ij = %scan(&&mat2_&mat.., 2, " ");
%let ji = %scan(&&mat2_&mat.., 3, " "); %let jj = %scan(&&mat2_&mat.., 4, " ");

%let x = %scan(&&vec2_&vec.., 1, " ");
%let y = %scan(&&vec2_&vec.., 2, " ");

%let vec2_&name. = %sysevalf(&ii. * &x. + &ij. * &y.) %sysevalf(&ji. * &x. + &jj. * &y.);

%mend mul;

/*
* mmul
* multiplies the matrix &mat1 by the matrix &mat2 creating a matrix &name
*/
%macro mmul(name = , mat1 = , mat2 = );

%global mat2_&name.;

%let m1ii = %scan(&&mat2_&mat1.., 1, " "); %let m1ij = %scan(&&mat2_&mat1.., 2, " ");
%let m1ji = %scan(&&mat2_&mat1.., 3, " "); %let m1jj = %scan(&&mat2_&mat1.., 4, " ");

%let m2ii = %scan(&&mat2_&mat2.., 1, " "); %let m2ij = %scan(&&mat2_&mat2.., 2, " ");
%let m2ji = %scan(&&mat2_&mat2.., 3, " "); %let m2jj = %scan(&&mat2_&mat2.., 4, " ");

%let mat2_&name. = %sysevalf(&m1ii. * &m2ii. + &m1ij. * &m2ji.) %sysevalf(&m1ii. * &m2ij. + &m1ij. * &m2jj.)
%sysevalf(&m1ji. * &m2ii. + &m1jj. * &m2ji.) %sysevalf(&m1ji. * &m2ij. + &m1jj. * &m2jj.);

%mend mmul;

/*
* dmul
* multiplies &xvar and &yvar in a dataset &in by the matrix &mat creating a new dataset &out
*/
%macro dmul(mat = , in = , out = , id = , xvar = x, yvar = y);

%let ii = %scan(&&mat2_&mat.., 1, " "); %let ij = %scan(&&mat2_&mat.., 2, " ");
%let ji = %scan(&&mat2_&mat.., 3, " "); %let jj = %scan(&&mat2_&mat.., 4, " ");

data &out. (drop = _&xvar. _&yvar.);
set &in. (rename = (&xvar. = _&xvar. &yvar. = _&yvar.));

id = "&id.";

&xvar. = &ii. * _&xvar. + &ij. * _&yvar.;
&yvar. = &ji. * _&xvar. + &jj. * _&yvar.;

run;
%mend dmul;

/*
* vadd
* adds the vector &vec to &xvar and &yvar in a dataset &in creating a new dataset &out
*/
%macro vadd(vec = , in = , out = , id = , xvar = , yvar = );

%let x = %scan(&&vec2_&vec.., 1, " ");
%let y = %scan(&&vec2_&vec.., 2, " ");

data &out. (drop = _&xvar. _&yvar.);
set &in. (rename = (&xvar. = _&xvar. &yvar. = _&yvar.));

id = "&id.";

&xvar. = &x. + _&xvar.;
&yvar. = &y. + _&yvar.;

run;
%mend vadd;

/*
* dvec2
* converts a list of vectors &vecs into a dataset &name
*/
%macro dvec2(name = , vecs = );

%let i = 1;
data &name.;
%do %until(%scan(&vecs., &i., " ") = );
%let vec = %scan(&vecs., &i., " ");
%let x = %scan(&&vec2_&vec.., 1, " ");
%let y = %scan(&&vec2_&vec.., 2, " ");
x = &x.;
y = &y.;
output;
%let i = %eval(&i. + 1);
%end;;
run;

%mend dvec2;

/*
* Generate koch snowflake fractals
*/

%let res_high = 4;
%let res_med = 2;
%let res_low = 1;

/*
* iterate
* takes in &data and iterates on the Koch snowflake, applying id &poly
*/
%macro iterate( data = , poly = );

data temp1;
length id $20;
set &data.;
num + 1;
id = "&poly.";
run;

proc sql noprint;
select x, y
into :x_1-, :y_1-
from temp1;
quit;

data &data. ( keep = id x y);
length id $20;
set temp1;

output;
if num < &sqlobs then do;
x_0 = symget('x_'!!strip(put(num, 8.0)));
y_0 = symget('y_'!!strip(put(num, 8.0)));
x_1 = symget('x_'!!strip(put(num + 1, 8.0)));
y_1 = symget('y_'!!strip(put(num + 1, 8.0)));
x = (2/3) * x_0 + (1/3) * x_1;
y = (2/3) * y_0 + (1/3) * y_1;
id = "&poly.";
output;
x_0 = symget('x_'!!strip(put(num, 8.0)));
y_0 = symget('y_'!!strip(put(num, 8.0)));
x_1 = symget('x_'!!strip(put(num + 1, 8.0)));
y_1 = symget('y_'!!strip(put(num + 1, 8.0)));
x = x_0/2 + x_1/2 +
((y_0 - y_1)*(((x_1 - x_0)**2 + (y_1 - y_0)**2)**0.5)*sin(&pi /3))/(3*(((x_1 - x_0)**2 + (y_0 - y_1)**2)**0.5));
y = y_0/2 + y_1/2 +
((x_1 - x_0)*(((x_1 - x_0)**2 + (y_1 - y_0)**2)**0.5)*sin(&pi /3))/(3*(((x_1 - x_0)**2 + (y_0 - y_1)**2)**0.5));
id = "&poly.";
output;
x_0 = symget('x_'!!strip(put(num, 8.0)));
y_0 = symget('y_'!!strip(put(num, 8.0)));
x_1 = symget('x_'!!strip(put(num + 1, 8.0)));
y_1 = symget('y_'!!strip(put(num + 1, 8.0)));
x = (1/3) * x_0 + (2/3) * x_1;
y = (1/3) * y_0 + (2/3) * y_1;
id = "&poly.";
output;
end;
else do;
x_0 = symget('x_'!!strip(put(num, 8.0)));
y_0 = symget('y_'!!strip(put(num, 8.0)));
x_1 = symget('x_'!!strip(put(1, 8.0)));
y_1 = symget('y_'!!strip(put(1, 8.0)));
x = (2/3) * x_0 + (1/3) * x_1;
y = (2/3) * y_0 + (1/3) * y_1;
id = "&poly.";
output;
x_0 = symget('x_'!!strip(put(num, 8.0)));
y_0 = symget('y_'!!strip(put(num, 8.0)));
x_1 = symget('x_'!!strip(put(1, 8.0)));
y_1 = symget('y_'!!strip(put(1, 8.0)));
x = x_0/2 + x_1/2 +
((y_0 - y_1)*(((x_1 - x_0)**2 + (y_1 - y_0)**2)**0.5)*sin(&pi /3))/(3*(((x_1 - x_0)**2 + (y_0 - y_1)**2)**0.5));
y = y_0/2 + y_1/2 +
((x_1 - x_0)*(((x_1 - x_0)**2 + (y_1 - y_0)**2)**0.5)*sin(&pi /3))/(3*(((x_1 - x_0)**2 + (y_0 - y_1)**2)**0.5));
id = "&poly.";
output;
x_0 = symget('x_'!!strip(put(num, 8.0)));
y_0 = symget('y_'!!strip(put(num, 8.0)));
x_1 = symget('x_'!!strip(put(1, 8.0)));
y_1 = symget('y_'!!strip(put(1, 8.0)));
x = (1/3) * x_0 + (2/3) * x_1;
y = (1/3) * y_0 + (2/3) * y_1;
id = "&poly.";
output;
end;

run;

proc datasets library=work nolist;
delete temp1;
quit;

%mend;

/*
* koch
* Set up initial conditions and iterate to create a koch snowflake of size &size id &poly and
* iterates &it times, adding to the dataset &data.
*/
%macro koch(data = , it = , size = , poly = );

%vec2(name = center, x = 0, y = 0);
%vec2(name = v1, x = 0, y = %sysevalf(-(&size.)));
%vec2(name = v2, x = %sysevalf(-(&size.) * cos(&pi.*30/180)), y = %sysevalf((&size.) * sin(&pi.*30/180)));
%vec2(name = v3, x = %sysevalf((&size.) * cos(&pi.*30/180)), y = %sysevalf((&size.) * sin(&pi.*30/180)));
%dvec2(name = &data., vecs = v1 v2 v3);

%do koch = 1 %to &it;

%iterate(data = &data., poly = &poly.);

%end;

%mend koch;

/*
* lmh
* Create different resolutions &low &med and &high of the koch snowflake
*/
%macro lmh(low = , med = , high = );

%vec2(name = center, x = 0, y = 0);
%vec2(name = v1, x = 0, y = %sysevalf(-(10)));
%vec2(name = v2, x = %sysevalf(-(10) * %sysfunc(cos(&pi.*30/180))), y = %sysevalf((10) * %sysfunc(sin(&pi.*30/180))));
%vec2(name = v3, x = %sysevalf( (10) * %sysfunc(cos(&pi.*30/180))), y = %sysevalf((10) * %sysfunc(sin(&pi.*30/180))));
%dvec2(name = polyhigh, vecs = v1 v2 v3);

%put &=vec2_v1.;
%put &=vec2_v2.;
%put &=vec2_v3.;

%do koch = 1 %to &high.;

%iterate(data = polyhigh, poly = polyhigh);

%if &koch. = &low. %then %do;
data polylow;
length id $20;
set polyhigh;
id = "polylow";
run;
%end;
%else %if &koch. = &med. %then %do;
data polymed;
length id $20;
set polyhigh;
id = "polymed";
run;
%end;

%end;
%mend lmh;

%*Generate a low medium and high resolution fractal;
%lmh(low = &res_low., med = &res_med., high = &res_high.);

%*Halve the size of the low polytope;
%mat2(name = half, ii = 0.5, ij = 0, ji = 0, jj = 0.5);
%dmul(mat = half, in = polylow, out = polylow, id = polylow, xvar = x, yvar = y);

%*Increase the size of the high polytope;
%mat2(name = x10, ii = 10, ij = 0, ji = 0, jj = 10);
%dmul(mat = x10, in = polyhigh, out = polyhigh, id = polyhigh, xvar = x, yvar = y);

/*
* Simulation
*/

%let g = 9.8;
%let wind_min = 2;
%let wind_max = 5;
%let wind_fore = 330;
%let wind_back = 280;
%let size_max = 20;
%let bound_x = %sysevalf( &width. + (&size_max. * 0.5));
%let bound_y = %sysevalf( &height. + (&size_max. * 0.5));

%put &=width;
%put &=height;

%put &=bound_x;
%put &=bound_y;

/*
* simulate
* Simulates movement over &num many snowflakes for &frames many frames
*/
%macro simulate(num = , frames = );

/*
* v = v_0;
* y = v_0*t + p_0;
*/

/*
* a_x = v_w * cos(t)
* a_y = v_w * sin(t) - g
*
* v_x = v_w * sin(t) + v_0_x
* v_y = v_w * -cos(t) - g * t + v_0_y
*
* p_x = v_w * -cos(t) + v_0_x * t + p_0_x
* p_y = v_w * -sin(t) - 0.5 * g * t^2 + v_0_y * t + p_0_y
*/

%*Temporarily suppress log output (too many datasteps and operations);
options nosource nonotes errors=0;

data polydata;
set _null_;
run;

%*Initialize positions, velocities and rotations;
data initdata;
call streaminit(&rand_seed.);

do i = 1 to &num.;

if i <= &numsmall. then do;
/*pos_x = 8 * &width. * rand("UNIFORM") - 4 * &width.;
pos_y = 2 * &height. * rand("UNIFORM") + 1.5 * &height.;*/

pos_x = 2 * &width. * rand("UNIFORM") - &width.;
pos_y = 2 * &height. * rand("UNIFORM") - &height.;

vel_x = ((&wind_max. - &wind_min.) * rand("UNIFORM") + &wind_min.) * cos(constant("pi") * &wind_back. / 180);
vel_y = ((&wind_max. - &wind_min.) * rand("UNIFORM") + &wind_min.) * sin(constant("pi") * &wind_back. / 180);
end;
else if i <= &numsmall. + &nummedium. then do;
pos_x = - &width.;
pos_y = 0.75 * &height.;

vel_x = ((&wind_max. - &wind_min.) * rand("UNIFORM") + &wind_min.) * cos(constant("pi") * &wind_fore. / 180);
vel_y = ((&wind_max. - &wind_min.) * rand("UNIFORM") + &wind_min.) * sin(constant("pi") * &wind_fore. / 180);
end;
else do;
if i = &numsmall. + &nummedium. + 1 then do;
pos_x = &width. + rand("UNIFORM") * 10;
pos_y = &height. + rand("UNIFORM") * 10;
end;
else do;
pos_x = - &width. - rand("UNIFORM") * 10;
pos_y = - &height. - rand("UNIFORM") * 10;
end;

vel_x = 0;
vel_y = 0;
end;

scale = 0.25 * rand("UNIFORM") + 0.75;

acc_x = 0;
acc_y = 0;

rot_x = 360 * rand("UNIFORM");
rot_y = 360 * rand("UNIFORM");
rot_z = 360 * rand("UNIFORM");

output;

end;
drop i;
run;

%*Put into macro variables;
proc sql noprint;
select pos_x, pos_y, vel_x, vel_y, acc_x, acc_y, rot_x, rot_y, rot_z, scale
into :pos_x_1-,
:pos_y_1-,
:vel_x_1-,
:vel_y_1-,
:acc_x_1-,
:acc_y_1-,
:rot_x_1-,
:rot_y_1-,
:rot_z_1-,
:scale_1-
from initdata;
quit;

%*Update position vectors and create separate dataset per frame;
%do i = 1 %to &frames.;

%do j = 1 %to &num.;

/*%put Frame: &=i. of &=frames.;
%put Polygon: &=j. of &=num.;
%put Small: &=numsmall.;
%put Medium: &=nummedium.;*/

data updated;
%if &j. <= &numsmall. %then %do;
set polylow;
%end;
%else %if &j. <= &numsmall. + &nummedium. %then %do;
set polymed;
%end;
%else %do;
set polyhigh;
%end;
frame = &i.;
%if &i. < 40 %then %do;
text = "Season's Greetings";
%end;
%else %if &i. < 70 %then %do;
text = "from Phastar";
%end;
%else %do;
text = " ";
%end;
tx = 0;
ty = (&i. < 20) * 2 * &height.;
%*tsize = (100 * (&frames. - (&frames. - &i.))) / &frames.; /*lerp*/
tsize = 100;

text2 = "from Phastar";
t2x = 0;
t2y = (&i. < 40) * 2 * &height.;
t2size = 100;
run;

%*Kinematics;
/*
* p_x = v_w * -cos(t) + v_0_x * t + p_0_x
* p_y = v_w * -sin(t) - 0.5 * g * t^2 + v_0_y * t + p_0_y
*/
%if &j. <= &numsmall. %then %do;
%let pos_x_final = %sysevalf(&&vel_x_&j.. * (&i. - 1) + &&pos_x_&j..);
%let pos_y_final = %sysevalf(&&vel_y_&j.. * (&i. - 1) + &&pos_y_&j..);
%if %sysevalf(&pos_x_final. > &bound_x.) %then %do;
%let pos_x_&j. = %sysevalf(&&pos_x_&j.. - 2 * &bound_x.);
%end;
%if %sysevalf(&pos_y_final. < -&bound_y.) %then %do;
%let pos_y_&j. = %sysevalf(&&pos_y_&j.. + 2 * &bound_y.);
%end;
%end;
%else %if &j. <= &numsmall. + &nummedium. %then %do;
%let exp = 100 * %sysfunc(exp(-0.025 * (&i. - 20)**2));

%let pos_x_final = %sysevalf(&exp. * -%sysfunc(cos(&pi. * &i. / 20)) + &&vel_x_&j.. * (&i. - 1) + &&pos_x_&j..);
%let pos_y_final = %sysevalf(&exp. * -%sysfunc(sin(&pi. * &i. / 20)) + &&vel_y_&j.. * (&i. - 1) + &&pos_y_&j..);
%end;
%else %do;
%let pos_x_final = %sysevalf(&&pos_x_&j..);
%let pos_y_final = %sysevalf(&&pos_y_&j..);
%end;

%*Rotations;
%if &j. <= &numsmall. + &nummedium. %then %do;
%let rot_x_final = %sysevalf(&&rot_x_&j. * &i. * 0.01);
%let rot_y_final = %sysevalf(&&rot_y_&j. * &i. * 0.01);
%let rot_z_final = %sysevalf(&&rot_z_&j. * &i. * 0.01);
%end;
%else %do;
%let rot_x_final = 0;
%let rot_y_final = 0;
%let rot_z_final = %sysevalf(&&rot_z_&j. * &i. * 0.0025);
%end;

%*Apply transformations;
%mat2(name = scale, ii = &&scale_&j.., jj = &&scale_&j..);

%mat2_rot(name = rot_z, theta = &rot_z_final.);
%mat2(name = rot_y, ii = %sysfunc(cos(&pi. * &rot_y_final. / 180)), jj = 1);
%mat2(name = rot_x, ii = 1, jj = %sysfunc(cos(&pi. * &rot_x_final. / 180)));

%mmul(name = rot_final, mat1 = rot_y, mat2 = rot_z);
%mmul(name = rot_final, mat1 = rot_x, mat2 = rot_final);

%mmul(name = rot_final, mat1 = rot_final, mat2 = scale);

%vec2(name = move, x = &pos_x_final., y = &pos_y_final.);

%*Update rotations and scale;
%dmul(mat = rot_final, in = updated, out = updated, id = polylow_&j._&i., xvar = x, yvar = y);

%*Update positions;
%vadd(vec = move, in = updated, out = updated, id = polylow_&j._&i., xvar = x, yvar = y);

%*Add all polydata;
data polydata;
set polydata
updated;
run;

%end;

%end;

options source notes errors = 1;

%mend simulate;

%simulate(num = &num., frames = &frames.);

/*
* Gif output
*/

%*Check ods destination;
/*proc print data=sashelp.vdest; run;*/

%*Close ods listing destination;
ods listing close;

ods graphics / imagefmt=GIF width=5in height=5in; /* each image is 5in x 5in GIF */
options papersize=('5 in', '5 in') /* set size for images */
nodate nonumber /* hide date, time, or frame number */
animduration=0.1 animloop=yes noanimoverlay /* animation details */
printerpath=gif animation=start; /* start recording images to GIF */
ods printer file="&out_path.\&out_name..gif"; /* images saved into animated GIF */
proc sgplot data=polydata nowall noborder noautolegend aspect=1;
%*title "&out_title.";
by frame;
%*inset "Merry Christmas" / TEXTATTRS=(Color=Green Family=Arial Size=8 Style=Italic Weight=Bold)
VALUEALIGN=CENTER;
styleattrs backcolor = cx000000;
polygon X=x Y=y id=id / fill
lineattrs=GraphData3 fillattrs=(color=cxfefeff);
text X=tx Y=ty text=text / POSITION = CENTER
TEXTATTRS = (COLOR='cx73002E' family="Monotype Corsiva")
SIZERESPONSE = tsize
SIZEMAXRESPONSE=100
SIZEMIN=0
SIZEMAX=40;
%*text X=t2x Y=t2y text=text2 / POSITION = CENTER
TEXTATTRS = (COLOR='cx73002E' family="Monotype Corsiva")
SIZERESPONSE = t2size
SIZEMAXRESPONSE=100
SIZEMIN=0
SIZEMAX=40;
xaxis offsetmin=0.02 offsetmax=0.02 min=-&width. max=&width. display=NONE;
yaxis offsetmin=0.02 offsetmax=0.02 min=-&height. max=&height. display=NONE;
run;
options printerpath=gif animation=stop; /* stop recording images */
ods printer close; /* close the animated GIF file */

/*
* Image output
*/

/*
ods graphics / noborder;
ods pdf file="P:\Training Area\TrainingCode_UK_US_2020\James P\Graphics\test.pdf";

proc sgplot data=polydata nowall noborder;
styleattrs backcolor = cx000000;
polygon X=x Y=y id=id / fill
lineattrs=GraphData3 fillattrs=(color=cxfefeff);
xaxis min=-100 max=100 display=NONE;
yaxis min=-100 max=100 display=NONE;
run;

ods pdf close;
ods graphics off;
*/

/*
* Check out the following link for more information on animating with proc sgplot in sas
* https://blogs.sas.com/content/iml/2016/08/22/animation-by-statement-proc-sgplot.html#:~:text=You%20can%20control%20the%20animation%20by%20using%20SAS,PRINTERPATH%3DGIF%20ANIMATION%3DSTOP%20after%20you%20generate%20the%20image%20files.
*/