function imagemodel(source_position_index)
%%
% /*
%  *      cref.c is for use in generating room-simulation data for plotting.
%  *              It calculates Cartesian coordinates for each reflection
%  *              as well as the delay and dB attenuation on each.
%  *	In contrast to reflect, refplot, image18, etc., it takes source
%  *		position in Cartesian coordinates.
%  *
%  *      input parameters:
%  *              xl          - listener's x coordinate (meters)
%  *              yl          - listener's y coordinate
%  *              zl          - listener's z coordinate
%  *              xr          - width of simulated room (meters)
%  *              yr          - length of simulated room
%  *              zr          - height of simulated room
%  *              xs          - source's x coordinate (meters)
%  *              ys          - source's y coordinate
%  *              zs          - source's z coordinate
%  *          optional:
%  *              repeated sets of az, el and r for a repeated sequence
%  *                  of locations
%  *
%  *      output parameters:
%  *              x       - x coordinate of image-model reflection
%  *              y       - y coordinate of image-model reflection
%  *              z       - z coordinate of image-model reflection
%  *              az      - azimuth incidence angles (degrees)
%  *              el      - elevation incidence angles (degrees)
%  *              delay   - reflection latencies (sec)
%  *              db      - dB down from reference 1 m away

%%
    %load basic values for image model
    load 'basic_values.mat';
    
    %load room data: Receiving and Source position, Room size
    load 'room_data.mat';

%     if (narg < 10) {
%         fprintf(stderr,"narg = %d\n", narg);
%         fprintf(stderr,"USAGE: %s xl yl zl xr yr zr xs ys zs\n", argv[0]);
%         exit(-1);
%     }

%%
%  setup the time window */

    sr = basic_values.SR;      % Sampling rate
    clength = basic_values.CLENGTH;
    window = clength / (2. * sr);
    spm = basic_values.SPM;
    
%     lc = cell2mat(room_data.l_coord(1,:));
%     sc = cell2mat(room_data.s_coord(1,:));
%     rs = cell2mat(room_data.room_size(1,:));
    lc = room_data(5,:);
    sc = room_data(source_position_index,:);
    rs = room_data(4,:);
        
%  Get listener location in meters and convert it to seconds	*/
%  atof function : Convert string to double
    xl = lc(1) * spm;        % SPM converts meters to seconds */
    yl = lc(2) * spm;
    zl = lc(3) * spm;

%  Get room size in meters and convert it to seconds	*/
    
    xr = rs(1) * spm;
    yr = rs(2) * spm;
    zr = rs(3) * spm;
    
%  Get source location in meters and convert it to seconds	*/

    xref(1) = sc(1) * spm;
    yref(1) = sc(2) * spm;
    zref(1) = sc(3) * spm;
    
    xs = xref(1);
    ys = yref(1);
    zs = zref(1);
  
%%    
%      Calculate coordinates of the image model virtual sources      */

    ir = 1;
    for ic = -basic_values.N_DEEP:basic_values.N_DEEP
        xc(ir) = cvs(ic,xs,xr);
        yc(ir) = cvs(ic,ys,yr);
        zc(ir) = cvs(ic,zs,zr);
        
        ir = ir+1;
    end
    
    nr = ir;
    ir = 1;
     
    for ix = 1:nr-1
        for iy = 1:nr-1
            for iz = 1:nr-1
                xref(ir) = xc(ix);
                yref(ir) = yc(iy);
                zref(ir) = zc(iz);
        	    reford(ir) = abs(ix-basic_values.N_DEEP) + abs(iy-basic_values.N_DEEP) + abs(iz-basic_values.N_DEEP);
        
                % Get the az, el, and distance by using x,y,z
                [az(ir) el(ir) delay(ir)] = ctos( xref(ir) - xl, yref(ir) - yl, zref(ir) - zl);
                %ctos(xref(ir)-xl,yref(ir)-yl,zref(ir)-zl,az(ir),el(ir),delay(ir));
                
                db(ir) = 20.*log10(basic_values.REFDIST/delay(ir))-delay(ir)*basic_values.ABSORB-(reford(ir))*basic_values.BOUNCE;
                
                if(delay(ir) < basic_values.MAXDEL) 
%                     if confirm == 1
%                         result_IM_cref2d = zeros (ir, 8, 3);
%                         confirm = 0;
%                     end
                    %result_IM_cref2d(ir,:, source_position_index) = [basic_values.MPS*xref(ir), basic_values.MPS*yref(ir), basic_values.MPS*zref(ir),az(ir),el(ir),delay(ir),db(ir),reford(ir)];
                    temp(ir,:) = [basic_values.MPS*xref(ir), basic_values.MPS*yref(ir), basic_values.MPS*zref(ir),az(ir),el(ir),delay(ir),db(ir),reford(ir)];
                    ir = ir+1;
                end
            end
        end       
    end
    
%                         % Save all data to result_ISM
%                     result_ISM {source_position_index} = temp;
%     save result_ISM result_ISM;
    if source_position_index == 1
        result_ISM{source_position_index} = {temp};
        save result_ISM result_ISM;
    elseif source_position_index == 2
        load result_ISM;
        %temp_1 = result_ISM{source_position_index-1};
        result_ISM{source_position_index} = {temp};
%         result_IM_cref2d(:,:,source_position_index-1) = temp(:,:,:);
%         result_IM_cref2d(:,:,source_position_index) = result_IM_cref2d(:,:,source_position_index-1);
        save result_ISM result_ISM;
    elseif source_position_index == 3
        load result_ISM;
        result_ISM{source_position_index} = {temp};
        save result_ISM result_ISM;
    end

