Turns out I’m not done. I decided to add bilinear increment into the code:
New Code
#@cli rep_samila: _samila_formula,_result_x,_result_y,_result_ang,_result_scale,_result_coordinate_mode={ 0=cartesian_coordinates | 1=polar_coordinates},0<_ref_ratio_a<=1,0<_ref_ratio_b<=1,_ref_surface_x,_ref_surface_y,_ref_surface_ang,_ref_scale,_ref_surface_density>0,_custom_formula_f1,_custom_formula_f2,_custom_formula_x,_custom_formula_y,_output_mapped_coordinates={ 0=coordinates_map | 1=mapped_image },_const_A,_const_B,_const_C,_const_D,...
#@cli : Creates Samila mapping image.
#@cli : Default values: '_samila_formula=0','_result_x=0','_result_y=0','_result_ang=0','_result_scale=75%','_result_coordinate_mode=1','_ref_ratio_a=0','_ref_ratio_b=0','_ref_surface_x=0','_ref_surface_y=0',_ref_surface_ang=0','_ref_surface_density=10','_custom_formula_f1=n/a','_custom_formula_f2=n/a','_custom_formula_x=n/a','_custom_formula_y=n/a','_output_mapped_coordinates=1','_const_A=u','_const_B=u','_const_C=u','_const_D=u'
rep_samila:
skip ${1=0},${2=0},${3=0},${4=0},${5=75%},${6=0},${7=1},${8=0},${9=0},${10=0},${11=0},${12=1},${13=10},${14=},${15=},${16=},${17=},${18=1},${19=u(.01,10)},${20=u(.01,10)},${21=u(.01,10)},${22=u(.01,10)},${23=u(.01,10)},${24=u(.01,10)},${25=u(.01,10)},${26=u(.01,10)},${27=u(.01,10)},${28=u(.01,10)},${29=u(.01,10)},${30=u(.01,10)},${31=u(.01,10)},${32=u(.01,10)},${33=u(.01,10)},${34=u(.01,10)},${35=u(.01,10)},${36=u(.01,10)},${37=u(.01,10)},${38=u(.01,10)},${39=u(.01,10)}
check "isfinite($1) && isnum($1) &&
isfinite($2) && isnum($2) &&
isfinite($3) && isnum($3) &&
isfinite($4) && isnum($4) &&
((isfinite($5) && isnum($5) ) || ispercentage($5) ) && $5 &&
isbool($6) &&
isfinite($7) && isnum($7) &&
isfinite($8) && isnum($8) &&
isfinite($9) && isnum($9) &&
isfinite($10) && isnum($10) &&
isfinite($11) && isnum($11) &&
((isfinite($12) && isnum($12) ) || ispercentage($12) ) && $12 &&
((isfinite($13) && isnum($13) ) || ispercentage($13) ) && $13 &&
(narg($14)?is_expr($14):1) &&
(narg($15)?is_expr($15):1) &&
(narg($16)?is_expr($16):1) &&
(narg($17)?is_expr($17):1) &&
isbool($18) &&
isfinite($19) &&
isfinite($20) &&
isfinite($21) &&
isfinite($22)"
samila_mode,result_x,result_y,result_ang,\
result_coordinate_mode,\
ref_ratio_a,ref_ratio_b,ref_surface_x,ref_surface_y,ref_surface_ang,ref_scale=${1-4},${6-12}
result_scale,ref_surface_density:=$5,$13
custom_formula_f1=$14
custom_formula_f2=$15
custom_formula_x=$16
custom_formula_y=$17
output_mode=$18
A,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R,S,T:=${19-38}
use_custom_ratio:=$ref_ratio_a&&$ref_ratio_b
foreach {
if (s>1)||(stats(#-1)!=[0,0,0,0,0,0,0,0,0,0,0,0,0,0])
100%,100% rm..
elif
!s continue
fi
{round((($use_custom_ratio?(or=w/h;nr=$ref_ratio_a/$ref_ratio_b;or!=nr))?(nr>or?[h*nr,h]:[w,w/nr]):[w,h])*$ref_surface_density)},1,2,:"
begin(
rot_x(a,b,sin_ang,cos_ang)=a*cos_ang-b*sin_ang;
rot_y(a,b,sin_ang,cos_ang)=a*sin_ang+b*cos_ang;
const samila_mode=$samila_mode;
const icx=w#-1/2;
const icy=h#-1/2;
const itx=icx+icx*$result_x;
const ity=icy-icy*$result_y;
const scale_factor=(min(icx,icy)/(2*pi))*$result_scale;
const use_polar_mode=$result_coordinate_mode;
const result_ang=$result_ang°;
const rotate_result_mode=$result_ang%360;
if(rotate_result_mode,
const result_sin_ang=sin(result_ang);
const result_cos_ang=cos(result_ang);
);
const max_ref_x=w-1;
const max_ref_y=h-1;
const ref_cx=max_ref_x/2;
const ref_cy=max_ref_y/2;
const ref_shift_x=$ref_surface_x;
const ref_shift_y=$ref_surface_y;
const ref_sd=max(w,h)/min(w,h);
const ref_scale=$ref_scale;
const ref_ang=$ref_surface_ang°;
const ref_sx=1/((w>h?ref_cx/ref_sd:ref_cx)/(pi*ref_scale));
const ref_sy=1/((w>h?ref_cy:ref_cy/ref_sd)/(pi*ref_scale));
const ref_rotate_mode=$ref_surface_ang%360;
const ref_shift_mode=ref_shift_x||ref_shift_y;
if(ref_rotate_mode,
const ref_sin_ang=sin(ref_ang);
const ref_cos_ang=cos(ref_ang);
);
ref_rotate_mode&&ref_shift_mode?(
create_ref_surface()=(
x=(x-ref_cx)*ref_sx;
y=(y-ref_cy)*ref_sy;
tx=x;
x=rot_x(x,y,ref_sin_ang,ref_cos_ang);
x=rot_y(tx,y,ref_sin_ang,ref_cos_ang);
x-=ref_shift_x;
y-=ref_shift_y
);
):(
create_ref_surface()=(
x=(x-ref_cx)*ref_sx-ref_shift_x;
y=(y-ref_cy)*ref_sy-ref_shift_y;
if(ref_rotate_mode,
tx=x;
x=rot_x(x,y,ref_sin_ang,ref_cos_ang);
y=rot_y(tx,y,ref_sin_ang,ref_cos_ang);
);
);
);
const A=$A;
const B=$B;
const C=$C;
const D=$D;
const E=$F;
const F=$F;
const G=$G;
const H=$H;
const I=$I;
const J=$J;
const K=$K;
const L=$L;
const M=$M;
const N=$N;
const O=$O;
const P=$P;
const Q=$Q;
const R=$R;
const S=$S;
const T=$T;
samila_mode>=0?(
rx(x,y)=f1(x,y);
ry(x,y)=f2(x,y);
);
!samila_mode?(
f1(x,y)=(A*.01)*x^C-sin(y^E)+abs(y-x);
f2(x,y)=(B*.01)*y^D-cos(x^F)+G*x;
):
samila_mode==1?(
f1(x,y)=cos(B*tan(x^A))-sin(C*y)+x*y;
f2(x,y)=sin(D*y^A)+cos(B*x)-x;
):
samila_mode==2?(
f1(x,y)=sin(cos(x^A))^C-sin(y^E)*x+abs(x-y);
f2(x,y)=cos(sin(x^B))^D+cos(y^F)*y+abs(y+x)*G;
):
samila_mode==3?(
f1(x,y)=sin(A*y)+C*cos(E*x)*y;
f2(x,y)=sin(B*x)+D*cos(F*y)*x;
):
samila_mode==4?(
f1(x,y)=sin(log(1+abs(x+A))-y)*cos(x^C)*E+abs(y-x);
f2(x,y)=cos(log(1+abs(y+B))+x)*sin(y^D)*F-G*x;
):
samila_mode==5?(
f1(x,y)=A+log(1+abs(x^C))-sin(gamma(y))+abs(y-x);
f2(x,y)=B-log(1+abs(y^D))-cos(gamma(x))-abs(x+y);
):
samila_mode==6?(
f1(x,y)=sin(cos(x*A)*y-atan(x))+x*C;
f2(x,y)=cos(sin(y*B)*x+atan(y))-y*D;
):(
f1(x,y)="$custom_formula_f1";
f2(x,y)="$custom_formula_f2";
rx(x,y)="$custom_formula_x";
ry(x,y)="$custom_formula_y";
);
);
create_ref_surface();
nx=rx(x,y);
ny=ry(x,y);
if(use_polar_mode,
tx=nx;
nx*=cos(ny);
ny=tx*sin(ny);
);
if(rotate_result_mode,
tx=nx;
nx=rot_x(nx,ny,result_sin_ang,result_cos_ang);
ny=rot_y(tx,ny,result_sin_ang,result_cos_ang);
);
[nx,ny]*scale_factor+[itx,ity];
"
if !$output_mode continue fi
eval. ":
inrange(i0,-1,w#-2,1,1)&&inrange(i1,-1,h#-2,1,1)?(
A=isint(i0);B=isint(i1);
A&&B?(
++i(#-2,i0,i1);
):
A?(
Y=floor(i1);
top=Y+1-i1;
bottom=1-top;
i(#-2,i0,Y++)+=top;
i(#-2,i0,Y)+=bottom;
):
B?(
X=floor(i0);
left=X+1-i0;
right=1-left;
i(#-2,X++,i1)+=left;
i(#-2,X,i1)+=right;
):(
X=floor(i0);
Y=floor(i1);
cover_x0=X+1-i0;
cover_y0=Y+1-i1;
cover_x1=1-cover_x0;
cover_y1=1-cover_y0;
i(#-2,X,Y)+=cover_x0*cover_y0;
i(#-2,X+1,Y)+=cover_x1*cover_y0;
i(#-2,X,Y+1)+=cover_x0*cover_y1;
i(#-2,X+1,Y+1)+=cover_x1*cover_y1;
);
);I;"
rm.
}