SHARE
TWEET

Untitled

a guest Feb 19th, 2019 102 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. ' frequency modulator & demodulator for fm_radio_template.reason
  2. ' by veezay
  3. '
  4. ' FModem.exe (fmdemod9.exe)
  5. ' latest version 2018-11-27 (prev 2018-07-31)
  6. '
  7. ' new in version:
  8. ' - doesn't overwrite output file fmodem_out.wav, instead uses running numbering from 001 to 999
  9. '
  10. '
  11. ' stuff todo:
  12. ' - add timestamp of start & finishing ETA of processing
  13. ' - select new RANDOMIZE seed from input audio
  14. ' - rewrite WAV header as a TYPE
  15. '
  16. '
  17. ' stuff to try out:
  18. ' - figure out a nice carrier frequency? or dunno? maybe the current is just fine? depends on stuff.
  19. ' - try different carrier waveforms such as triangle (cleanest zero crossing calculation?)
  20.  
  21.  
  22. RANDOMIZE TIMER
  23.  
  24.  
  25. DIM lutsize AS INTEGER, lutsizeoffset AS INTEGER, lfolutsize AS INTEGER, audiolutsize AS INTEGER, lutindex AS _UNSIGNED INTEGER, lutindextemp AS DOUBLE
  26. DIM palutratio AS _UNSIGNED LONG, lfopalutratio AS _UNSIGNED LONG, audiopalutratio AS _UNSIGNED LONG
  27. DIM outputsample(2) AS INTEGER, previousoutputsample(2) AS INTEGER
  28. DIM sinelookup(1024, 1024) AS INTEGER
  29. DIM modulatoroutput AS INTEGER
  30. DIM samplerate AS _UNSIGNED LONG, newsample AS INTEGER, samplex AS _UNSIGNED LONG
  31. DIM carrier AS _UNSIGNED LONG, carrieros AS _UNSIGNED LONG
  32. DIM currentsample AS _UNSIGNED _INTEGER64, s AS _UNSIGNED _INTEGER64, position AS _UNSIGNED _INTEGER64, lastposition AS _UNSIGNED _INTEGER64
  33. DIM sampleoffset AS _UNSIGNED _INTEGER64, outputoffset AS _UNSIGNED _INTEGER64, samplebuffer AS _UNSIGNED LONG, outputbuffer AS _UNSIGNED LONG, n AS _UNSIGNED LONG
  34. DIM dummy AS _BYTE, dummyint AS INTEGER
  35. DIM downsampled AS LONG, writedownsampled AS _BIT, lastchunk AS _BIT
  36. lastchunk = 0
  37.  
  38. oversampling = 4 ' 8
  39. samplerate = 88200 * oversampling
  40. carrier = 45059 '35059 '* 2 + 1
  41. carrieros = carrier \ oversampling
  42. writedownsampled = 1
  43.  
  44. ' lookup table sizes
  45. lfolutsize = 64 ' not in use currently
  46. audiolutsize = 256
  47. lfopalutratio = 2 ^ 32 / lfolutsize ' phase accumulator size / lookup table size ratio
  48.  
  49. ' how many input samples to load per chunk + same for output
  50. samplebuffer = 2 ^ 16
  51. outputbuffer = 2 ^ 16
  52. outputfileextra = 256 ' compensation for me being a shitty coder
  53.  
  54. ' make lookup tables of oscillator
  55. GOSUB generatelookuptable
  56.  
  57.  
  58. ' define input file
  59. IF COMMAND$ = "" THEN
  60.     inputfiles$ = "C:\files\reason\exports\fm_in.wav" ' probably should replace this with some error message telling user to plz drag&drop .wav on .exe
  61. ELSE
  62.     inputfiles$ = COMMAND$
  63.     inputfile = 1
  64.     FOR i = 1 TO LEN(inputfiles$)
  65.         char$ = MID$(inputfiles$, i, 1)
  66.         IF char$ = " " THEN
  67.             inputfile = inputfile + 1
  68.         ELSE
  69.             inputfile$(inputfile) = inputfile$(inputfile) + char$
  70.         END IF
  71.     NEXT i
  72.     numberofinputs = inputfile
  73.     FOR i = 1 TO numberofinputs
  74.         LOCATE i: PRINT " input file "; trim$(i); ": "; inputfile$(i)
  75.     NEXT i
  76. END IF
  77.  
  78. ' start timing the execution
  79. DIM time AS DOUBLE
  80. time = TIMER(0.001)
  81.  
  82. ' get input file specifics
  83. DIM maxchannels AS INTEGER, numberofchannels(numberofinputs) AS INTEGER, samplelength(numberofinputs) AS _UNSIGNED _INTEGER64, sampletime(numberofinputs)
  84. DIM inputfilenumber(numberofinputs)
  85.  
  86. FOR i = 1 TO numberofinputs
  87.     inputfilenumber(i) = FREEFILE
  88.     OPEN inputfile$(i) FOR BINARY AS inputfilenumber(i)
  89.     GET inputfilenumber(i), 23, numberofchannels(i)
  90.  
  91.     IF maxchannels < numberofchannels(i) THEN maxchannels = numberofchannels(i)
  92.  
  93.     samplelength(i) = (LOF(inputfilenumber(i)) - 44) / 2 / numberofchannels(i)
  94.     sampletime(i) = samplelength(i) / 88200 ' we assume our input sample rate is always 88200 Hz, which it of course always is
  95. NEXT i
  96.  
  97. ' initialize output file
  98. DIM zcrossprev(maxchannels) AS _FLOAT, zcrosscurr(maxchannels) AS _FLOAT, samplecount(maxchannels) AS _UNSIGNED LONG, writtensamples(maxchannels) AS _UNSIGNED LONG
  99. DIM sample(numberofinputs, samplebuffer * oversampling, maxchannels) AS INTEGER
  100. IF writedownsampled THEN
  101.     DIM outputfile(outputbuffer * maxchannels / oversampling) AS INTEGER ' the actual buffer to be written to file
  102. ELSE
  103.     DIM outputfile(outputbuffer * maxchannels) AS INTEGER ' the actual buffer to be written to file
  104. END IF
  105. DIM outputfilecopy((outputbuffer + outputfileextra) * maxchannels) AS INTEGER ' the buffer that is being filled + shifted
  106. DIM newinputsample(i, maxchannels) AS INTEGER, oldinputsample(i, maxchannels) AS INTEGER
  107. DIM bufferfull(maxchannels) AS INTEGER
  108.  
  109. ' write output to the same directory as the first input file
  110. FOR i = LEN(inputfile$(1)) TO 1 STEP -1
  111.     IF MID$(inputfile$(1), i, 1) = "\" THEN EXIT FOR
  112. NEXT i
  113. outputfile$ = MID$(inputfile$(1), 1, i) + "fmodem_out_001.wav"
  114.  
  115. filenameid = 1
  116.  
  117. ' let's see if fmodem_out_001.wav exists and increment the number in case that happens to be the case
  118. WHILE _FILEEXISTS(outputfile$)
  119.     'KILL outputfile$
  120.  
  121.     filenameoffset = LEN(outputfile$) - 6
  122.     outputfile$ = MID$(outputfile$, 0, filenameoffset)
  123.     filenameid = filenameid + 1
  124.     IF filenameid = 1000 THEN
  125.         PRINT " too many outputs, i'm broken. remove fmodem_out_xxx.wavs" ' :D
  126.         END
  127.     END IF
  128.  
  129.     filenameidstr$ = ""
  130.     IF filenameid < 10 THEN filenameidstr$ = "0"
  131.     IF filenameid < 100 THEN filenameidstr$ = filenameidstr$ + "0"
  132.     filenameidstr$ = filenameidstr$ + LTRIM$(RTRIM$(STR$(filenameid)))
  133.     outputfile$ = outputfile$ + filenameidstr$ + ".wav"
  134. WEND
  135.  
  136. outputfilenumber = FREEFILE
  137. OPEN outputfile$ FOR BINARY AS outputfilenumber
  138.  
  139. ' rewrite this shit into TYPE
  140. DIM foolong AS LONG, fooint AS INTEGER
  141. foostr$ = "RIFF": PUT outputfilenumber, 1, foostr$
  142. foostr$ = "WAVE": PUT outputfilenumber, 9, foostr$
  143. foostr$ = "fmt ": PUT outputfilenumber, 13, foostr$
  144. foolong = 16: PUT outputfilenumber, 17, foolong
  145. fooint = 1: PUT outputfilenumber, 21, fooint
  146. PUT outputfilenumber, 23, maxchannels
  147. foolong = samplerate
  148. IF writedownsampled THEN foolong = foolong / oversampling
  149. PUT outputfilenumber, 25, foolong
  150. foolong = samplerate * 2 * maxchannels
  151. IF writedownsampled THEN foolong = foolong / oversampling
  152. PUT outputfilenumber, 29, foolong
  153. fooint = 2 * maxchannels: PUT outputfilenumber, 33, fooint
  154. fooint = 16: PUT outputfilenumber, 35, fooint
  155. foostr$ = "data": PUT outputfilenumber, 37, foostr$
  156.  
  157.  
  158. LOCATE 2 + numberofinputs: PRINT " making progress"
  159. progresstimer = _FREETIMER
  160. ON TIMER(progresstimer, 0.1) GOSUB printprogress
  161. TIMER(progresstimer) ON
  162.  
  163.  
  164. ' initialize carrier gain specifics
  165. DIM carriergain(numberofinputs) AS DOUBLE, initial_carriergain(numberofinputs) AS DOUBLE, carrierlimitmin(numberofinputs) AS DOUBLE, carrierlimitmax(numberofinputs) AS DOUBLE
  166. ' range 0 -> 1 (logarithmic), 0 being no signal at all, 1 good strong signal
  167. ' magic numbers defining main signal strength
  168. carrierlimitmin(1) = 0.00035
  169. carrierlimitmax(1) = 0.01120 '0.01280
  170. carriergain(1) = RND * (carrierlimitmax(1) - carrierlimitmin(1)) + carrierlimitmin(1)
  171. ' give rest of input signals lower gain limits
  172. FOR i = 2 TO numberofinputs
  173.     carrierlimitmin(i) = 0.00001
  174.     carrierlimitmax(i) = 0.00060
  175.     carriergain(i) = RND * (carrierlimitmax(i) - carrierlimitmin(i)) + carrierlimitmin(i)
  176. NEXT i
  177.  
  178. ' update the carrier signal gains drift every 2^nth sample
  179. DIM carrierdrift_counter AS _UNSIGNED _BIT * 6
  180. 'initial_carriergain = carriergain
  181.  
  182.  
  183. DIM phaseacc(numberofinputs, maxchannels) AS _UNSIGNED LONG, frequency(numberofinputs) AS _UNSIGNED LONG
  184.  
  185. FOR i = 1 TO numberofinputs
  186.     frequency(i) = carrier / samplerate * 2 ^ 32
  187. NEXT i
  188.  
  189. DIM gain AS LONG, phasenoisegain(numberofinputs) AS LONG, phasenoisegainhalf(numberofinputs) AS LONG, noisegain AS SINGLE, noisegainhalf AS SINGLE, demodulationgain AS DOUBLE
  190. gain = 8192 ' so this this is obviously the bandwidth of the carrier signal)
  191. phasenoisegain(1) = 150 ^ 3 ' 200 ^ 3 seems to give a fair bit of noise
  192. phasenoisegainhalf(1) = phasenoisegain(1) / 2
  193. noisegain = 10: noisegainhalf = noisegain / 2
  194. demodulationgain = 32768 / gain ' plz
  195.  
  196. ' lfos not used, also in case need lfos: rewrite using arrays
  197. lfo1 = 2
  198. lfo2 = 3
  199.  
  200. lfo1rate = 10000
  201. 'lfo1gain = 1
  202. lfo2rate = 10000
  203. 'lfo1gain = 1
  204.  
  205. 'phaseacc(lfo1) = 0
  206.  
  207.  
  208. sampleoffset = 0
  209.  
  210. ' initialize zero crossings
  211. FOR channel = 0 TO maxchannels - 1
  212.     samplecount(channel) = 0
  213.     writtensamples(channel) = 0
  214.     zcrosscurr(channel) = 0 ' initial zero crossing "happened" at the first sample
  215. NEXT channel
  216.  
  217.  
  218. ' we'll just figure these out here in case we don't need to recalculate them in the loop, which we don't
  219. lutsize = audiolutsize '+ modulatoroutput ' lfo modulating audio lookup table size (should probably interpolate if ever use it)
  220. audiopalutratio = 4294967296 / lutsize ' phase accumulator size / lookup table size ratio
  221.  
  222. outputmodulatedsignal = 0
  223. IF outputmodulatedsignal = 1 THEN
  224.     IF _FILEEXISTS("signal.fm") THEN KILL "signal.fm"
  225.     OPEN "signal.fm" FOR BINARY AS 10
  226. END IF
  227.  
  228. carrierdrift_counter = 1 ' start from 1 in case we don't want to mess about with the signal gain
  229.  
  230. DIM random_drift(numberofinputs) AS DOUBLE, random_drift_multiplier AS DOUBLE
  231. random_drift_multiplier = 0.0005
  232. ' the following calculation keeps the "speed" of the occasional interference noise somewhat constant
  233. ' regardless of oversampling setting or number of input signals used by scaling the multiplier
  234. ' according to how many RND calls are being made
  235. random_drift_multiplier = random_drift_multiplier / oversampling / (2 * numberofinputs + 1)
  236.  
  237.  
  238. ' main loop
  239. FOR s = 0 TO samplelength(1) * oversampling - 1
  240.     IF s MOD samplebuffer = 0 THEN GOSUB loadsamplechunk
  241.  
  242.     carrierdrift_counter = carrierdrift_counter + 1 ' comment this out for testing without drift
  243.     IF carrierdrift_counter = 0 THEN
  244.         FOR i = 1 TO numberofinputs
  245.             random_drift(i) = (RND - 0.5) * random_drift_multiplier
  246.             carriergain(i) = carriergain(i) + random_drift(i) ' linear gain drift seems to work fine (also it's less expensive than logarithmic so win/win)
  247.             'carriergain = carriergain + (initial_carriergain - carriergain) * 0.000015 ' autotune (no, not that kind)
  248.             IF carriergain(i) < carrierlimitmin(i) THEN carriergain(i) = carrierlimitmin(i) ' hard clip the gain values
  249.             IF carriergain(i) > carrierlimitmax(i) THEN carriergain(i) = carrierlimitmax(i) ' -||-
  250.         NEXT i
  251.     END IF
  252.  
  253.     FOR channel = 0 TO maxchannels - 1
  254.         samplecount(channel) = samplecount(channel) + 1
  255.  
  256.         currentsample = s - sampleoffset + samplebuffer
  257.  
  258.         outputsample(channel) = 0
  259.         FOR i = 1 TO numberofinputs
  260.             ' following two lines are not necessary to compute in the loop since NOTHING'S MODULATING
  261.             'lutsize = audiolutsize '+ modulatoroutput ' lfo modulating audio lookup table size (should probably interpolate if ever use it)
  262.             'audiopalutratio = 4294967296 / lutsize ' phase accumulator size / lookup table size ratio
  263.             phaseacc(i, channel) = phaseacc(i, channel) + frequency(i) + sample(i, currentsample, channel) * gain + phasenoisegain(i) * RND - phasenoisegainhalf(i)
  264.             audiolutindex(i) = phaseacc(i, channel) / audiopalutratio
  265.  
  266.             outputsample(channel) = outputsample(channel) + sinelookup(lutsize, audiolutindex(i)) * carriergain(i) '+ noisegain(i) * RND - noisegainhalf(i)
  267.         NEXT i
  268.  
  269.         ' add noise to the fm'd signal
  270.         outputsample(channel) = outputsample(channel) + noisegain * RND - noisegainhalf
  271.  
  272.         IF SGN(previousoutputsample(channel)) <> SGN(outputsample(channel)) THEN GOSUB countzerocrossing
  273.  
  274.         previousoutputsample(channel) = outputsample(channel)
  275.  
  276.     NEXT channel
  277.  
  278.     ' debug, output fm'd signal
  279.     IF outputmodulatedsignal THEN PUT #10, , outputsample(0)
  280.  
  281.     IF INP(96) = 1 THEN
  282.         prematureexit = 1
  283.         EXIT FOR
  284.     END IF
  285. NEXT s
  286.  
  287.  
  288. IF prematureexit = 0 THEN
  289.  
  290.     lastchunk = 1
  291.     GOSUB writeoutputchunk
  292.  
  293.     foolong = LOF(outputfilenumber) - 8: PUT outputfilenumber, 5, foolong
  294.     foolong = LOF(outputfilenumber) - 44: PUT outputfilenumber, 41, foolong
  295.     CLOSE ' both input and output files
  296.     CLS
  297.     TIMER(progresstimer) OFF
  298.  
  299.     FOR i = 1 TO numberofinputs
  300.         LOCATE i: PRINT " input file "; trim$(i); ": "; inputfile$(i)
  301.     NEXT i
  302.     LOCATE 2 + numberofinputs: PRINT " finished everything"
  303.     LOCATE 3 + numberofinputs: PRINT " 100.00 % ("; outputfile$; ")"
  304.  
  305.     IF TIMER(.001) < time THEN time = time - 86400
  306.     executiontime = INT((TIMER(.001) - time) * 100) / 100
  307.     executiontoinputratio = INT(executiontime / sampletime(1) * 1000) / 1000
  308.     LOCATE 5 + numberofinputs: PRINT " execution time:"; executiontime; "seconds (";
  309.     IF executiontoinputratio < 1 THEN PRINT "0";
  310.     PRINT RTRIM$(LTRIM$(STR$(executiontoinputratio))); "x)"
  311.  
  312.     FOR i = 1 TO numberofinputs
  313.         LOCATE 7 + numberofinputs + i: PRINT " gain "; trim$(i); ": ";
  314.         PRINT USING "#.#########"; carriergain(i)
  315.     NEXT i
  316. ELSE
  317.     CLOSE ' both/all input and output files
  318.     'IF _FILEEXISTS(outputfile$) THEN KILL outputfile$
  319.     CLS
  320.     TIMER(progresstimer) OFF
  321.     LOCATE 1: PRINT " exited."
  322. END IF
  323.  
  324.  
  325. END
  326.  
  327.  
  328. countzerocrossing:
  329. ' calculate the instantaneous frequency from two adjacent zero crossings with linear interpolation
  330. zcrossprev(channel) = zcrosscurr(channel)
  331. zcrosscurr(channel) = s + (previousoutputsample(channel) / (previousoutputsample(channel) - outputsample(channel)))
  332. newsample = (samplerate / (2 * oversampling) / (zcrosscurr(channel) - zcrossprev(channel)) - carrieros) * demodulationgain
  333.  
  334. FOR n = 1 TO samplecount(channel)
  335.     position = writtensamples(channel) * maxchannels + channel - outputoffset
  336.     outputfilecopy(position) = newsample
  337.     writtensamples(channel) = writtensamples(channel) + 1
  338. NEXT n
  339.  
  340. IF position >= outputbuffer THEN bufferfull(channel) = 1
  341. bcount = 0
  342. FOR n = 0 TO maxchannels - 1
  343.     bcount = bcount + bufferfull(n)
  344. NEXT n
  345. IF bcount = maxchannels THEN
  346.     GOSUB writeoutputchunk
  347.     'PRINT outputoffset;
  348. END IF
  349.  
  350. samplecount(channel) = 0
  351. RETURN
  352.  
  353.  
  354. printprogress:
  355. progress$ = LTRIM$(RTRIM$(STR$(FIX(writtensamples(0) / (samplelength(1) * oversampling) * 10000))))
  356. WHILE LEN(progress$) <= 2
  357.     progress$ = "0" + progress$
  358. WEND
  359. progress$ = STRING$(5 - LEN(progress$), " ") + MID$(progress$, 1, LEN(progress$) - 2) + "." + MID$(progress$, LEN(progress$) - 1, 2) + " %"
  360. LOCATE 3 + numberofinputs: PRINT progress$
  361. RETURN
  362.  
  363.  
  364.  
  365. loadsamplechunk:
  366. FOR i = 1 TO numberofinputs
  367.     GET inputfilenumber(i), 45 + sampleoffset * numberofchannels(i) * 2 \ oversampling, dummybyte
  368.  
  369.     FOR samplex = 0 TO samplebuffer \ oversampling - 1
  370.         FOR channel = 0 TO numberofchannels(i) - 1
  371.             GET inputfilenumber(i), , newinputsample(i, channel)
  372.             FOR os = 0 TO oversampling - 1
  373.                 ' linear interpolation
  374.                 sample(i, samplex * oversampling + os, channel) = oldinputsample(i, channel) - (oldinputsample(i, channel) - newinputsample(i, channel)) * ((os + 1) / oversampling)
  375.             NEXT os
  376.             oldinputsample(i, channel) = newinputsample(i, channel)
  377.         NEXT channel
  378.     NEXT samplex
  379. NEXT i
  380.  
  381. sampleoffset = sampleoffset + samplebuffer
  382. RETURN
  383.  
  384.  
  385. writeoutputchunk:
  386. IF writedownsampled THEN
  387.     ' sample averaging downsampling method
  388.     FOR n = 0 TO outputbuffer - 1 STEP (oversampling * maxchannels)
  389.         FOR c = 0 TO maxchannels - 1
  390.             downsampled = 0
  391.             FOR m = 0 TO oversampling - 1
  392.                 ' since sample data is interlaced for all channels, we need to skip in the array to get the samples of a given channel
  393.                 downsampled = downsampled + outputfilecopy(n + c + m * maxchannels)
  394.             NEXT m
  395.             downsampled = downsampled / oversampling
  396.             outputfile(n / oversampling + c) = downsampled 'outputfilecopy(n)
  397.         NEXT c
  398.     NEXT n
  399.     PUT outputfilenumber, 45 + (outputoffset * 2 / oversampling), outputfile()
  400. ELSE
  401.     ' output oversampled data as is
  402.     FOR n = 0 TO outputbuffer - 1
  403.         outputfile(n) = outputfilecopy(n)
  404.     NEXT n
  405.     PUT outputfilenumber, 45 + (outputoffset * 2), outputfile()
  406. END IF
  407. IF lastchunk THEN RETURN
  408.  
  409. FOR n = 0 TO outputfileextra - 1
  410.     outputfilecopy(n) = outputfilecopy(outputbuffer + n)
  411. NEXT n
  412. FOR n = outputfileextra TO outputbuffer
  413.     outputfilecopy(n) = 0
  414. NEXT n
  415.  
  416. FOR n = 0 TO maxchannels - 1
  417.     bufferfull(n) = 0
  418. NEXT n
  419.  
  420. outputoffset = outputoffset + outputbuffer
  421. RETURN
  422.  
  423.  
  424.  
  425. generatelookuptable:
  426. ' look-up tables containing the sinewave period
  427. FOR size = 1 TO 1023
  428.     FOR i = 0 TO size - 1
  429.         sinelookup(size, i) = -SIN(i / size * _PI * 2) * 16384
  430.     NEXT i
  431. NEXT size
  432. RETURN
  433.  
  434. FUNCTION trim$ (i)
  435. trim$ = LTRIM$(RTRIM$(STR$(i)))
  436. END FUNCTION
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top