clear();
clc();
mode(-1);
Path = get_absolute_file_path("PolarPattern_example.sce");
chdir(Path);
funcprot(0);
getd(Path)
files="";
// data are placed in a text file in csv format
// Phi, from 0 to 360°, is on the first line  
// Theta,from 0 to +180° is the first column  
// the other data is a matrix corresponding to the modulus[Theta,Phi]

function [Phi,Theta,data]=ReadData(fil,sep)
    // read data in csv file
    M=csvRead(fil,sep,'.','double'); 
    Theta=M(2:$,1);  // Theta 1st column
    Phi=M(1,2:$);    // Phi 1st line
    data=M(2:$,2:$); // matrix of modulus[Theta,Phi]
endfunction;

//files=uigetfile(["*.csv"],Path, "Choose a file name", %t);
//[Ph,Thet,Gph]= ReadData(files(1),",");
//[Ph,Thet,Gth]= ReadData(files(2),",");
 [Ph,Thet,Gph]= ReadData(Path+'Gphi.csv',",");
 [Ph,Thet,Gth]= ReadData(Path+'Gtheta.csv',",");
 
  Gtot=10*log10(Gth.^2+Gph.^2);
  Gph=20*log10(Gph);
  Gth=20*log10(Gth);
  scmax=max(Gtot);
  scale=25;
  fg=gcf();
  fg.color_map       = jetcolormap(128);
  fg.figure_size=[700,900]
  subplot(3,2,1);
  polar3D(Ph,Thet,Gph,scale,scmax,'G Phi');
  subplot(3,2,2);
  polar3D(Ph,Thet,Gth,scale,scmax,'G Theta');
  subplot(3,2,3);
  polar3D(Ph,Thet,Gtot,scale,scmax,'G Total');
  
  // Etheta Ephi2 @ Theta=90°  horizontal plane
  i_Thet=find(round(Thet)==90); // search index line corresponding to Theta=90°
  subplot(3,2,4);
  polar2D(Ph, [Gth(i_Thet,:);Gph(i_Thet,:);Gtot(i_Thet,:)], ...
          scale,scmax, 'cut Theta=90° ',...
          ["V pol","H pol","Total"],"xy");

    // Etheta1 Ephi1 @ phi=0° vertical plane
    //retreive two columns in matrix corresponding to the Phi cut
    Et=CutPhi(Gth,0);       // H polar
    Ep=CutPhi(Gph,0);       // V polar
    Ett=CutPhi(Gtot,0);     // V+H  polar
    Thx=[Thet; 180+Thet(2:$)]'; // built Theta over 360°
  subplot(3,2,5);
  polar2D(Thx, [Et;Ep;Ett], ...
          scale,scmax, 'cut Phi=0° ',...
          ["V pol","H pol","Total"],"xz");

    // Etheta1 Ephi1 @ phi=90° vertical plane
    //retreive two columns in matrix corresponding to the Phi cut
    Et=CutPhi(Gth,90);   // H polar
    Ep=CutPhi(Gph,90);   // V polar
    Ett=CutPhi(Gtot,90); // V+H  polar
    Thx=[Thet; 180+Thet(2:$)]';
  subplot(3,2,6);
  polar2D(Thx,[Et;Ep;Ett], scale,scmax, 'cut Phi=90° ',...
          ["V pol","H pol","Total"],"yz");
// draw 3axis reference
  triax(.42,.63);
