function [R,Z,f] = plotRZGauss(N,beta,a,b,delta,tau0)
%function [R,Z,f] = plotRZGauss(N,beta,a,b,tau0,L) nos va a generar una carta de puntos (R,Z) en el instante que yo le diga
%utilizando como parametros:

%   N:         Numero de puntos a utilizar
%   beta:      Parametro beta que modela el comportamiento gaussiano en la direccion radial que vamos a utilizar.
%   a:         Hace referencia a la altura maxima del cilindro
%   b:         Hace referencia a la distancia radial maxima del sistema que vamos a estudiar
%   delta:     Parametro beta que modela el comportamiento gaussiano en la direccion Z que vamos a utilizar.
%   tau0:      Instante en el que voy a considerar el campo escalar de temperaturas.

%En este programa estamos suponieno que el perfil va a ser gaussiano en R, controlado a traves del parametro beta, y tambien gaussiano en Z,
%controlado a traves del parametro delta. Beta y delta se especificaran mas adelante, una vez estudiado el comportamiento del laser y como la
%optica podria modificar el resultado.

    [Am,betam] = coefAs(N,beta,b);
    
%Temperaturas o niveles de temperatura que estamos buscando

    %v = [1.0,0.9,0.8,0.7,0.6,0.5,0.4,0.3,0.2,0.1,0];
    v = linspace(0,1.0,20);
    
%Veamos los parametros asociados al experimento. Las temperaturas se miden en ºC

%    Tint = 60;
%    Text = -150;
%    Ts = 0;
%Veamos el dato referente al radio de la celula

%    R0 = 20*10^-6;  %Globulo rojo humano
%    d = 5*10^-6;      %Anchura o espesor de un portaobjetos

%Los parametros que requiere el modelo para considerar la evolucion lineal en la variable Z sera:

%    ns = (Ts - Text)/(Tint - Text);
%    H1 = d/R0;
    
%Queremos obtener los valores de la funcion de temperatura para distintos valores de R y Z, barriendo un amplio
%margen del disco

    B=0;i=1;k=1;
    
    Btau = 1/sqrt(1+4*tau0/delta);
    
    for Z=0:0.2:a
        for R=-b:0.2:b
        
            for n=1:1:N
                B = B + Am(n)*besselJ(0,betam(n)*R)*exp(-1*tau0*betam(n)^2);
            end
        
            f(k,i) = B*(erf(Z/(4*tau0)) -Btau*exp(-(Z^2)*(Btau^2)/delta)*erf(Btau*Z/sqrt(4*tau0)));
            i=i+1;
            B = 0;
        end
        k = k + 1;
        i=1;
    end
 
 %Una vez generada la carta devolvamos tambien los datos de R y Z por si se requiere alguna representacion adimensional.
 
    R=-b:0.2:b;
    Z=0:0.2:a;
    
    %plot(R,H1*ones(size(R)));hold on;
    %h = contour(R,Z,f,v);
    %clabel(h);title('Carta R,Z en el el instante tau=0.07 supuesto perfil gaussiano en Z');
    
    title('Curvas Asociadas a un perfil gaussiano tanto en R como en Z');
    subplot(1,2,1);
    pcolor(R,Z,real(f))
    shading interp  %Para borrar las lineas de la rejilla
    title('Curva con pseudocolor de la Temperatura adimensional');
    caxis([0,1])
    colorbar
    
    subplot(2,2,2);
    h = contour(R,Z,f,v);
    clabel(h);title('Carta de contornos R,Z')
    
    subplot(2,2,4);
    h = surf(R,Z,real(f));
    shading interp;
    xlabel('R');ylabel('Z');zlabel('Temperatura adimensional');title('Curva tridimensional de la Tª');
    
    