(*  The following are examples of curvature analysis  on  individual tubercles in seeds of Silene species*)






	(*******************************  S. spinescens 1 to 6  *)

Needs["Splines`"]

spinescens[1]={{0.2806,0.2806},{1.4029,1.1223},{2.5252,2.2447},{3.6476,2.8058},{4.7699,3.0864},{5.8922,3.6476},{7.0146,4.2087},{8.1369,4.7699},{9.2592,5.0505},{10.3815,5.3311},{11.5039,5.8922},{12.6262,6.7340},{13.7485,7.2951},{14.8708,7.8563},{15.9932,8.4175},{17.1155,8.9786},{18.2378,9.5398},{19.3602,10.6621},{20.4825,11.2233},{21.6048,12.0650},{22.7271,12.6262},{23.8495,13.1874},{24.9718,13.1874},{26.0941,13.7485},{27.2165,13.7485},{28.3388,13.7485},{29.4611,13.7485},{30.5834,13.7485},{31.7058,13.1874},{32.8281,13.1874},{33.9504,12.9068},{35.0728,12.6262},{36.1951,12.6262},{37.3174,12.3456},{38.4397,11.7844},{39.5621,11.2233},{40.6844,10.6621},{41.8067,10.1010},{42.9291,9.5398},{44.0514,8.9786},{45.1737,8.6980},{46.2960,7.8563},{47.4184,7.2951},{48.5407,6.7340},{49.6630,6.4534}}
spinescens[2]={{0.2806,8.9786},{1.4029,9.2592},{2.5252,9.2592},{3.6476,9.8204},{4.7699,9.8204},{5.8922,9.8204},{7.0146,9.8204},{8.1369,10.3815},{9.2592,11.2233},{10.3815,12.3456},{11.5039,13.4679},{12.6262,14.5903},{13.7485,15.9932},{14.8708,16.8349},{15.9932,18.5184},{17.1155,19.9213},{18.2378,21.3242},{19.3602,22.4466},{20.4825,23.2883},{21.6048,24.1301},{22.7271,24.9718},{23.8495,25.5330},{24.9718,25.5330},{26.0941,24.9718},{27.2165,24.4106},{28.3388,23.8495},{29.4611,23.2883},{30.5834,22.7271},{31.7058,22.1660},{32.8281,21.3242},{33.9504,20.4825},{35.0728,19.6407},{36.1951,19.0796},{37.3174,17.9572},{38.4397,16.8349},{39.5621,15.7126},{40.6844,14.0291},{41.8067,12.3456},{42.9291,11.2233},{44.0514,9.8204},{45.1737,8.6980},{46.2960,8.1369}}
spinescens[3]={{4.7699,12.6262},{5.8922,12.9068},{7.0146,13.4679},{8.1369,14.5903},{9.2592,15.7126},{10.3815,16.2738},{11.5039,17.6767},{12.6262,19.0796},{13.7485,20.2019},{14.8708,21.3242},{15.9932,22.4466},{17.1155,23.5689},{18.2378,24.6912},{19.3602,25.8135},{20.4825,26.6553},{21.6048,27.4970},{22.7271,28.3388},{23.8495,29.1805},{24.9718,29.4611},{26.0941,30.0223},{27.2165,30.8640},{28.3388,31.1446},{29.4611,31.7058},{30.5834,31.9864},{31.7058,32.2669},{32.8281,32.2669},{33.9504,31.9864},{35.0728,31.7058},{36.1951,31.1446},{37.3174,30.3029},{38.4397,29.7417},{39.5621,28.6194},{40.6844,28.0582},{41.8067,27.4970},{42.9291,26.3747},{44.0514,25.8135},{45.1737,24.6912},{46.2960,24.1301},{47.4184,23.0077},{48.5407,21.8854},{49.6630,20.7631},{50.7853,19.3602},{51.9077,17.9572},{53.0300,15.9932},{54.1523,15.1514},{55.2747,14.5903},{56.3970,14.0291},{57.5193,13.1874},{58.6416,12.6262}}
spinescens[4]={{1.4029,10.3815},{2.5252,10.6621},{3.6476,11.2233},{4.7699,12.3456},{5.8922,13.4679},{7.0146,14.0291},{8.1369,14.5903},{9.2592,15.4320},{10.3815,15.9932},{11.5039,17.3961},{12.6262,18.5184},{13.7485,19.0796},{14.8708,20.2019},{15.9932,21.3242},{17.1155,22.4466},{18.2378,23.0077},{19.3602,23.8495},{20.4825,24.4106},{21.6048,24.9718},{22.7271,25.2524},{23.8495,25.5330},{24.9718,26.0941},{26.0941,25.5330},{27.2165,25.5330},{28.3388,25.5330},{29.4611,25.5330},{30.5834,26.0941},{31.7058,26.0941},{32.8281,25.5330},{33.9504,25.2524},{35.0728,24.6912},{36.1951,23.5689},{37.3174,22.4466},{38.4397,21.3242},{39.5621,20.7631},{40.6844,19.6407},{41.8067,18.5184},{42.9291,16.8349},{44.0514,15.7126},{45.1737,14.0291},{46.2960,12.9068},{47.4184,11.7844},{48.5407,10.9427},{49.6630,10.3815}}
spinescens[5]={{8.1369,5.8922},{9.2592,6.7340},{10.3815,7.2951},{11.5039,8.4175},{12.6262,8.9786},{13.7485,10.1010},{14.8708,11.2233},{15.9932,12.3456},{17.1155,14.0291},{18.2378,14.3097},{19.3602,14.5903},{20.4825,14.8708},{21.6048,14.8708},{22.7271,15.7126},{23.8495,16.8349},{24.9718,17.6767},{26.0941,18.5184},{27.2165,19.6407},{28.3388,20.2019},{29.4611,20.7631},{30.5834,21.3242},{31.7058,21.6048},{32.8281,22.1660},{33.9504,22.1660},{35.0728,22.7271},{36.1951,22.7271},{37.3174,22.7271},{38.4397,22.7271},{39.5621,22.7271},{40.6844,22.4466},{41.8067,22.1660},{42.9291,21.6048},{44.0514,21.3242},{45.1737,21.0437},{46.2960,20.4825},{47.4184,19.6407},{48.5407,19.0796},{49.6630,18.2378},{50.7853,17.6767},{51.9077,16.8349},{53.0300,15.7126},{54.1523,14.5903},{55.2747,13.4679},{56.3970,12.3456},{57.5193,11.2233},{58.6416,10.1010},{59.7640,8.9786},{60.8863,8.4175},{62.0086,8.1369},{63.1310,7.5757}}
spinescens[6]={{3.6476,7.5757},{4.7699,7.8563},{5.8922,8.6980},{7.0146,9.5398},{8.1369,10.6621},{9.2592,11.7844},{10.3815,12.9068},{11.5039,14.0291},{12.6262,14.5903},{13.7485,15.7126},{14.8708,16.2738},{15.9932,17.3961},{17.1155,18.5184},{18.2378,18.5184},{19.3602,18.2378},{20.4825,18.2378},{21.6048,18.2378},{22.7271,18.2378},{23.8495,18.7990},{24.9718,18.7990},{26.0941,19.3602},{27.2165,20.2019},{28.3388,20.7631},{29.4611,21.0437},{30.5834,21.6048},{31.7058,21.6048},{32.8281,21.6048},{33.9504,21.6048},{35.0728,21.6048},{36.1951,21.6048},{37.3174,21.6048},{38.4397,21.6048},{39.5621,21.0437},{40.6844,21.0437},{41.8067,20.4825},{42.9291,20.2019},{44.0514,19.9213},{45.1737,19.6407},{46.2960,19.0796},{47.4184,18.5184},{48.5407,17.3961},{49.6630,16.2738},{50.7853,15.7126},{51.9077,14.5903},{53.0300,13.4679},{54.1523,12.3456},{55.2747,11.5039},{56.3970,10.3815},{57.5193,9.8204},{58.6416,9.2592},{59.7640,8.9786},{60.8863,8.6980}}

(*Poner aqui la familia*)
var=spinescens;
(* el número de curvas estudiadas/semillas de la familia*)
m=6;
For[j=1,j≤m,j++,
puntos=var[j];
(* pru es el spline cubico
d es el vector que contiene la información de pru 
n es el número de trozos del spline, es decir, el número de puntos menos 1
*)
pru=SplineFit[puntos,Cubic];
d=InputForm[pru];
n=IntegerPart[d[[1,2,2]]];
spli=ParametricPlot[pru[t],{t,0,n},PlotRange->All, DisplayFunction->Identity];
(* Con los coeficientes se calculan los polinomios cúbicos que dan el spline. Habrá n polinomios para la primera coordenada (X) y otros tantos para la segunda
*)
For[i=0,i<n,i++,
SX[i,t_]=d[[1,4,i+1,1,1]]+d[[1,4,i+1,1,2]]*(t-i)+d[[1,4,i+1,1,3]]*(t-i)^2+d[[1,4,i+1,1,4]]*(t-i)^3];
For[i=0,i<n,i++,SY[i,t_]=d[[1,4,i+1,2,1]]+d[[1,4,i+1,2,2]]*(t-i)+d[[1,4,i+1,2,3]]*(t-i)^2+d[[1,4,i+1,2,4]]*(t-i)^3];
(* Derivada primera y segunda de cada trozo del spline
*)
For[i=0,i<n,i++,DX[i,t_]=d[[1,4,i+1,1,2]]+2*d[[1,4,i+1,1,3]]*(t-i)+3*d[[1,4,i+1,1,4]]*(t-i)^2];
For[i=0,i<n,i++,DY[i,t_]=d[[1,4,i+1,2,2]]+2*d[[1,4,i+1,2,3]]*(t-i)+3*d[[1,4,i+1,2,4]]*(t-i)^2];
For[i=0,i<n,i++,DDX[i,t_]=2*d[[1,4,i+1,1,3]]+6*d[[1,4,i+1,1,4]]*(t-i)];
For[i=0,i<n,i++,DDY[i,t_]=2*d[[1,4,i+1,2,3]]+6*d[[1,4,i+1,2,4]]*(t-i)];
(*Cálculo de la curvatura a trozos
*)
For[i=0,i<n,i++,cur[i,t_]=(DX[i,t]DDY[i,t]-DDX[i,t]DY[i,t])/((DX[i,t])^2+(DY[i,t])^2)^(3/2)];
(* Dibujos de los trozos de curvatura 
*)
For[i=0,i<n,i++,dib[i]=Plot[cur[i,t],{t,i,i+1},DisplayFunction->Identity]];
pp=dib[0];
For[i=1,i<n,i++,pp=Show[pp,dib[i],PlotRange->All]];
(* bez es el spline de Bezier
b es el vector que contiene la información de bez
nb es el número de trozos del spline, es decir, el número de puntos menos 1 (Como nb=n, no se define nuevo)
*)
bez=SplineFit[puntos,Bezier];
b=InputForm[bez];
(*Calculo de las funciones de Bezier para la coordenada X y la Y; la curva es la imagen del intervalo[0,1]*)
For[i=1,i≤n+1,i++,BX[t_]=bez[[4,1,1]]/.#1->t];
For[i=1,i≤n+1,i++,BY[t_]=bez[[4,1,2]]/.#1->t];
dibez[j]=ParametricPlot[{BX[t],BY[t]},{t,0,1},PlotRange->All,DisplayFunction->Identity];
curbe[j,t_]=(BX'[t]BY''[t]-BX''[t]BY'[t])/((BX'[t])^2+(BY'[t])^2)^(3/2);
cur[j]=Abs[NMinValue[{curbe[j,t],0.30<t<0.7},t]];
vmed[j]=Abs[NIntegrate[curbe[j,t],{t,0.3,0.7}]/0.4];
dibcurbe[j]=Plot[{curbe[j,t],-vmed[j]},{t,0,1},PlotRange->{-0.2,0.1},DisplayFunction->Identity,PlotStyle->{RGBColor[0,1,0],{Dotted,Dashed}}];
total=GraphicsGrid[{
(*{spli,pp},*)
{dibez[j],dibcurbe[j]}},Spacings->0.3];
Print[var,j, ":", "Máximo valor de la curvatura (en valor absoluto)",cur[j]];
Print[total]
]
Print["Curvaturas de cada ", var ,":  ", Table[cur[j],{j,1,6}]];
Print["Valores medios de curvatura de cada ", var ," en [0.3,0.7]:  ", Table[vmed[j],{j,1,6}]];

Show[Table[dibcurbe[j],{j,1,6}]]

Table[NMinValue[{curbe[j,t],0.25<t<0.75},t],{j,1,6}]

dibcurbe[j],{j,1,6}

curbe[j,t_]
dibcurbe[j]