int dotranslate = chi("dotranslate");
int dotrotate = chi("dotrotate");
int doscale = chi("doscale");
int uniformscale = chi("uniformscale");
vector min0,max0,min1,max1;
getbbox(0,min0,max0);
getbbox(1,min1,max1);
vector size0 = max0-min0;
vector center0 = (max0+min0)/2;
vector size1 = max1 -min1;
vector center1 = (max1+min1)/2;

matrix T = 1 ;
matrix R = 1 ;
matrix S = 1 ;

if ( dotranslate){
    vector t = center1 -center0;
    translate(T,t);
    }

if (dotrotate){
    float size0_array[]= set(size0);
    float size1_array[]= set(size1);
    
    int size0_argsort[]= argsort(size0_array);
    int size1_argsort[]= argsort(size1_array);
    vector newsize = size0;
    
    vector mr[];
    for(int i;i<3; i++){
        vector axis;
        axis[size1_argsort[i] ]=1;
        mr[size0_argsort[i] ]=axis;
        
        newsize[size0_argsort[i]]= size0[size0_argsort[i]];
        }

        
    size0 = newsize;
    //mr[2]= cross(mr[0],mr[1]);
    R = set(mr);

}
if(doscale){
    //@P = fit(@P,min0,max0,min1,max1);
    vector s = size1/size0;
    if (uniformscale){
        s = min(s);
        }
        scale(S,s);
}

vector p = center0;matrix X = 1;translate(X,p);

matrix xform = invert(X)*R*S*T*X;

string prim_type_name = primintrinsic(0,"typename",@primnum);
if(match('Packed*',prim_type_name)){
    matrix xform0 = getpackedtransform(0,@primnum);
    setpackedtransform(0,@primnum,xform0*xform);
    }else{
        @P *= xform;
    }
    

int differ_angle = chi("differ_angle");

int q[] = {};
int angle[] = {};

int indegree[] = {};
int flag[] = {};
int cell[] = {};

for (int i = 0; i < npoints(0); i++){
    int _angle = point(0, "angleY", i);
    append(angle, _angle);
    int _indegree = point(0, "indegree", i);
    append(indegree, _indegree);
    append(flag, 0);
    
    append(cell, -1);
    if (_indegree == 0){
        append(q, i);
        flag[i] = 1;
        cell[i] = 0;
    }
}

while (len(q)){
    int next_pts[] = point(0, "next_pts", q[0]);
    foreach (int next_pt; next_pts){
        if (flag[next_pt])  continue;
        indegree[next_pt]--;
        if (indegree[next_pt] == 0){
            int pre_pts[] = point(0, "pre_pts", next_pt);
            
            int fault_angle[] = {};
            foreach (int pre_pt; pre_pts){
                int start_angle = angle[pre_pt] - differ_angle;
                int end_angle = angle[pre_pt] + differ_angle;
                for (int cur_angle = start_angle; cur_angle <= end_angle; cur_angle++)
                    append(fault_angle, cur_angle < 0 ? 360 + cur_angle : cur_angle > 360 ? cur_angle - 360 : cur_angle);
            }
            int perhaps_angle[] = {};
            for (int i = 0; i < 360; i++)
                if (find(fault_angle, i) < 0)
                    append(perhaps_angle, i);
            if (find(fault_angle, angle[next_pt]) >= 0)
                if (len(perhaps_angle) > 0){
                    int index = floor(fit01(rand(next_pt),0, len(perhaps_angle)));
                    angle[next_pt] = perhaps_angle[index];
                 }
            flag[next_pt] = 1;
            cell[next_pt] = cell[q[0]] + 1;
            append(q, next_pt);
        }
    
    }
    removeindex(q, 0);
}

for (int i = 0; i < npoints(0); i++){
    setpointattrib(0, "angleY", i, angle[i], "set");
    setpointattrib(0, "cell", i, cell[i], "set");
}

f@angleY += fit01(rand(@ptnum), 5, -5);
v@angle.y = f@angleY;
matrix3 m = maketransform(2, v@angle);
p@orient = quaternion(m);