Ok I dont think I provided enough data last I asked this and as a result I am still stuck so here goes;

I have data in 4 matrices describing a wave going from one medium to another, the incident wave transmission and reflection coefficients (F1), the exiting waves transmittance and reflection coefficients (F2) and the propagating waves in each medium (P1 and P2). All this boils down to 4 matrices with data.

The propagating wave and incident wave vectors are summed and square rooted and the Real parts defined as cos_phi[i] for increasing k which the 4 matrices before are defined by.

I plotted the results as cos_bloch_phase1 using pointplot and the sequence of values of k versus cos_phi.

Then displayed cos_bloch_phase1,plot... etc

What I am struggling to figure out is how to export the data that is plotted in cos_bloch_phase1. I have copied and pasted my code in its entirety of which a further plot is included but can be ignored as they are constructed similarly;

restart;
with(plots): with(linalg):

matrix elements for the first interface

incident parameters
theta0:=89*Pi/180; k_min:=0.0; k_max:=25;

number_of_points:=1000;
step:=(k_max-k_min)/number_of_points;
                              1000
                         0.02500000000

parameters of a binary crystal
n1:=1.46; n2:=4.6; d1:=0.1; d2:=0.1;
                              1.46
                              4.6
                              0.1
                              0.1


for i from 0 by 1 to number_of_points do
k[i]:=k_min+i*step;

theta1:=evalf(arcsin(n0*sin(theta0)/n1)); theta2:=evalf(arcsin(n0*sin(theta0)/n2));
q1:=k[i]*n1*cos(theta1); q2:=k[i]*n2*cos(theta2);


F1:=Matrix([[(q1+q2)/(2*q1),(q1-q2)/(2*q1)],[(q1-q2)/(2*q1),(q1+q2)/(2*q1)]]);
F2:=Matrix([[(q2+q1)/(2*q2),(q2-q1)/(2*q2)],[(q2-q1)/(2*q2),(q2+q1)/(2*q2)]]);
P1:=Matrix([[exp(-I*q1*d1),0],[0,exp(I*q1*d1)]]);
P2:=Matrix([[exp(-I*q2*d2),0],[0,exp(I*q2*d2)]]);
M_period:=multiply(P1, F1, P2, F2);

cos_phi[i]:=Re(1/2*(M_period[1,1]+M_period[2,2]));
end:


cos_bloch_phase1:=pointplot([seq([k[i],cos_phi[i]],i=0..number_of_points)],connect=true):

display(cos_bloch_phase1,plot([1,-1],x=0..k_max,color=[red,red]));



for i from 0 by 1 to number_of_points do
k[i]:=k_min+i*step;

theta1:=evalf(arcsin(n0*sin(theta0)/n1)); theta2:=evalf(arcsin(n0*sin(theta0)/n2));
q1:=k[i]*n1*cos(theta1); q2:=k[i]*n2*cos(theta2);

F1p:=Matrix([[(q1/n1^2+q2/n2^2)/(2*q1/n1^2),(q1/n1^2-q2/n2^2)/(2*q1/n1^2)],[(q1/n1^2-q2/n2^2)/(2*q1/n1^2),(q1/n1^2+q2/n2^2)/(2*q1/n1^2)]]);
F2p:=Matrix([[(q2/n2^2+q1/n1^2)/(2*q2/n2^2),(q2/n2^2-q1/n1^2)/(2*q2/n2^2)],[(q2/n2^2-q1/n1^2)/(2*q2/n2^2),(q2/n2^2+q1/n1^2)/(2*q2/n2^2)]]);
P1p:=Matrix([[exp(-I*q1*d1),0],[0,exp(I*q1*d1)]]);
P2p:=Matrix([[exp(-I*q2*d2),0],[0,exp(I*q2*d2)]]);
Mp_period:=multiply(P1p, F1p, P2p, F2p);


cos_phi2[i]:=Re(1/2*(Mp_period[1,1]+Mp_period[2,2]));
end:



cos_bloch_phase2:=pointplot([seq([k[i],cos_phi2[i]],i=0..number_of_points)],connect=true):

display(cos_bloch_phase2,plot([1,-1],x=0..k_max,color=[red,red]));


display(cos_bloch_phase2,cos_bloch_phase1,plot([1,-1],x=0..k_max+5,color=[red,red]));

I have a plot of a sinusoidal function(cos_phi) which changes with incrementing values of variable k. When this function goes above 1 or below -1 I would like to have the range of k for which it occurs outputed somehow. Would anybody know how to do this?

I have theoretically 3(could eventually be more) layers with an incident wave with a wave equation for that wave.

It refracts into the 2nd layer from the first and now has a 2nd wave equation, then from the 2nd into the 3rd layer with a 3rd wave equation.

All the wave equations are of the form, Psi(z) = A_1psi_1(z) + B_1psi_2(z); this is just a general solution where psi_1&2 are linearly independant solutions that make up the general equation above and A_1 and B_1 are constant coefficients that would be A_2,B_2 and A_3,B_3 for the 2nd and 3rd layers respectively.

Transfer matrix method gives A_1,B_1 in terms of A_2,B_2(as it transfers from layer 1 to 2 they equate under boundary conditions so you can solve the simultaneous equations for results). You create a matrix of these results and multiply it with the respective matrix of the 2nd layer to 3rd layer to give you the overall transfer matrix from one side of the system to the other.

I think something to do with transfer function but not sure how to use it or set up the problem. 

Thanks in advance for any pointers.

Hey all new to Maplesoft my question is this;