ft_read_virgo_map.m

用matlab做的Frequency domain INterferomEter Simulation SoftwarE,在国外网站下的,还在研究中,不会用,有会用的回复我

%
%--------------------------------------------------------------------------
% function [map]=FT_read_virgo_map(filename)
%
% A function for Matlab which reads a mirror surface map in the standard
% virgo format.
%
% filename: name of virgo file
% map:      the surface map (variable must be in the standard virgo
%           structure)
%
% Charlotte Bond  10.07.2008
%--------------------------------------------------------------------------
%

function [map] = FT_read_virgo_map(filename)

  baseid='read_virgo_map';
    
  % Open file for reading
  [FID, result]=fopen(filename,'r');
    
  if (FID==-1)
    msgid=[baseid,':fileopen'];
    error(msgid,result);
  end
  
  name=filename;
  
  % Standard virgo format maps have type 'phase' (they store surface
  % heights)
  type=0;

  % In virgo fromat maps both (reflected and transmitted) light fields are 
  % affected 
  field=0;
  
  % In the standard virgo format the measurements are given in nanometers
  scaling=1e-09;
  
  % Now read the data
  [M,count]=fscanf(FID,'%g ');
  
  if (count<8)
    msgid=[baseid,':fileread'];
    result='could not read map header';
    error(msgid,result);
  end
  
  % The first two numbers are the grid dimensions in unit steps
  cols=M(1);
  rows=M(2);
  
  % The next two numbers give the optical centre of the grid in coordinates
  x0=M(3);
  y0=M(4);
  
  % The next two numbers are zeros and are not used
  
  % The next number is the zero point number (for all points not on the
  % optical surface) given in nanometers
  zero=M(7);
  
  % The next number is the step size given in meters.  The grid is square
  % so xstep and ystep are the same
  xstep=M(8);
  ystep=M(8);
  
  if(length(M)<9)
    % read dummy line of '=================='
    tmp=fgets(FID);
    [M_tmp,count_tmp]=fscanf(FID,'%g ');
    M=[M;M_tmp];
    count=count+count_tmp;
  end
  
  % Close the file
  fclose(FID);
    
  if (count<x0*y0)
    msgid=[baseid,':fileread'];
    result='number of data points not consistent with header information';
    error(msgid,result);
  end
  
  % Then the data follows.  Sometime the number before the data is 0.0 in 
  % which case we have to ignore this point.  If it is the first data point
  % it should be equal to the zero point.  In the case where the first 
  % number is the first data point, the count of the data differs for 
  % different maps.  If the last number is 0.0 then the last data point is
  % at M(count-1).  Otherwise it is M(count).
  
  if M(count)==0.0
      data_count=count-1;
  else
      data_count=count;
  end
  
  % We skip the first point if it is not data
  if M(9)==zero
      tmp_data=M(9:data_count);
  else
      tmp_data=M(10:count);
  end
  
  % Reshape the data into a 2D matrix:
  % - the reshape creates a matrix from the vector.
  % - however, it has x and y exchanged, therefor we take
  %   the transpose. 
  % - further, the virgo map format starts
  %   in the top left corner. More intuitive for using the
  %   data in Matglab or Finesse would be to start in the
  %   lower left corner. Thus we applu 'flipud'
  
  data=flipud(reshape(tmp_data,cols,rows)');
  
  
  % In the standard virgo format the data not on the surface is set to
  % an arbitrary 'zero point' value. We need to set all the points which 
  % are the same numerical value as the zero point to zero.
  zerodata=find(data==zero);
  data(zerodata)=0;
  
  % Create map structure from the data
  map=FT_new_surface_map(name,data,type,field,x0,y0,xstep,ystep,scaling);