*A set of macros for calculation of various geometries based on origin/destination data. *Any questions, email apwheele@gmail.com, Andy Wheeler *Value of Pi. DEFINE !PI ( ) (4*ARTAN(1)) !ENDDEFINE. *Converts from Radians to Degrees. DEFINE !DegToRad (!POSITIONAL !ENCLOSE("(",")")) (((!1)/180)*!PI) !ENDDEFINE. *Converts from Degrees to Radians. DEFINE !RadToDeg (!POSITIONAL !ENCLOSE("(",")")) (((!1)/!PI)*180) !ENDDEFINE. *Calculates the degree (from due north) given a set of from coordinates and to coordinates. *From and To take two tokens, x and then y. DEFINE !DirDeg (Orig = !TOKENS(2) /Dest = !TOKENS(2) /Out = !TOKENS(1) ) !LET !X1 = !HEAD(!Orig) !LET !Y1 = !TAIL(!Orig) !LET !X2 = !HEAD(!Dest) !LET !Y2 = !TAIL(!Dest) DO IF !X2 <> !X1. COMPUTE #deg = !RadToDeg(ARTAN((!Y2 - !Y1)/(!X2 - !X1))). ELSE. COMPUTE #deg = 0. END IF. DO IF !X2 = !X1 AND !Y1 <= !Y2. COMPUTE !Out = 0. ELSE IF !X2 = !X1 and !Y1 > !Y2. COMPUTE !Out = 180. ELSE IF !X1 > !X2. COMPUTE !Out = 270 - #deg. ELSE IF !X1 < !X2. COMPUTE !Out = 90 - #deg. END IF. FORMATS !Out (F3.0). !ENDDEFINE. *The arctan2 function. *Takes inputs in radians and returns radians. DEFINE !ATAN2 (Y = !TOKENS(1) /X = !TOKENS(1) /OUT = !TOKENS(1) ) DO IF !X > 0. COMPUTE !OUT = ARTAN(!Y/!X). ELSE IF !X < 0 AND !Y >= 0. COMPUTE !OUT = ARTAN(!Y/!X) + !PI. ELSE IF !X < 0 AND !Y < 0. COMPUTE !OUT = ARTAN(!Y/!X) - !PI. ELSE IF !X = 0 AND !Y > 0. COMPUTE !OUT = !PI/2. ELSE IF !X = 0 AND !Y < 0. COMPUTE !OUT = -1*(!PI/2). END IF. !ENDDEFINE. *Computes the distance between two points. DEFINE !Dist (Orig = !TOKENS(2) /Dest = !TOKENS(2) /Out = !TOKENS(1) ) COMPUTE !Out = (SQRT((!HEAD(!Orig) - !HEAD(!Dest))**2 + (!TAIL(!Orig) - !TAIL(!Dest))**2)). !ENDDEFINE. *to calculate new x y given seed, orientation(in radians), and distance. DEFINE !Move (Orig = !TOKENS(2) /Orient = !TOKENS(1) /Dist = !TOKENS(1) /Out = !TOKENS(2) ) COMPUTE !HEAD(!Out) = !HEAD(!Orig) + SIN(!Orient)*!Dist. COMPUTE !TAIL(!Out) = !TAIL(!Orig) + COS(!Orient)*!Dist. !ENDDEFINE. *Quadratic bezier curve between two points given a control point. DEFINE !QuadBez (Orig = !TOKENS(2) /Dest = !TOKENS(2) /Cont = !TOKENS(2) /Dens = !DEFAULT(15) !TOKENS(1) /Out = !TOKENS(1) ) !LET !origin_x = !HEAD(!Orig) !LET !origin_y = !TAIL(!Orig) !LET !control_x = !HEAD(!Cont) !LET !control_y = !TAIL(!Cont) !LET !dest_x = !HEAD(!Dest) !LET !dest_y = !TAIL(!Dest) !LET !OutX = !CONCAT(!Out,"_x") !LET !OutY = !CONCAT(!Out,"_y") VECTOR !OutX(!Dens). VECTOR !OutY(!Dens). LOOP #i = 1 to !Dens. COMPUTE #t = #i/(!Dens+1). COMPUTE !OutX(#i) = ((1-#t)**2)*!origin_x + (2*(1-#t)*#t*!control_x) + ((#t**2)*!dest_x). COMPUTE !OutY(#i) = ((1-#t)**2)*!origin_y + (2*(1-#t)*#t*!control_y) + ((#t**2)*!dest_y). END LOOP. !ENDDEFINE. *******************************************************************. *Some examples of use. * DATA LIST FREE / ID X1 Y1 X2 Y2. * BEGIN DATA 1 1 3 2 2 2 5 2 3 4 3 1 1 2 1 4 4 1 3 1 5 1 4 1 5 6 2 5 2 4 END DATA. * FORMATS ALL (F1.0). * VARIABLE LEVEL ID (NOMINAL) X1 TO Y2 (SCALE). * DATASET NAME Test. *Plot of the directions. * GGRAPH /GRAPHDATASET NAME="graphdataset" VARIABLES=X1 Y1 X2 Y2 ID /GRAPHSPEC SOURCE=INLINE. * BEGIN GPL SOURCE: s=userSource(id("graphdataset")) DATA: X1=col(source(s), name("X1")) DATA: Y1=col(source(s), name("Y1")) DATA: X2=col(source(s), name("X2")) DATA: Y2=col(source(s), name("Y2")) DATA: ID=col(source(s), name("ID"), unit.category()) GUIDE: axis(dim(1), label("X")) GUIDE: axis(dim(2), label("Y")) SCALE: linear(dim(1), min(0), max(6)) SCALE: linear(dim(2), min(0), max(6)) ELEMENT: edge(position((X1*Y1)+(X2*Y2)), shape(shape.arrow)) ELEMENT: point(position(X1*Y1), color.interior(color.black), label(ID)) END GPL. *Direction in degrees from due north. * /*!DirDeg Orig = X1 Y1 Dest = X2 Y2 Out = Degrees.*/ EXECUTE. *Convert to radians. * /*COMPUTE Radians = !DegToRad(Degrees).*/ /*COMPUTE Degrees2 = !RadToDeg(Radians).*/ EXECUTE. *Distance between the two points. * /*!Dist Orig = X1 Y1 Dest = X2 Y2 Out = Distance.*/ EXECUTE. *Creating a set of control points given origin, distance and direction (in radians). * COMPUTE #or = !DegToRad(Degrees2+90). *Mid point. * COMPUTE #dis = Distance/2. * /*!Move Orig = X1 Y1 Orient = Radians Dist = #dis Out = MidX MidY.*/ *Control point perpindicular to midpoint at half the distance. * COMPUTE #or = !DegToRad(Degrees2-90). * /*!Move Orig = MidX MidY Orient = #or Dist = #dis Out = ControlX ControlY.*/ FORMATS MidX MidY ControlX ControlY (F1.0). * EXECUTE. *Plotting control points. * GGRAPH /GRAPHDATASET NAME="graphdataset" VARIABLES=X1 Y1 X2 Y2 ID ControlX ControlY /GRAPHSPEC SOURCE=INLINE. * BEGIN GPL SOURCE: s=userSource(id("graphdataset")) DATA: X1=col(source(s), name("X1")) DATA: Y1=col(source(s), name("Y1")) DATA: X2=col(source(s), name("X2")) DATA: Y2=col(source(s), name("Y2")) DATA: ID=col(source(s), name("ID"), unit.category()) DATA: ControlX=col(source(s), name("ControlX")) DATA: ControlY=col(source(s), name("ControlY")) GUIDE: axis(dim(1), label("X")) GUIDE: axis(dim(2), label("Y")) SCALE: linear(dim(1), min(0), max(6)) SCALE: linear(dim(2), min(0), max(6)) ELEMENT: edge(position((X1*Y1)+(X2*Y2)), shape(shape.arrow)) ELEMENT: point(position(X1*Y1), color.interior(color.black), label(ID)) ELEMENT: point(position(ControlX*ControlY), color.interior(color.green)) END GPL. *Creating Bezier curves - then reshaping to plot. * /*!QuadBez Orig = X1 Y1 Dest = X2 Y2 Cont = ControlX ControlY Dens = 15 Out = Bez.*/ EXECUTE. * VARSTOCASES /MAKE X FROM X1 Bez_x1 TO Bez_x15 X2 /MAKE Y FROM Y1 Bez_y1 TO Bez_y15 Y2 /INDEX Type /DROP Radians Degrees2 Distance MidX MidY ControlX ControlY. * RECODE Type (1 = 1)(2 THRU 16 = 2)(17 = 3). * VALUE LABELS Type 1 'Origin' 2 'Bezier' 3 'Destination'. * GGRAPH /GRAPHDATASET NAME="graphdataset" VARIABLES=X Y ID Type /GRAPHSPEC SOURCE=INLINE. * BEGIN GPL SOURCE: s=userSource(id("graphdataset")) DATA: X=col(source(s), name("X")) DATA: Y=col(source(s), name("Y")) DATA: ID=col(source(s), name("ID"), unit.category()) DATA: Type=col(source(s), name("Type"), unit.category()) GUIDE: axis(dim(1), label("X")) GUIDE: axis(dim(2), label("Y")) SCALE: linear(dim(1), min(0), max(6)) SCALE: linear(dim(2), min(0), max(6)) ELEMENT: path(position(X*Y), split(ID)) ELEMENT: point(position(X*Y), color.interior(Type), size(size."4")) END GPL. *Making order of arrows. * DO IF ID = LAG(ID). * COMPUTE XO = LAG(X). * COMPUTE YO = LAG(Y). * END IF. * GGRAPH /GRAPHDATASET NAME="graphdataset" VARIABLES=X Y XO YO ID Type MISSING=VARIABLEWISE /GRAPHSPEC SOURCE=INLINE. * BEGIN GPL SOURCE: s=userSource(id("graphdataset")) DATA: X=col(source(s), name("X")) DATA: Y=col(source(s), name("Y")) DATA: XO=col(source(s), name("XO")) DATA: YO=col(source(s), name("YO")) DATA: ID=col(source(s), name("ID"), unit.category()) DATA: Type=col(source(s), name("Type"), unit.category()) GUIDE: axis(dim(1), label("X")) GUIDE: axis(dim(2), label("Y")) SCALE: linear(dim(1), min(0), max(6)) SCALE: linear(dim(2), min(0), max(6)) ELEMENT: edge(position((XO*YO)+(X*Y)), split(ID), shape(shape.arrow)) END GPL. *To make continuous circular ramps you can map color to half direction. *And transparency (or saturation) to abs value of degrees. * DATASET CLOSE ALL. * INPUT PROGRAM. * COMPUTE #dens = 360. * LOOP #i = 1 TO #dens. * COMPUTE Deg = (#i-1)*(360/#dens). * END CASE. * END LOOP. * END FILE. * END INPUT PROGRAM. * DATASET NAME Circ. * /*COMPUTE #or = !DegToRad(Deg).* * /*!Move Orig = 0 0 Orient = #or Dist = 1 Out = LocX LocY.*/ * COMPUTE Col = (Deg >= 180). * DO IF Deg >= 180. * COMPUTE Sat = ABS(180-Deg). * ELSE. * COMPUTE Sat = ABS(Deg-180). * END IF. * FORMATS Sat (F3.0) Col (F1.0). * EXECUTE. * GGRAPH /GRAPHDATASET NAME="graphdataset" VARIABLES=LocX LocY Sat Col /GRAPHSPEC SOURCE=INLINE. * BEGIN GPL PAGE: begin(scale(500px,500px)) SOURCE: s=userSource(id("graphdataset")) DATA: LocX=col(source(s), name("LocX")) DATA: LocY=col(source(s), name("LocY")) DATA: Sat=col(source(s), name("Sat")) DATA: Col=col(source(s), name("Col"), unit.category()) TRANS: base=eval(0) GUIDE: axis(dim(1), null()) GUIDE: axis(dim(2), null()) GUIDE: legend(aesthetic(aesthetic.color), null()) GUIDE: legend(aesthetic(aesthetic.color.saturation), null()) GUIDE: legend(aesthetic(aesthetic.transparency), null()) COORD: rect(dim(1,2), sameRatio()) SCALE: linear(dim(1), min(-1), max(1)) SCALE: linear(dim(2), min(-1.05), max(.96)) SCALE: cat(aesthetic(aesthetic.color), map(("0",color.green),("1",color.purple))) SCALE: linear(aesthetic(aesthetic.color.saturation), aestheticMinimum(color.saturation."0"), aestheticMaximum(color.saturation."0.8"), min(0), max(180)) SCALE: linear(aesthetic(aesthetic.transparency), aestheticMinimum(transparency."0"), aestheticMaximum(transparency."0.97"), min(0), max(180)) ELEMENT: edge(position((base*base)+(LocX*LocY)), color.saturation(Sat), color(Col), size(size."4"), transparency(Sat)) PAGE: end() END GPL. *******************************************************************.