/* add_srtm_mask.c - takes an existing srtm mask and updates it using a sub grid i
                     1 - land
                     0 - water */

#include <stdio.h>
#include <math.h>
#include <stdlib.h>

unsigned char maskset[8] = { 128, 64, 32, 16, 8, 4, 2, 1};

main(argc, argv)
int	argc;
char	**argv;
{
	FILE	*fp1, *fp2;
	int	i, j, ntot, nln, nlt, nln0, nlt0, nmask = 116640000, wet_flag = 0;
        int	ioff, joff, kk, kkoff, ia;
        unsigned char *a;
        short int *b;

/* read the command line */

        if (argc != 8) {
                fprintf(stderr,"usage:  add_strm_mask <srtm.bit > <new_mask> nlt0 nln0 nln nlt wet_flag \n");
                fprintf(stderr,"                      srtm.bit  - existing srtm_bitmask \n");
                fprintf(stderr,"                      new_mask  - new 30 second pixel registered mask to add \n");
                fprintf(stderr,"                      nln0, nlt0 - upper left corner\n");
                fprintf(stderr,"                      nln, nlt   - number of columns, rows (4800 6000 or antatctica 7200 3600 \n");
                fprintf(stderr,"                      wetflag    - elevation value signifying wet area \n");
                exit(-1);
        }

        nln0  = atoi(argv[3]);
        nlt0  = atoi(argv[4]);
        nln   = atoi(argv[5]);
        nlt   = atoi(argv[6]);
        wet_flag = atoi(argv[7]);
        ntot  = nln*nlt;

/*  allocate memory */

        if(( a = (unsigned char *) malloc(nmask)) == NULL){
                printf("sorry, couldn't allocate memory for mask\n");
                exit(-1);
        }

        if(( b = (short int *) malloc(ntot*2)) == NULL){
                printf("sorry, couldn't allocate memory for bytemask\n");
                exit(-1);
        }

/*  open the files */

        if ( (fp1 = fopen(argv[1], "r+")) == NULL) {

/*  create a zero mask because the files does not exist */

                fprintf(stderr,"create a zero mask \n");
  		fp1 = fopen(argv[1], "w");
                for (i = 0; i < nmask; i++) {
        	        a[i]= 0;
                }
		if ((fwrite(a, 1, nmask, fp1)) != nmask) {
			fprintf(stderr,"add_srtm_mask:  ERROR during write of %s\n", argv[1]);
                }
		fclose(fp1);
		exit(0);
	}

        if ( (fp2 = fopen(argv[2], "r")) == NULL) {
                fprintf(stderr,"add_srtm_mask:  ERROR.  Cannot open w %s\n", argv[2]);
                exit(-1);
        }

/*  read the bit mask */

	if ((fread(a, 1, nmask, fp1)) != nmask) {
		fprintf(stderr,"add_srtm_mask:  ERROR during read of %s\n", argv[1]);
		exit(-1);
	}
        fprintf(stderr, "read the bit mask \n");

/*  read the topo array */

	if ((fread(b, 2, ntot, fp2)) != ntot) {
		fprintf(stderr,"add_srtm_mask:  ERROR during read of %s\n", argv[2]);
		exit(-1);
	}
        fclose(fp2);
        fprintf(stderr, "read the topo file \n");

/*  convert bits to bytes */

	ioff = 120*(90-nlt0);
        joff = -120*(-180-nln0);
    	for (i = 0; i < nlt; i++) {
		for(j = 0; j < nln; j++) {
			kk = nln*i + j;
			kkoff=120*360*(ioff+i)+(joff+j);
			ia = kkoff/8;
			if(b[kk] == wet_flag) {
				a[ia] = a[ia] - (a[ia] & maskset[kkoff%8]);
			}
			else {
				a[ia] = a[ia] | maskset[kkoff%8];
			}
		}
	}

/*  write the bit mask */

        rewind(fp1);
	if ((fwrite(a, 1, nmask, fp1)) != nmask) {
		fprintf(stderr,"add_srtm_mask:  ERROR during write of %s\n", argv[1]);
		exit(-1);
	}
	fclose(fp1);

	exit(0);
}