Data input

Fluxmagdata

Needs["Histograms`"]

fluxmagdata = Import["/Users/ccarter/People/Neri/MeshSubject_03/vertical_results/jmag.csv","CSV"];

Length[fluxmagdata]

flux scaled so that min - max is 0,1

fluxvalues = Rescale[Transpose[fluxmagdata][[3]]];

Length[fluxvalues]

Histogram[fluxvalues]

jmagmax =Max[fluxvalues]

jmagmin =Min[fluxvalues]

Use Log Flux

Histogram[Log[fluxvalues+.00000001]]

ChopMaxMin[vals_, min_, max_] := vals/.{a_/;a<min->min,b_/;b>max->max}

Histogram[ChopMaxMin[Log[fluxvalues+.00000001],-12,-0.5]]

scaledfluxmagvalues= Rescale[ChopMaxMin[Log[fluxvalues+.00000001],-12,-0.5]];

Histogram[scaledfluxmagvalues]

scaledfluxmagdata=Transpose[ Append[Most[Transpose[fluxmagdata]],scaledfluxmagvalues]];

fluxmagdata[[123]]
scaledfluxmagdata[[123]]

stressdata

stressdata = Import["/Users/ccarter/People/Neri/MeshSubject_03/vertical_results/stress.csv","CSV"];

Length[stressdata]

stressvalues = Rescale[Transpose[stressdata][[3]]];

Histogram[stressvalues]

strain scaled so that min - max is 0.5 , 4 (for line widths) (chopped data between 0.35 and 0.6)

scaledstressvalues=0.5 + 3.5Rescale[ChopMaxMin[stressvalues,0.35,0.6]];

Min[scaledstressvalues]
Max[scaledstressvalues]

Histogram[scaledstressvalues]

scaledstressdata=Transpose[ Append[Most[Transpose[stressdata]],scaledstressvalues]];

stressdata[[123]]
scaledstressdata[[123]]

straindata

straindata = Import["/Users/ccarter/People/Neri/MeshSubject_03/vertical_results/strainyy.csv","CSV"];

Length[straindata]

strainvalues = Rescale[Transpose[straindata][[3]]];

Histogram[strainvalues]

strain scaled so that min - max is 0.5 , 4 (for line widths) (chopped data between 0.05 and 0.7)

scaledstrainvalues=0.5 + 3.5Rescale[ChopMaxMin[strainvalues,0.05,0.7]];

Min[scaledstrainvalues]
Max[scaledstrainvalues]

Histogram[scaledstrainvalues]

Histogram[Log[scaledstrainvalues]]

scaledstraindata=Transpose[ Append[Most[Transpose[straindata]],scaledstrainvalues]];

straindata[[153]]
scaledstraindata[[153]]

energydata

energydata = Import["/Users/ccarter/People/Neri/MeshSubject_03/vertical_results/energy.csv","CSV"];

Length[energydata]

energyvalues = Rescale[Transpose[energydata][[3]]];

Histogram[energyvalues]

Use Log Energy

Histogram[Log[energyvalues+.00000001]]

energy scaled so that min - max is 0,100 (for numbers)

scaledenergyvalues= 100*Rescale[ChopMaxMin[Log[energyvalues+.00000001],-12,0]];

Histogram[scaledfluxmagvalues]

scaledenergydata=Transpose[ Append[Most[Transpose[energydata]],scaledenergyvalues]];

energydata[[153]]
scaledenergydata[[153]]

zdata

zdata = Import["/Users/ccarter/People/Neri/MeshSubject_03/vertical_results/z.csv","CSV"];

Length[zdata]

zvalues = Rescale[Transpose[zdata][[3]]];

Histogram[zvalues]

z scaled so that min - max is 0,3 (roughly 10%)(chopped 0.5,0.9)

scaledzvalues= 70*Rescale[ChopMaxMin[zvalues,0.5,0.9]];;

Histogram[scaledzvalues]

scaledzdata=Transpose[ Append[Most[Transpose[zdata]],scaledzvalues]];

zdata[[513]]
scaledzdata[[513]]

Interpolation of Data

fluxfunc

fluxfunc = Interpolation[scaledfluxmagdata/.{a_,b_,c_}->{{b,a},c}];

{{xmin,xmax},{ymin,ymax}}=fluxfunc[[1]]

Plot3D[fluxfunc[x,y],{x,xmin,xmax},{y,ymin,ymax},PlotRange-> All,PlotPoints->64]

stressfunc

stressfunc = Interpolation[scaledstressdata/.{a_,b_,c_}->{{b,a},c}];
Plot3D[stressfunc[x,y],{x,xmin,xmax},{y,ymin,ymax},PlotRange->All,PlotPoints-> 64]

strainfunc

strainfunc = Interpolation[scaledstraindata/.{a_,b_,c_}->{{b,a},c}];
Plot3D[strainfunc[x,y],{x,xmin,xmax},{y,ymin,ymax},PlotRange->All,PlotPoints-> 64]

energyfunc

energyfunc = Interpolation[scaledenergydata/.{a_,b_,c_}->{{b,a},c},InterpolationOrder-> 3];
Plot3D[energyfunc[x,y],{x,xmin,xmax},{y,ymin,ymax},PlotRange->All,PlotPoints->75]

zfunc

zfunc = Interpolation[scaledzdata/.{a_,b_,c_}->{{b,a},c},InterpolationOrder-> 3];
Plot3D[zfunc[x,y],{x,xmin,xmax},{y,ymin,ymax},PlotRange->All]

Skeleton Structure

read data

pts = Import["/Users/ccarter/People/Neri/MeshSubject_03/Medium_Fine_Mesh.points.csv","CSV"];

pgons = Import["/Users/ccarter/People/Neri/MeshSubject_03/Medium_Fine_Mesh.nodes.csv","CSV"]+1;

Length[pts]
Length[pgons]

pts[[123]]
pgons[[123]]

reconcile data structure: make edges

currentedge=1;
edgesOfVertex= Table[{},{Length[pts]}];
UniqueEdge[va_, vb_] :=Which[
ValueQ[edge[va,vb]],edge[va,vb],
ValueQ[edge[vb,va]],-edge[vb,va],
True,edge[va,vb] = currentedge;edgevertices[currentedge]={va,vb};AppendTo[edgesOfVertex[[va]],currentedge]; AppendTo[edgesOfVertex[[vb]],currentedge];currentedge++]

MakeEdges[ v_List] :=Module[{len=Length[v]},
Return[
Append[
 Table[UniqueEdge[v[[i]],v[[i+1]]],{i,1,len-1}],
UniqueEdge[v[[len]],v[[1]]]
]

]
]

Map[MakeEdges,pgons];

currentedge

edgevertices[1234]

edgevertexlist= Map[edgevertices[#]&,Range[currentedge-1]];

edgevertexlist[[1234]]

temp=Map[Union,edgesOfVertex];

temp[[1234]]

edgesOfVertex[[1234]]

Dimensions[edgesOfVertex]

Dimensions[temp]

Dimensions[Flatten@edgesOfVertex]
Dimensions[Flatten@edgesOfVertex]

edgesOfVertex//Short

Union[{66901,66902,66903,66904,66905,66907}]

Length[edgevertexlist]
currentedge-1

starEdgesVertices[v_] := Select[edgevertexlist,MemberQ[#,v]&] (*slow, don't use)*)

starEdges[v_] := Map[UniqueEdge[First[#],Last[#]]& ,starEdgesVertices[v]] (*slow don't use*)

currentedgerange= Range[currentedge-1];

starEdges[v_] := (*slow don't use*)
Pick[currentedgerange,Map[MemberQ[#,v]&,edgevertexlist]]
(*Map[UniqueEdge[First[#],Last[#]]& ,starEdgesVertices[v,edgevertices]]*)

starEdgesVertices[5228]

starEdges[5228]

edgesOfVertex[[5228]]

UniqueEdge[5228,5230]

pts[[234]]

Max[pts]

Min[pts]

scaledzdata[[234]]

zht[vtx_] :=Module[{pt=pts[[vtx]]},zfunc[pt[[2]],pt[[1]]]]

pts//Short

maxx=Max[ Map[First,pts]]
maxy=Max[ Map[Last,pts]]
minx=Min[ Map[First,pts]]
miny=Min[ Map[Last,pts]]

ptsscaled= pts*30/500;

ptsscaled//Short

zhts= Table[zht[i],{i,1,Length[pts]}];

zmax = Max[zhts]
zmin=Min[zhts]

ptszhtscaled=Table[
{
30(First[pts[[i]]]-xmin)/(xmax-xmin),
30(Last[pts[[i]]]-ymin)/(ymax-ymin),
 6(zht[i]-zmin)/(zmax-zmin)
} ,{i,1,Length[pts]}];

ptszhtscaled[[245]]

Max[Last/@ptszhtscaled]

Min[Last/@ptszhtscaled]

Max[First/@ptszhtscaled]

Min[First/@ptszhtscaled]

MidPoint[edge_] := Mean[pts[[edgevertices[edge]]]]

strainmag[edge_] :=strainmag[edge]=Module[{pt=MidPoint[edge]},strainfunc[pt[[2]],pt[[1]]]]

jmag[edge_] :=Module[{pt=MidPoint[edge]},fluxfunc[pt[[2]],pt[[1]]]]

temp = Table[strainmag[i],{i,1,currentedge-1}];
strainmin = Min[Table[strainmag[i],{i,1,currentedge-1}]]
strainmax = Max[Table[strainmag[i],{i,1,currentedge-1}]]

strainmax = Max[Table[strainmag[i],{i,1,currentedge-1}]]

Histogram[temp]

(*strainmagdata= Table[0.06 +  (.4 - 0.06)strainmag[i]/(strainmax-strainmin), {i,1,currentedge-1}];*)
strainmagdata= Table[0.06 +  (.8 - 0.06)strainmag[i]/(strainmax-strainmin), {i,1,currentedge-1}];

Histogram[strainmagdata]

jmagdata= Table[jmag[i], {i,1,currentedge-1}];

Histogram[temp]

maxedgevalue[vertex_, edgevaluelist_]:= Module[
{edges,vals,max},
edges = edgesOfVertex[[vertex]];
vals=Map[edgevaluelist[[#]]&, edges];
Flatten@{max = Max[vals],
Pick[edges,vals, max]}
]

avededgevalue[vertex_, edgevaluelist_]:= Module[
{edges,vals,max},
edges = edgesOfVertex[[vertex]];
vals=Map[edgevaluelist[[#]]&, edges];
Flatten@{max = Mean[vals],
Pick[edges,vals, max]}
]

ptszhtscaled[[451]]

Here are the dimensions (0, 7.5, 15, 22.5, 30)

Writing Panels

a=-.001; b=7.5; c=15.0;d=22.5; e=30.001

xl=d;
 xr=e; 
yd= d;
 yu=e;
AcceptedVertexCoords= Select[ptszhtscaled,((xl≤ #[[1]] ≤xr)&& (yd≤ #[[2]] ≤ yu))&];
Length[AcceptedVertexCoords]
AcceptedVertexCoords[[1;;10]]

filename  = "/Users/ccarter/People/Neri/Butterfly/Panels/Vertical/Vertical__"<>ToString[Round[10 xl]]<>"-"<>ToString[Round[10 xr]]<> "__"<> ToString[Round[10 yd]]<> "-"<> ToString[Round[10 yu]]<>".wrl"

Length[ptszhtscaled]/16//N

currentedge

Length[ptszhtscaled]

VertexIndices =
Pick[Range[Length[ptszhtscaled]],Map[((xl≤ #[[1]] ≤xr)&& (yd≤ #[[2]] ≤ yu))&,ptszhtscaled]];
(*Length[AcceptedVertexIndices]
AcceptedVertexIndices[[1;;10]]*)

xl 
xr
yd
yu

InteriorPointQ[x_,y_] := xl ≤ x ≤  xr && yd≤ y ≤ yu
ExteriorPointQ[x_,y_] := x < xl || x > xr || y < yd || y > yu
InteriorPointQ[x_List] := xl ≤ x[[1]] ≤  xr && yd≤ x[[2]] ≤ yu
ExteriorPointQ[x_List] := x[[1]] < xl || x[[1]] > xr || x[[2]] < yd || x[[2]] > yu

InteriorEdgeQ[edge_] := Module[{v1,v2,p1x,p1y,p1z,p2x,p2y,p2z},{v1,v2}= edgevertices[edge];
{p1x,p1y,p1z}= ptszhtscaled[[v1]];
{p2x,p2y,p2z}= ptszhtscaled[[v2]];
InteriorPointQ[p1x,p1y] && InteriorPointQ[p2x,p2y]
]

(*BorderEdgeQ[edge_] := Module[{v1,v2,p1x,p1y,p1z,p2x,p2y,p2z,t1,t2},{v1,v2}= edgevertices[edge];
{p1x,p1y,p1z}= ptszhtscaled[[v1]];
{p2x,p2y,p2z}= ptszhtscaled[[v2]];
t1= InteriorPointQ[p1x,p1y];
t2=InteriorPointQ[p2x,p2y];
(t1 && !t2) || (!t1 && t2)
(InteriorPointQ[p1x,p1y] && ExteriorPointQ[p2x,p2y]) || (InteriorPointQ[p2x,p2y] && ExteriorPointQ[p1x,p1y] )
]*)

ExteriorEdgeQ[edge_] :=  Module[{v1,v2,p1x,p1y,p1z,p2x,p2y,p2z},{v1,v2}= edgevertices[edge];
{p1x,p1y,p1z}= ptszhtscaled[[v1]];
{p2x,p2y,p2z}= ptszhtscaled[[v2]];
ExteriorPointQ[p1x,p1y] && ExteriorPointQ[p2x,p2y]
]

BorderEdgeQ[edge_] := And[Not[InteriorEdgeQ[edge]],Not[ExteriorEdgeQ[edge]]]

InteriorEdges= Select[Range[currentedge-1],
InteriorEdgeQ[#]&];
InteriorEdges//Short
InteriorEdges//Length

ptszhtscaled[[23,1;;2]]

InteriorVertexIndices= Select[Range[Length[ptszhtscaled]],
(InteriorPointQ@@ptszhtscaled[[#,1;;2]])&];
InteriorVertexIndices//Short
InteriorVertexIndices//Length

ExteriorEdges= Select[Range[currentedge-1],
ExteriorEdgeQ[#]&];
ExteriorEdges//Short
ExteriorEdges//Length

currentedge-1 - Length[ExteriorEdges] - Length[InteriorEdges]

BorderEdges= Select[Range[currentedge-1],
BorderEdgeQ[#]&];
BorderEdges//Short
BorderEdges//Length

Map[ptszhtscaled[[#]]&,(edgevertices[InteriorEdges[[1]]])]
(*Map[ptszhtscaled[[#]]&,(edgevertices[ExteriorEdges[[1]]])]
Map[ptszhtscaled[[#]]&,(edgevertices[BorderEdges[[1]]])]*)

edgevertices[12]
ptszhtscaled[[4734]]

CutEdge[borderedge_] := Module[{ e1,e2,trx,tlx,tdy,tuy,tmin,pin,pout,x1,x2,dx,dy},
{e1,e2}= edgevertices[borderedge];
x1 = ptszhtscaled[[e1]];
x2 = ptszhtscaled[[e2]];

(*Print[{InteriorPointQ[x1],InteriorPointQ[x2]}];
Print[{ExteriorPointQ[x1],ExteriorPointQ[x2]}];*)
If[InteriorPointQ[x1] && InteriorPointQ[x2],Print["Two interior cut points", {x1,x2}]];
If[ExteriorPointQ[x1] &&ExteriorPointQ[x2],Print["Two exterior cut points", {x1,x2}]];

If[InteriorPointQ[x1],pin=x1; pout=x2,pin=x2;pout=x1];



(*Print[{pin,pout}];*)
dx=(pin[[1]]-pout[[1]]);
dy = (pin[[2]]-pout[[2]]);

(*Print[{pin,pout,dx,dy}];*)
If[PossibleZeroQ[dx],Print["possible dx zero", {N[pin[[1]],15],N[pout[[1]],15]},{pin,pout}];dx=10^(-12)];
If[PossibleZeroQ[dy],Print["possible dy zero", {N[pin[[2]],15],N[pout[[2]],15]}];dy=10^(-12)];
trx= If[PossibleZeroQ[(xr-pout[[1]])],10^16,(xr-pout[[1]])/dx];

tlx= If[PossibleZeroQ[(xl-pout[[1]])],10^16,(xl-pout[[1]])/dx];

tuy= If[PossibleZeroQ[(yu-pout[[2]])],10^16,(yu-pout[[2]])/dy];

tdy= If[PossibleZeroQ[(yd-pout[[2]])],10^16,(yd-pout[[2]])/dy];

tmin = Min[Select[{trx,tlx,tuy,tdy},(0< #<1)&]];
(*Print[{tmin,tlx,tlx,tuy,tdy}];*)
{pin, tmin pin + (1-tmin)pout}
]

(*SelectedEdgeVertices[i_] :=Map[PlaceInAcceptedList[#,AcceptedVertexIndices]&,edgevertices[i]]*)

interiorlistEdges =
Table[
{
Hue[0.66jmag[i]],
AbsoluteThickness[strainmagdata[[i]]],
Line[edgevertices[i]]
}
,
{i,InteriorEdges}];
interiorlistEdges//Short

BorderEdges[[2]]

ProgressIndicator[Dynamic[n],{1,currentedge-1}]

InteriorEdgeQ[798]

InteriorPointQ[5069]

InteriorPointQ[9.059605618953016`,0.`]
InteriorPointQ[8.905509599718043`,0.`]

edgevertices[798]

ptszhtscaled[[5069]]
ptszhtscaled[[5070]]

edgevertices[2]

Map[ptszhtscaled[[#]]&,edgevertices[2]]

(*borderlistCylinders*)

(*Length[clistSelected]
clistSelected[[1;;23]]*)

(*slist=Table[n=i;em = maxedgevalue[i,strainmagdata];
{Hue[0.66First[avededgevalue[i,jmagdata]]],EdgeForm[],Sphere[i,First[em]]},{i,1,Length[pts]}];*)

(*slistSelected=
Table[{Hue[0.66First[avededgevalue[i,jmagdata]]],Sphere[i,First[maxedgevalue[i,strainmagdata]]]},
{i,InteriorVertexIndices}];
slistSelected//Short
slistSelected//Length*)

(*slistSelected=
Table[{Hue[0.66First[avededgevalue[i,jmagdata]]],Sphere[ptszhtscaled[[i]],First[maxedgevalue[i,strainmagdata]]]},
{i,InteriorVertexIndices}];
slistSelected//Short
slistSelected//Length*)

Min[Table[strainmagdata[[i]],{i,1,currentedge-1}]]

Max[Table[strainmagdata[[i]],{i,1,currentedge-1}]]

colfunc[val_] :=
N@ Blend[
{
RGBColor[0,0.5 ,0.99],
RGBColor[0,0.67,0.99],
RGBColor[0,0.825,1],
RGBColor[0,1,1],
RGBColor[0,1,0.5],
RGBColor[0,1,0.33],
RGBColor[0,1,.165],
RGBColor[0,1,0],
RGBColor[0.5,1,0],
RGBColor[0.66,1,0],
RGBColor[0.825,1,0],
RGBColor[1,1,0],
RGBColor[1,0.5,0],
RGBColor[1,0.33,0],
RGBColor[1,.167,0],
RGBColor[1,0,0]
},1-val]

42/255//N
84/255//N
168/255//N
211/255//N

colfunc[0]

WriteCylinder[i_,{v1_,v2_}, rad_, hue_]:=
Module[{axis, angle, midpoint,color,height,red,green,blue},
angle = Round[VectorAngle[{0,1,0},v2-v1],.001];
axis = Map[Round[#,.01]&,Normalize[Cross[{0,1,0},v2-v1]]];
midpoint =Map[ Round[#,.001]&,(v1+v2)/2];
height = Round[Norm[v2-v1],.001];
color  = colfunc[hue];
red=Clip[color[[1]],{0.001,1.0}];
green = Clip[color[[2]],{0.001,1.0}];
blue =Clip[ color[[3]],{0.001,1.0}];

(*color=Blend[{Red,Green,Blue},hue];*)(*ColorData["Rainbow"][1-hue];*)(* Blend[{Red,Blue},hue];*)
(*color=Blend[{Red,Green,Blue},hue];*)
WriteString[outstream,"#ic ",i, "\nTransform {\n translation ", midpoint[[1]]," ", midpoint[[2]], " ",midpoint[[3]], "\n rotation ", axis[[1]], " " , axis[[2]], " ", axis[[3]], " ", angle,"\n children\n  Shape {\n   appearance Appearance {material Material { diffuseColor " , red, " ", green," ", blue," }}\n geometry Cylinder {radius ",rad, " height ", height, "  }\n}\n}\n"];
]

WriteSphere[v1_, rad_, hue_]:=
Module[{color,red,green,blue,p1,p2,p3},
color = colfunc[hue];(*Blend[{Red,Green,Blue},hue];*)(*ColorData["Rainbow"][1-hue];*)(*color=Blend[{Red,Blue},hue];*)
red=Clip[color[[1]],{0.001,1.0}];
green = Clip[color[[2]],{0.001,1.0}];
blue =Clip[ color[[3]],{0.001,1.0}];
p1 = Round[v1[[1]],.01];
p2 = Round[v1[[2]],.01];
p3 = Round[v1[[3]],.01];
WriteString[outstream,"Transform {\n translation ", p1," ", p2, " ",p3, "\n children\n  Shape {\n appearance Appearance {material Material { diffuseColor " , red, " ", green," ", blue," }}\n geometry Sphere {radius ",rad, " }\n }\n}\n"];
]

Clear[outstream]

outstream =OpenWrite[filename]
WriteString[outstream,"#VRML V2.0 utf8\n\n Viewpoint {\n
        position 15 15 50\n
        orientation 0 0 -1 0\n
    }\n\n
"];

icount = 1;
Table[WriteCylinder[icount++,
Map[ptszhtscaled[[#]]&,edgevertices[i]],0.8*strainmagdata[[i]],jmag[i]],{i,InteriorEdges}];

Table[WriteCylinder[i,
CutEdge[i],0.8*strainmagdata[[i]],jmag[i]],{i,BorderEdges}];

Table[WriteSphere[ptszhtscaled[[i]],0.8*First[maxedgevalue[i,strainmagdata]],
First[avededgevalue[i,jmagdata]]],{i,InteriorVertexIndices}];

Close[outstream]

BorderEdges

Table[i,{i,BorderEdges}]