Monday, February 21, 2011

perspective opcv



#include "stdafx.h"
#include <cv.h>
#include <highgui.h>
#include <cxcore.h>
#include <math.h>
using namespace cv;
int main(int argc, char** argv) {
  CvPoint2D32f srcQuad[4], dstQuad[4];
  CvMat*       warp_matrix = cvCreateMat(3,3,CV_32FC1);
  IplImage     *src1, *dst1, *dst2, *dst3, *dst4, *dst5, *dst7, *dst8;
 src1=cvLoadImage("C:\\ENEL\\OpenCV2.1\\samples\\c\\airplane.jpg");
IplImage *src = cvCreateImage( cvSize(800 ,700), src1->depth, src1->nChannels );
IplImage *dst6 = cvCreateImage( cvSize(800 ,700), src1->depth, src1->nChannels );
IplImage *checkerfloor = cvCreateImage( cvSize(200,600), src1->depth, src1->nChannels );
IplImage *endwall = cvCreateImage( cvSize(800,700), src1->depth, src1->nChannels );
//Draw a checkerfloor using OpenCV drawing tools, question 1a
//Mostly using cvFillPoly()
cvZero(checkerfloor); //make image all black
cvNamedWindow( "checkerfloor", CV_WINDOW_AUTOSIZE);
//generate coodinates for white squares
CvPoint  square1[]={0,0,  0,100,  100,100,  100,0};//draw topleft white square, array of coordinates
CvPoint  square2[]={100,100,  100,200,  200,200,  200,100};//draw second from top right white square, and so on etc...
CvPoint  square3[]={0,200,  0,300,  100,300,  100,200};
CvPoint  square4[]={100,300,  100,400,  200,400,  200,300};
CvPoint  square5[]={0,400,  0,500,  100,500,  100,400};
CvPoint  square6[]={100,500,  100,600,  200,600,  200,500};
CvPoint* squareArr[6]={square1,square2,square3,square4,square5,square6};// array of array
int      nSquarePts[6]={4,4,4,4,4,4};//number of points for each square array
int      nSquares=6;// number of squares
int      isPolyLineClosed=1;//use this for cvPolyLine
int      lineWidth=1; //use this for cvPolyLine
//cvPolyLine(checkerfloor,squareArr,nSquarePts,nSquares,isPolyLineClosed,cvScalar(0,255,255),lineWidth);
//Draw a set of filled polygons, squares in this case:
cvFillPoly(checkerfloor,squareArr,nSquarePts,nSquares,cvScalar(255,255,255), 8);
//cvLine(checkerfloor, cvPoint(100,100), cvPoint(200,200), CV_RGB(255,255,255),1,8);
//cvLine(checkerfloor, cvPoint(200,200), cvPoint(500,300), CV_RGB(255,255,255),1,8);
cvShowImage("checkerfloor", checkerfloor);
//ends question 1a
//Draw a endwall using OpenCV drawing tools, question 1
//Mostly using cvFillPoly()
CvPoint  wall[]={0,0,  0,700,  800,700,  800,0};
CvPoint* wallArr[1]={wall};// array of array
int      nWallPts[4]={4};//number of points for each square array
int      nWalls=1;// number of squares
//Draw a set of filled polygons, squares in this case:
cvFillPoly(endwall,wallArr,nWallPts,nWalls,cvScalar(123,56,78), 8);
cvNamedWindow( "endwall", CV_WINDOW_AUTOSIZE);
cvShowImage("endwall", endwall);
//end endwall

//This is part is for question 1c)
cvResize(src1, src);
cvResize(checkerfloor, dst6); //clone dst* to intialize
    dst1 = cvCloneImage(src);
    dst1->origin = src->origin;
    dst2 = cvCloneImage(src);
    dst2->origin = src->origin;
    dst3 = cvCloneImage(src);
    dst3->origin = src->origin;
    dst4 = cvCloneImage(src);
    dst4->origin = src->origin;
 dst5 = cvCloneImage(src);
 cvZero(dst5);//This  sets all values in all channels of the array to zero
 dst7 = cvCloneImage(src);
 dst7->origin = src->origin;
 cvZero(dst7);
//dst for for wall
 dst8 = cvCloneImage(src);
 dst8->origin = src->origin;
 cvZero(dst8);
//srcQuad is parameter for cvGetPerspectiveTransform()
    //----------------------------------------------
    srcQuad[0].x = 0;               //src Top left
    srcQuad[0].y = 0;
    srcQuad[1].x = src->width - 1;  //src Top right
    srcQuad[1].y = 0;
    srcQuad[2].x = 0;               //src Bottom left
    srcQuad[2].y = src->height - 1;
    srcQuad[3].x = src->width - 1;  //src Bot right
    srcQuad[3].y = src->height - 1;
//pic1: left side nearest
    dstQuad[0].x = 0; //dst Top left
    dstQuad[0].y = 0;
    dstQuad[1].x = 100;  //dst Top right
    dstQuad[1].y = 100;
    dstQuad[2].x = 0;  //dst Bottom left
    dstQuad[2].y = 700;
    dstQuad[3].x = 100;  //dst Bot right
    dstQuad[3].y = 700-400/3;
    cvGetPerspectiveTransform(
      srcQuad,
      dstQuad,
      warp_matrix
    );
    cvWarpPerspective( src, dst1, warp_matrix );
//pic2: left side farest
    dstQuad[0].x = 100; //dst Top left
    dstQuad[0].y = 100;
    dstQuad[1].x = 300*6/7;  //dst Top right
    dstQuad[1].y = 300*6/7;
    dstQuad[2].x = 100;  //dst Bottom left
    dstQuad[2].y = 700-400/3;
    dstQuad[3].x = 300*6/7;  //dst Bot right
    dstQuad[3].y = 300*6/7+100;
    cvGetPerspectiveTransform(
      srcQuad,
      dstQuad,
      warp_matrix
    );
    cvWarpPerspective( src, dst2, warp_matrix );
//pic3 right side farest
    dstQuad[0].x = 300+500/7; //dst Top left
    dstQuad[0].y = 300*6/7;
    dstQuad[1].x = 300+500*2/3;  //dst Top right
    dstQuad[1].y = 100;
    dstQuad[2].x = 300+500/7;  //dst Bottom left
    dstQuad[2].y = 300*6/7+100;
    dstQuad[3].x = 300+500*2/3;  //dst Bot right
    dstQuad[3].y = 700-400/3;
    cvGetPerspectiveTransform(
      srcQuad,
      dstQuad,
      warp_matrix
    );
    cvWarpPerspective( src, dst3, warp_matrix );
//pic4 right side nearest
    dstQuad[0].x = 300+500*2/3; //dst Top left
    dstQuad[0].y = 100;
    dstQuad[1].x = 800;  //dst Top right
    dstQuad[1].y = 0;
    dstQuad[2].x = 300+500*2/3;  //dst Bottom left
    dstQuad[2].y = 700-400/3;
    dstQuad[3].x = 800;  //dst Bot right
    dstQuad[3].y = 700;
    cvGetPerspectiveTransform(
      srcQuad,
      dstQuad,
      warp_matrix
    );
    cvWarpPerspective( src, dst4, warp_matrix );
//-------------------------------------------------
//added checker floor, Added march 6, 2011
 dstQuad[0].x = 300*6/7; //dst Top left
    dstQuad[0].y = 300*6/7+100;
    dstQuad[1].x =300+500/7 ;  //dst Top right
    dstQuad[1].y = 300*6/7+100;
    dstQuad[2].x = 0;  //dst Bottom left
    dstQuad[2].y = 700;
    dstQuad[3].x = 800;  //dst Bot right
    dstQuad[3].y = 700;
    cvGetPerspectiveTransform(
      srcQuad,
      dstQuad,
      warp_matrix
    );
    cvWarpPerspective( dst6, dst7, warp_matrix );

//added endwall, Added march 8, 2011
 dstQuad[0].x = 300*6/7; //dst Top left
    dstQuad[0].y = 300*6/7;
    dstQuad[1].x =300+500/7 ;  //dst Top right
    dstQuad[1].y = 300*6/7;
    dstQuad[2].x = 300*6/7;  //dst Bottom left
    dstQuad[2].y = 300*6/7+100;
    dstQuad[3].x = 300+500/7;  //dst Bot right
    dstQuad[3].y =300*6/7+100 ;
    cvGetPerspectiveTransform(
      srcQuad,
      dstQuad,
      warp_matrix
    );
    cvWarpPerspective( endwall, dst8, warp_matrix );
//add all pics on dst5
// end question 1c
 cvAddWeighted(dst1,1,dst2,1,0,dst5);
 cvAddWeighted(dst3,1,dst5,1,0,dst5);
 cvAddWeighted(dst4,1,dst5,1,0,dst5);
 cvAddWeighted(dst7,1,dst5,1,0,dst5);
 cvAddWeighted(dst8,1,dst5,1,0,dst5);
    cvNamedWindow( "Perspective_Warp", 1 );
      cvShowImage( "Perspective_Warp", dst5 );
int p[3]; p[0] = CV_IMWRITE_JPEG_QUALITY; p[1] = 100; p[2] = 0; cvSaveImage("out2.jpg", dst5, p);
     cvWaitKey(0);
   
 cvReleaseImage(&dst1);
 cvReleaseImage(&dst2);
 cvReleaseImage(&dst3);
 cvReleaseImage(&dst4);
 cvReleaseImage(&dst5);
 cvReleaseImage(&dst6);
 cvReleaseImage(&dst7);
 cvReleaseImage(&dst8);
 cvReleaseImage(&src);
 cvReleaseImage(&src1);
 cvReleaseImage(&checkerfloor);
 cvDestroyWindow("checkerfloor");
 cvDestroyWindow("endwall");
    cvReleaseMat(&warp_matrix);
 cvDestroyWindow( "Perspective_Warp" );
    return 0;
  }
---------------------------------------------------------------------------------

Sunday, February 13, 2011

cvSobel

#include "stdafx.h"
#include <cv.h>
#include <highgui.h>
#include <cxcore.h>
#include <math.h>
using namespace cv;
int g_slider_position=0,high=0,low=0;
CvCapture* capture         = NULL;
void onTrackbarSlide(int pos) {}
//------------------------------------------
int main( int argc, char* argv[] ) {
    capture = cvCreateFileCapture( "C:\\Users\\abc\\Desktop\\MVI_0035.avi" );
    if(!capture){
        return -1;
    }
 CvSize size = cvSize(
       (int)cvGetCaptureProperty( capture, CV_CAP_PROP_FRAME_WIDTH),
       (int)cvGetCaptureProperty( capture, CV_CAP_PROP_FRAME_HEIGHT)
    );
    IplImage* frame = cvCreateImage(size,IPL_DEPTH_8U,3);
    IplImage* grey = cvCreateImage(size,IPL_DEPTH_8U,1);
 IplImage* edge = cvCreateImage(size, IPL_DEPTH_8U, 1);
    IplImage* dirsx = cvCreateImage(size, IPL_DEPTH_16S, 1);
    IplImage* dirsy = cvCreateImage(size, IPL_DEPTH_16S, 1);
    IplImage* dirsxy = cvCreateImage(size, IPL_DEPTH_16S, 1);
    IplImage* dirx = cvCreateImage(size, IPL_DEPTH_8U, 1);
 IplImage* diry = cvCreateImage(size, IPL_DEPTH_8U, 1);
    IplImage* dirxy = cvCreateImage(size, IPL_DEPTH_8U, 1);
    while( (frame=cvQueryFrame(capture)) != NULL ) {
  cvCvtColor(frame,grey,CV_RGB2GRAY); 
  cvNamedWindow( "gray", CV_WINDOW_AUTOSIZE );
  cvNamedWindow( "canny", CV_WINDOW_AUTOSIZE );
  cvNamedWindow( "sobel", CV_WINDOW_AUTOSIZE );
  cvCreateTrackbar("T","gray",&g_slider_position,1000,onTrackbarSlide);
  cvCreateTrackbar("high","canny",&high,1000,onTrackbarSlide);
  cvCreateTrackbar("low","canny",&low,1000,onTrackbarSlide);
  cvSmooth(grey,grey,CV_GAUSSIAN,13);
  cvCanny(grey, edge, low, high);
  //-----------------------
  cvSobel(grey,dirsx,1,0, 3);
  cvConvertScaleAbs(dirsx,dirx,1,0);
  cvSobel(grey,dirsy,0,1, 3);
  cvConvertScaleAbs(dirsy,diry,1,0);
  cvAddWeighted(dirsx,0.707,dirsy,0.707,0,dirsxy);
     cvConvertScaleAbs(dirsxy,dirxy,1,0);
  //-------------------
  cvShowImage( "gray", grey );
  cvShowImage( "canny", edge );
  cvShowImage("sobel", dirxy);
  cvWaitKey(1001-g_slider_position);
    }
    cvReleaseImage( &frame );
 cvReleaseImage( &grey );
    cvReleaseCapture( &capture );
 cvDestroyWindow( "gray" );
 cvDestroyWindow( "canny" );
 cvDestroyWindow( "sobel" );
    return(0);
}

Saturday, February 12, 2011

opcv header

http://www710.univ-lyon1.fr/~bouakaz/OpenCV-0.9.5/docs/ref/OpenCVRef_ImageProcessing.htm
http://www710.univ-lyon1.fr/~bouakaz/OpenCV-0.9.5/docs/ref/OpenCVRef_BasicFuncs.htm#decl_CvPoint2D32f

polyfit polyval matlab

polyfit - Polynomial curve fitting

Syntax

p = polyfit(x,y,n)
[p,S] = polyfit(x,y,n)
[p,S,mu] = polyfit(x,y,n)

Description

p = polyfit(x,y,n) finds the coefficients of a polynomial p(x) of degree n that fits the data, p(x(i)) to y(i), in a least squares sense. The result p is a row vector of length n+1 containing the polynomial coefficients in descending powers:
[p,S] = polyfit(x,y,n) returns the polynomial coefficients p and a structure S for use with polyval to obtain error estimates or predictions. Structure S contains fields R, df, and normr, for the triangular factor from a QR decomposition of the Vandermonde matrix of x, the degrees of freedom, and the norm of the residuals, respectively. If the data y are random, an estimate of the covariance matrix of p is (Rinv*Rinv')*normr^2/df, where Rinv is the inverse of R. If the errors in the data y are independent normal with constant variance, polyval produces error bounds that contain at least 50% of the predictions.
[p,S,mu] = polyfit(x,y,n) finds the coefficients of a polynomial in
where and . mu is the two-element vector [μ12]. This centering and scaling transformation improves the numerical properties of both the polynomial and the fitting algorithm.

Examples

This example involves fitting the error function, erf(x), by a polynomial in x. This is a risky project because erf(x) is a bounded function, while polynomials are unbounded, so the fit might not be very good.
First generate a vector of x points, equally spaced in the interval [0, 2.5]; then evaluate erf(x) at those points.
x = (0: 0.1: 2.5)';
y = erf(x); 
The coefficients in the approximating polynomial of degree 6 are
p = polyfit(x,y,6)

p =

 0.0084  -0.0983   0.4217   -0.7435  0.1471   1.1064  0.0004
There are seven coefficients and the polynomial is . To see how good the fit is, evaluate the polynomial at the data points with:
f = polyval(p,x);
A table showing the data, fit, and error is
table = [x y f y-f]
 
 table =

    0          0          0.0004    -0.0004
    0.1000     0.1125     0.1119     0.0006
    0.2000     0.2227     0.2223     0.0004
    0.3000     0.3286     0.3287    -0.0001
    0.4000     0.4284     0.4288    -0.0004
    ...
    2.1000     0.9970     0.9969     0.0001
    2.2000     0.9981     0.9982    -0.0001
    2.3000     0.9989     0.9991    -0.0003
    2.4000     0.9993     0.9995    -0.0002
    2.5000     0.9996     0.9994     0.0002
So, on this interval, the fit is good to between three and four digits. Beyond this interval the graph shows that the polynomial behavior takes over and the approximation quickly deteriorates.
x = (0: 0.1: 5)';
y = erf(x);
f = polyval(p,x);
plot(x,y,'o',x,f,'-')
axis([0  5  0  2])

Algorithm

The polyfit MATLAB file forms the Vandermonde matrix, V, whose elements are powers of x.
It then uses the backslash operator, \, to solve the least squares problem Vpy.
You can modify the MATLAB file to use other functions of x as the basis functions.

See Also

poly, polyval, roots, lscov, cov

Friday, February 11, 2011

Laws for thermocouples

Law of homogeneous material
A thermoelectric current cannot be sustained in a circuit of a single homogeneous material by the application of heat alone, regardless of how it might vary in cross section. In other words, temperature changes in the wiring between the input and output do not affect the output voltage, provided all wires are made of the same materials as the thermocouple. No current flows in the circuit made of a single metal by the application of heat alone.

[edit] Law of intermediate materials

The algebraic sum of the thermoelectric emfs in a circuit composed of any number of dissimilar materials is zero if all of the junctions are at a uniform temperature. So If a third metal is inserted in either wire and if the two new junctions are at the same temperature, there will be no net voltage generated by the new metal.

[edit] Law of successive or intermediate temperatures

If two dissimilar homogeneous materials produce thermal emf1 when the junctions are at T1 and T2 and produce thermal emf2 when the junctions are at T2 and T3 , the emf generated when the junctions are at T1 and T3 will be emf1 + emf2,provided T1<T2<T3.

Tuesday, February 8, 2011

ghdl tb_df2j

ghdl -s df2j.vhd fracmult.vhd lprc1.vhd sregisterbn.vhd tb_df2j.vhd tb_lprc1.vhd
ghdl -i df2j.vhd fracmult.vhd lprc1.vhd sregisterbn.vhd tb_df2j.vhd tb_lprc1.vhd
ghdl -m --warn-error tb_df2j
ghdl -r tb_df2j --stop-time=400ns --vcd=tb_df2j.vcd
C:\ENEL\gtkwave\gtkwave -S script.txt tb_df2j.vcd

-----------------------------------------------------------------------
library ieee;
use ieee.std_logic_1164.all;
use ieee.numeric_std.all;
entity tb_df2j is
end tb_df2j;
architecture Behavioral of tb_df2j is
   component df2j
      Port(
         clk :  in  std_logic;
         enable :  in  std_logic;
         reset :  in  std_logic;
         a1m1 :  in  signed(7 downto 0);
         a2 :  in  signed(7 downto 0);
         fin :  in  signed(15 downto 0);
         fout :  out  signed(15 downto 0)
      );
   end component;
   constant clkPERIOD : time := 10 ns;
   signal W_clk : std_logic;
   signal W_enable : std_logic;
   signal W_reset : std_logic;
   signal W_a1m1 : signed(7 downto 0);
   signal W_a2 : signed(7 downto 0);
   signal W_fin : signed(15 downto 0);
   signal W_fout : signed(15 downto 0);
   signal Trp : std_logic;
begin
   DUT: df2j
      port map(
         clk => W_clk,
         enable => W_enable,
         reset => W_reset,
         a1m1 => W_a1m1,
         a2 => W_a2,
         fin => W_fin,
         fout => W_fout
      );
   defclk: process
      begin
         W_clk <= '0';
         wait for clkPERIOD/2;
         W_clk <= '1';
         wait for clkPERIOD/2;
      end process;
   deftrp: process
     begin
       Trp <= '0';
       wait for 3*clkPERIOD;
       Trp <= '1';
       wait for clkPERIOD;
     end process;
     W_enable <= Trp;
    
   STIMULI: process
      begin
        W_reset <= '0';
         W_a2 <= to_signed(123,8);
         W_a1m1 <= to_signed(111,8);
         W_fin <= to_signed(0,16);
       
        wait for clkPERIOD + 3 ns;
        W_reset <= '1';
        wait for clkPERIOD;        
        W_reset <= '0';
        wait for 4*clkPERIOD;                
        W_fin <= to_signed(1000,16);       
        wait for 2*clkPERIOD;
        W_fin <= to_signed(0,16); 
         assert false report "end of test" severity error;
         wait;
      end process STIMULI;
end Behavioral;
------------------------------------------------------------

Sunday, February 6, 2011

opcv code 4

higui function

http://www.cognotics.com/opencv/docs/1.0/ref/opencvref_highgui.htm#decl_CvCapture
------------------------------------------------------------
grey->edge detector
#include "stdafx.h"
#include <cv.h>
#include <highgui.h>
#include <cxcore.h>
#include <math.h>
using namespace cv;
int main()
{
 IplImage* img = cvLoadImage("C:\\ENEL\\OpenCV2.1\\samples\\c\\airplane.jpg");
 IplImage* edge = cvCreateImage(cvGetSize(img), IPL_DEPTH_8U, 1);
 cvCvtColor(img, edge, CV_BGR2GRAY);
 cvCanny(edge, edge, 0.3, 0.8);
 cvShowImage( "Example1", edge );
 cvWaitKey(0);
 cvReleaseImage(&img);
 cvReleaseImage(&edge);
}
----------------------------------------------------------
video writer
int main( int argc, char* argv[] ) {
    CvCapture* capture = 0;
    capture = cvCreateFileCapture( "C:\\ENEL\\OpenCV2.1\\samples\\c\\tree.avi" );
    if(!capture){
        return -1;
    }
    double fps = cvGetCaptureProperty (
        capture,
        CV_CAP_PROP_FPS
    );
 CvSize size = cvSize(
       (int)cvGetCaptureProperty( capture, CV_CAP_PROP_FRAME_WIDTH),
       (int)cvGetCaptureProperty( capture, CV_CAP_PROP_FRAME_HEIGHT)
    );
    CvVideoWriter *writer = cvCreateVideoWriter(
        "C:\\Users\\abc\\Desktop\\a.avi",
        CV_FOURCC('M','J','P','G'),
        fps,
        size
    );
    IplImage* frame = cvCreateImage(
        size,
        IPL_DEPTH_8U,
        3
    );
    while( (frame=cvQueryFrame(capture)) != NULL ) {
        cvWriteFrame( writer, frame );
    }
    cvReleaseVideoWriter( &writer );
    cvReleaseImage( &frame );
    cvReleaseCapture( &capture );
    return(0);
}
--------------------------------------------------------------
better edge detector

http://dasl.mem.drexel.edu/~noahKuntz/openCVTut5.html

sim city 4 mod

http://www.youtube.com/watch?v=Un0LJLLh_ms&feature=related

http://www.simtropolis.com/stex/

http://sc4devotion.com/

login: sanitoryplayer   

EA keygen

http://d.namipan.com/downfile/EAKG-_www.jinghua8.com.rar/1e1a1553b27851b63c0cf09c5cd9ce89556eeb00df360500

sanitorysocial8

Thursday, February 3, 2011

Standard deviation

Consider a population consisting of the following eight values:

    2,\  4,\  4,\  4,\  5,\  5,\  7,\  9
These eight data points have the mean (average) of 5:

    \frac{2 + 4 + 4 + 4 + 5 + 5 + 7 + 9}{8} = 5
To calculate the population standard deviation, first compute the difference of each data point from the mean, and square the result of each:

    \begin{array}{lll}
    (2-5)^2 = (-3)^2 = 9   &&  (5-5)^2 = 0^2 = 0 \\
    (4-5)^2 = (-1)^2 = 1  &&  (5-5)^2 = 0^2 = 0 \\
    (4-5)^2 = (-1)^2 = 1  &&  (7-5)^2 = 2^2 = 4 \\
    (4-5)^2 = (-1)^2 = 1  &&  (9-5)^2 = 4^2 = 16
    \end{array}
Next compute the average of these values, and take the square root:

    \sqrt{ \frac{9 + 1 + 1 + 1 + 0 + 0 + 4 + 16}{8} } = 2
This quantity is the population standard deviation; it is equal to the square root of the variance. The formula is valid only if the eight values we began with form the complete population. If they instead were a random sample, drawn from some larger,  parent” population, then we should have used 7 (n − 1) instead of 8 (n) in the denominator of the last formula, and then the quantity thus obtained would have been called the sample standard deviation. See the section Estimation below for more details

------------------------------------------------------------------------
#include <math.h>
#include <string>
#include <iostream>
#include <fstream>
using namespace std;
int main()
{
 int xn=0, N = 30;
 double x[N],m=0,s=0;
 string filename,outfilename;
   ifstream instream;
  
   cout<<"file"<<endl;
 cin>>filename;                        
 instream.open(filename.c_str());
 cout <<"successful"<<endl;

while(xn != N)
{
 instream >> x[xn];
// cout << "x"<<xn << " = " << x[xn] << endl;
 xn++;
 }

 for(xn=0;xn!=N;xn++)
 {
  m=m+x[xn];
  }
 
 m=m/N;
 cout<<"median = "<<m<<endl;

 for(xn=0;xn!=N;xn++)
 {
  s=s+(x[xn]-m)*(x[xn]-m);
  }

 s=sqrt(s/N);
 cout<<"standard deviation = "<<s<<endl;

/*cout<<"outfile"<<endl;
cin>>outfilename;
ofstream outfile(outfilename.c_str());
for(xn=0;xn!=N;xn++)
{
 outfile<<xn<<" "<<x[xn];
 if(xn%2 != 0){outfile<<endl;}
 else{outfile<<" ";}
 }*/
cout<<"done"<<endl;
return 0;
 }

 
-------------------------------------------