```    //  NYD.C
//  Author: Peter Meyer
//  This program calculates the statistics for NYD occurrences
//  (dates in March/April) over a 4000-year period.
//  Last mod.: 1999-03-29

#include <stdio.h>

#include <mplc.h>

#define N 15

long freq[64][3];
int n[3] = { 201, 1001, 4001 };
double sum_pc[3] = { 0.0, 0.0, 0.0 };

//  Common Era    ------ Meyer-Palmen Solilunar Calendar ------     Julian day
//    0-03-22    069-06-01-01 Wednesday, Aristarchus  1, 069-06      1,721,141

/*-----------*/
void main(void)
{
int i, j, diy, month, day;
long year, ce_yr;
double pc;
MPDate mpdt;

for ( i=0; i<64; i++ )
for ( j=0; j<3; j++ )
freq[i][j] = 0;

mpdt.jdn = 1721141L;
//  This is the JD of the MP NYD in Year 0 CE.
jdn_to_mpdate(&mpdt);
//  mpdt now holds the MP era date correspondng to the MP NYD in Year 0 CE.

//  For each year in the range 0 through 4000 CE
//  get the date CE of the MP NYD.
for ( ce_yr=0; ce_yr<=4000; ce_yr++ )
{
//  Convert MP era date to CE date.
mpdate_to_ce_date(&mpdt,&year,&month,&day);

diy = day + 31*(month==4);
//  diy = days since last day of February.

freq[diy][2]++;
if ( ce_yr >=1500 && ce_yr <= 2500 )
freq[diy][1]++;
if ( ce_yr >= 1900 && ce_yr <= 2100 )
freq[diy][0]++;

//  Advance mpdt to next NYD.
if ( ++mpdt.year > 6840 )
{
mpdt.year = 1;
mpdt.cycle++;
}
//  mpdt.month and mpdt.day both remain 1.
}

//  Now print the results.

printf("\n%*s%*s%*s%*s",
N,"",N,"2000CE+/-100",N,"2000CE+/-500",N,"2000CE+/-2000");
printf("\n%*s%*s%*s%*s\n\n",
N,"",N,"201 NYDs",N,"1001 NYDs",N,"4001 NYDs");

for ( diy=1; diy<=61; diy++ )
{
if ( freq[diy][2] )
{
printf("%*s",N-3,( diy <= 31 ? "March" : "April" ));
printf("%3d", ( diy <= 31 ? diy : diy-31 ));
for ( i=0; i<3; i++ )
{
pc = (((double)freq[diy][i])*100)/n[i];
sum_pc[i] += pc;
printf("%*ld, %4.*f%%",N-7,freq[diy][i],(i?2:1),pc);
}
printf("\n");
}
}

printf("%*s%*s%*s%*s\n",
N,"",N,"-----",N,"-----",N,"-----");

printf("%*s",N,"");
for ( i=0; i<3; i++ )
printf("%*.2f%%",N-1,sum_pc[i]);
printf("\n");
}
```