C snippet for WSPR broadcasts

This is a short update about my WSPR beacon experiments. I am trying to set up a small standalone micro controller based WSPR beacon. As there are no details available about the WSPR protocol, I have analysed the original Fortran code of the original WSPR application written by K1JT.

The following code snippet generates an .au file of a WSPR broadcast containing the specified CALL, GRID and EIRP:

#include <stdio.h>

void main(){

   // pack call in n1
   static char c[6]=" K1JT "; //call
   // all must be in uppercase
   // char at position 3 must be last digit of call prefix
   // left and right idented with space characters
   unsigned int n1;
   n1=(c[0]>='0'&&c[0]<='9'?c[0]-'0':c[0]==' '?36:c[0]-'A'+10);
   n1=36*n1+(c[1]>='0'&&c[1]<='9'?c[1]-'0':c[1]==' '?36:c[1]-'A'+10);
   n1=10*n1+c[2]-'0';
   n1=27*n1+(c[3]==' '?26:c[3]-'A');
   n1=27*n1+(c[4]==' '?26:c[4]-'A');
   n1=27*n1+(c[5]==' '?26:c[5]-'A');

   // pack grid, dbm in n2
   static int lat  = 40;   //degrees latitude north
   static int lon = -74;   //degrees longitude east
   int ng =(-lon+180+1)*90 + lat;
   static int dbm = 33;    //EIRP dBm
   // must be in set {0,3,7,10,13,17,20,23,27,30,33,37,40,43,47,50,53,57,60}
   int n2=ng<<7|0x40|dbm;

   // pack n1,n2 into 50 bits
   char packed[11] = {n1>>20, n1>>12, n1>>4, (n1&0x0f)<<4|(n2>>18)&0x0f, n2>>10, n2>>2, (n2&0x03)<<6, 0, 0, 0, 0};
//   printf("packed=%d %d %d %d %d %d %d %d %d %d %d\n",packed[0],packed[1],packed[2],packed[3],packed[4],packed[5],packed[6],packed[7],packed[8],packed[9],packed[10]);

   // Convolutional encoder for a K=32, r=1/2 code.
   // Layland-Lushbaugh polynomials for a K=32, r=1/2 convolutional code,
   // and 8-bit parity lookup table.
   const unsigned char partab[256] = { 
      0,1,1,0,1,0,0,1,1,0,0,1,0,1,1,0,
      1,0,0,1,0,1,1,0,0,1,1,0,1,0,0,1,
      1,0,0,1,0,1,1,0,0,1,1,0,1,0,0,1,
      0,1,1,0,1,0,0,1,1,0,0,1,0,1,1,0,
      1,0,0,1,0,1,1,0,0,1,1,0,1,0,0,1,
      0,1,1,0,1,0,0,1,1,0,0,1,0,1,1,0,
      0,1,1,0,1,0,0,1,1,0,0,1,0,1,1,0,
      1,0,0,1,0,1,1,0,0,1,1,0,1,0,0,1,
      1,0,0,1,0,1,1,0,0,1,1,0,1,0,0,1,
      0,1,1,0,1,0,0,1,1,0,0,1,0,1,1,0,
      0,1,1,0,1,0,0,1,1,0,0,1,0,1,1,0,
      1,0,0,1,0,1,1,0,0,1,1,0,1,0,0,1,
      0,1,1,0,1,0,0,1,1,0,0,1,0,1,1,0,
      1,0,0,1,0,1,1,0,0,1,1,0,1,0,0,1,
      1,0,0,1,0,1,1,0,0,1,1,0,1,0,0,1,
      0,1,1,0,1,0,0,1,1,0,0,1,0,1,1,0 };
   unsigned char symbol[176];
   int nstate = 0;
   int k = 0;
   int i,j;
   for(j=0;j!=sizeof(packed);j++){
      for(i=7;i>=0;i--){
         int i4 = packed[j];
         if(i4 < 0) i4 = i4 + 0x100;
         nstate = (nstate<<1) | ((i4>>i) & 0x01);
         int n = nstate & 0xf2d05351;
         n=n^(n>>16);
         symbol[k] = partab[(n^(n>>8)) & 0xff];
         k++;
         n=nstate & 0xe4613c47;
         n=n^(n>>16);
         symbol[k] = partab[(n^(n>>8)) & 0xff];
         k++;
      }
   }

   // Interleave symbols
   // j0 is iteration of bit swapped bytes 0 to 255 with value < 162
   const unsigned char j0[162] = {
      0x00,0x80,0x40,0x20,0xa0,0x60,0x10,0x90,
      0x50,0x30,0x70,0x08,0x88,0x48,0x28,0x68,
      0x18,0x98,0x58,0x38,0x78,0x04,0x84,0x44,
      0x24,0x64,0x14,0x94,0x54,0x34,0x74,0x0c,
      0x8c,0x4c,0x2c,0x6c,0x1c,0x9c,0x5c,0x3c,
      0x7c,0x02,0x82,0x42,0x22,0x62,0x12,0x92,
      0x52,0x32,0x72,0x0a,0x8a,0x4a,0x2a,0x6a,
      0x1a,0x9a,0x5a,0x3a,0x7a,0x06,0x86,0x46,
      0x26,0x66,0x16,0x96,0x56,0x36,0x76,0x0e,
      0x8e,0x4e,0x2e,0x6e,0x1e,0x9e,0x5e,0x3e,
      0x7e,0x01,0x81,0x41,0x21,0xa1,0x61,0x11,
      0x91,0x51,0x31,0x71,0x09,0x89,0x49,0x29,
      0x69,0x19,0x99,0x59,0x39,0x79,0x05,0x85,
      0x45,0x25,0x65,0x15,0x95,0x55,0x35,0x75,
      0x0d,0x8d,0x4d,0x2d,0x6d,0x1d,0x9d,0x5d,
      0x3d,0x7d,0x03,0x83,0x43,0x23,0x63,0x13,
      0x93,0x53,0x33,0x73,0x0b,0x8b,0x4b,0x2b,
      0x6b,0x1b,0x9b,0x5b,0x3b,0x7b,0x07,0x87,
      0x47,0x27,0x67,0x17,0x97,0x57,0x37,0x77,
      0x0f,0x8f,0x4f,0x2f,0x6f,0x1f,0x9f,0x5f,
      0x3f,0x7f };
   unsigned char symbol2[162];
   for(i=0;i!=162;i++)
      symbol2[j0[i]]=symbol[i]; //interleave

   // Sync vector 162 bits
   const unsigned char npr3[162] = {
      1,1,0,0,0,0,0,0,1,0,0,0,1,1,1,0,
      0,0,1,0,0,1,0,1,1,1,1,0,0,0,0,0,
      0,0,1,0,0,1,0,1,0,0,0,0,0,0,1,0,
      1,1,0,0,1,1,0,1,0,0,0,1,1,0,1,0,
      0,0,0,1,1,0,1,0,1,0,1,0,1,0,0,1,
      0,0,1,0,1,1,0,0,0,1,1,0,1,0,1,0,
      0,0,1,0,0,0,0,0,1,0,0,1,0,0,1,1,
      1,0,1,1,0,0,1,1,0,1,0,0,0,1,1,1,
      0,0,0,0,0,1,0,1,0,0,1,1,0,0,0,0,
      0,0,0,1,1,0,1,0,1,1,0,0,0,1,1,0,
      0,0 };
   for(i=0;i!=162;i++){
      symbol2[i] = npr3[i] | symbol2[i]<<1;
   }

   for(i=0;i!=162;i++){
//      printf("%u", symbol2[i]&0x02);
//      printf("symbol[%u]=%u\n", i, symbol[i]);
   }

   // Modulate
   FILE* pf = fopen("a.au", "wb");
   fputc(0x2e,pf); //.au magic number
   fputc(0x73,pf);
   fputc(0x6e,pf);
   fputc(0x64,pf);
   fputc(0x00,pf); //offset
   fputc(0x00,pf);
   fputc(0x00,pf);
   fputc(24,pf);
   fputc(0xff,pf); //size
   fputc(0xff,pf);
   fputc(0xff,pf);
   fputc(0xff,pf);
   fputc(0x00,pf); //encoding
   fputc(0x00,pf);
   fputc(0x00,pf);
   fputc(0x03,pf);
   fputc(0x00,pf); //samplerate
   fputc(0x00,pf);
   fputc(0x2e,pf);
   fputc(0xe0,pf);
   fputc(0x00,pf); //channels
   fputc(0x00,pf);
   fputc(0x00,pf);
   fputc(0x01,pf);

   const int nmax=120*12000;
   double tsymbol=8192.0/12000.0;
   double dt=1.0/12000.0;
   double f0=1500;
   double dfgen=12000.0/8192.0;
   double t=0;
   double dphi;
   double phi=0;
   double dj0=0;
   for(i=0; i!=nmax; i++){
      t=t+dt;
      j=(int)(t/tsymbol)+1;
      if(j!=dj0){
         double f=f0+dfgen*(symbol2[j]-1.5);
         dj0=j;
         dphi=2*3.14*dt*f;
      }
      phi=phi+dphi;
      int wave=32767.0*sin(phi);
      fputc((wave>>8)&0xff,pf);
      fputc(wave&0xff,pf);
   }
   fclose(pf);
}